solve under-determined linear system by adding a constraint ?

Question

Hi all,

I have A a nxn matrix with rank(A) = n-1. Hence to solve A.u = b, I need a constraint, eg:

u1 = 0 

I would like to use a direct solver (eg: mumps). How can I tell mumps to solve for:

(A.u = b, u1 = 0)

?

Thanks a lot for you help.

Answer 1

The mathematically cleanest way is to use a lagrangian multiplier and impose that the solution has zero average.

See e.g. the singular Poisson demo:
https://fenicsproject.org/documentation/dolfin/1.6.0/python/demo/documented/singular-poisson/python/documentation.html

As an alternative, you can use a pointwise type of boundary condition. For example the script below will ensure that $u(x_0) = 0$, where $x_0 = (0,0)$.

import dolfin as dl
mesh = dl.UnitSquareMesh(10,10)
V = dl.FunctionSpace(mesh, "CG",1)
def point0(x,on_boundary):
    return x[0] < dl.DOLFIN_EPS and x[1] < dl.DOLFIN_EPS
bc_point0 = dl.DirichletBC(V, dl.Constant(0), point0, 'pointwise')

Note: $x_0$ must be a vertex of the mesh and V must be a "CG" space (of any order) for the above trick to work correctly.