How to solve linear system with constaint?

Question

I want to solve a singular linear system, but it is solvable with a constraint, like zero mean value.

I know this. But the code can not run in my computer, something run with solver.set_nullspace(null_space).
What's more, if I solve the singular poisson problem, how do I know basis of null_space?

My code is below:

from fenics import *

mesh = UnitCubeMesh(10, 10, 10)
DG = FunctionSpace(mesh, 'DG', 0)
RT = FunctionSpace(mesh, 'RT', 1)

p = TrialFunction(DG)
v = TestFunction(RT)

bc = DirichletBC(RT, Constant([0.0, 0.0, 0.0]), DomainBoundary())

a = p * div(v) * dx
L = inner(f, v) *dx

A, b = assemble_system(a, L, bc)

Then how do I solve the problem? And how do I use zero mean value( i.e. p*dx = 0)?

I can solve Au=b out of fenics, like using numpy, but it's strange.

Answer 1

To use the zero mean value constraint, try to add a function space for Lagrange multipliers

R = FunctionSpace(mesh, 'R', 0)

then

r  = TrialFunction(R)
r_ = TestFunction(R)

and

a += r_*p*dx + p_*r*dx

where p_ is the corresponding test function.