I am trying to see if I can replace COMSOL with Fenics for doing multiplysics simulations using FEM. A strong point of COMSOL is the ability to include constraints using Lagrange multipliers
Using Fenics, I have been able to impose global constraint (integral form), but not local (point-wise) constraint.
Let us consider a classical problem (solved with COMSOL)
Find u in [0,1]x[0,1] such that -Delta u = f_0 where f_0 is such that u(0.5,0.5)=1
f_0 is the the Lagrange multiplier associated with the point-wise constraint u(0.5,0.5)=1.
If I impose f_0 (as a Constant("2.67")) , I can solve the problem with fenics using a point-wise integration (with Measure('dP')) assuming (0,5,0,5) is a mesh vertex, but it does not work if a use a mixed formulation with Lagange multiplier (using FunctionSpace(mesh,"R",0)). Here is the Python code
mesh=UnitSquareMesh(4,4)
class MidlePoint(SubDomain):
def inside(self,x,on_boundary):
return near(x[0],0.5) and near(x[1],0.5)
# boundary
gamma0=CompiledSubDomain("on_boundary")
midlepoint=MidlePoint()
# points wise integration
points=VertexFunction("uint",mesh)
points.set_all(0)
midlepoint.mark(points,5)
dp = Measure('dP')[points]
# FE space
V = FunctionSpace(mesh,"Lagrange",1)
R = FunctionSpace(mesh,"R",0)
W = MixedFunctionSpace([V,R])
# FEM formulation
u,lb = TrialFunction(W)
v,lb_test = TestFunction(W)
a = inner(nabla_grad(u),nabla_grad(v))*dx() + lb*v*dp(5) + lb_test*u*dp(5)
L = lb_test*dp(5)
# solve
bcs = [DirichletBC(W.sub(0),0.0,gamma0)]
wh=Function(W)
solve(a==L,wh,bcs)
running the code on fenics 1.5 I get the following error on solve
<ipython-input-197-1e95aec09e42> in <module>()
wh=Function(W)
----> solve(a==L,wh,bcs)
/usr/lib/python2.7/dist-packages/dolfin/fem/solving.pyc in solve(*args, **kwargs)
267 # Call variational problem solver if we get an equation (but not a tolerance)
268 elif isinstance(args[0], ufl.classes.Equation):
--> 269 _solve_varproblem(*args, **kwargs)
/usr/lib/python2.7/dist-packages/dolfin/fem/solving.pyc in _solve_varproblem(*args, **kwargs)
296 solver = LinearVariationalSolver(problem)
297 solver.parameters.update(solver_parameters)
--> 298 solver.solve()
RuntimeError
*** -------------------------------------------------------------------------
*** Error: Unable to assemble form over vertices.
*** Reason: Expecting test and trial spaces to only have dofs on vertices for point integrals.
*** Where: This error was encountered inside Assembler.cpp.
*** Process: unknown
***
*** DOLFIN version: 1.5.0
*** Git changeset: unknown
*** -------------------------------------------------------------------------
Looking at the code in Assembler.cpp, the error occurs in a test to Check that the test and trial space as dofs on the vertices with num_entity_dofs(0)#0
But lb, the lagrange multiplier is in the space R (the real constant) and has a num_entity_dofs(0)=0
But as lb is constant it can be evaluated at any vertices,and thus the assembly should be possible.
Am I doing wrong ?
Can I solve problem with pointwise constraint using lagrange multiplier with Fenics ?
I know that for this particular problem I can impose u(0.5,0.5)=1 using DirichletBC pointwise
but this is not a solution for a general point-wise constraint problem
Best Regards
Marc
