How to solve a nonlinear problem with global constraint using derivative?


365
views
0
11 months ago by

Hello
I would like to solve a nonlinear equation with global constraint.
As an example I am trying to solve Poisson equation with contraint that the average of the solution is 1.
F is the energy I want to minimize.
u is the variable field.
m is the global Lagrange multiplier.
I want to use derivative to compute the weak form as in the example Hyperelastic

from dolfin import *

mesh = UnitSquareMesh(50, 50)

Velt = FiniteElement("CG", mesh.ufl_cell(), 1)
Relt = FiniteElement("Real", mesh.ufl_cell(), 0) # global Lagrange multiplier
X = FunctionSpace(mesh, Velt * Relt)

f = Constant(1) # target average of u
x = Function(X)
u, m = x.split()
bc = DirichletBC(X.sub(0), Constant(0), 'on_boundary')

F = dot(grad(u),grad(u))/2*dx + m*(u-f)*dx
h = derivative(F, x, TestFunction(X))
solve(h == 0, x, bc)
print(assemble(u*dx)) # print 0, should be 1


If instead I define the weak form h by hand, using TrialFunction, it works, the average u is 1 as it should.
But in more complex cases I wish to use derivative.
Can you please tell me what I am doing wrong?
Alain

Community: FEniCS Project

2 Answers


6
11 months ago by
This works

from dolfin import *

mesh = UnitSquareMesh(50, 50)

Velt = FiniteElement("CG", mesh.ufl_cell(), 1)
Relt = FiniteElement("Real", mesh.ufl_cell(), 0) # global Lagrange multiplier
X = FunctionSpace(mesh, Velt * Relt)

f = Constant(1) # target average of u
x = Function(X)
bc = DirichletBC(X.sub(0), Constant(0), 'on_boundary')

F = dot(grad(x[0]),grad(x[0]))/2*dx + x[1]*(x[0]-f)*dx
h = derivative(F, x, TestFunction(X))
solve(h == 0, x, bc)
print(assemble(x[0]*dx)) # prints 1
1
11 months ago by
Using
u, m = split(x)​

instead of

u, m = x.split()

also works.

Please login to add an answer/comment or follow this question.

Similar posts:
Search »