Apply a geometric constraint on the boundary

variational formulation
1

Dear Fenics community,

I would like to apply a geometric constraint on a PDE defined over a 2D domain [0, R]x[0,L], namely u_1(L, x_2) = u_1(L, 0) for all x_2 . I introduced a lagrange multiplier p and wrote the following code:

  • definition of the mesh and the functionspaces:

    mesh = RectangleMesh(Point(0.0, 0.0), Point(L, R), 300, 8)

    V_element = VectorElement(‘CG’,triangle, 1)
    V0_element = FiniteElement(‘CG’,triangle, 1)
    V = FunctionSpace(mesh, V_element)
    V0 = FunctionSpace(mesh, V0_element)
    W_element = MixedElement([V_element,V0_element])
    W = FunctionSpace(mesh, W_element)

    w = Function(W)
    (u,p) = split(w)

  • element of area to include informations about the boundary condition:
    class Boundary(SubDomain):
    def inside(self, x, on_boundary):
    return x[0]>L-R/20. and on_boundary

    boundary_parts = MeshFunction(“size_t”,mesh, mesh.topology().dim()-1, 0)
    my_boundary = Boundary()
    my_boundary.mark(boundary_parts, 1)
    my_ds = ds(subdomain_data=boundary_parts)

  • expression assigned to 0 in the interior of the domain (the Lagrange multiplier p acts only on the surface):
    Exp = Expression(‘x[0] < L ? 1. : 0.0001’, degree = 1)

  • energy with Lagrange multiplier:
    bulk_energy = bulk_energy_density(u) + Exp*inner(p,p))*dx
    surface_energy = inner(p,(u[0] - u[0](L, 0)))*my_ds(3)

However, the syntax of the last line is wrong, and I wonder what the right syntax would be ?

thank you in advance for your help,

Kind regards,

Claire

2

If I understand correctly, an equivalent local constraint would be

\partial_{x_2}u_1(\mathbf{x}) = 0~~~\forall\mathbf{x}\in\{L\}\times(0,R)\text{ ,}

i.e., that u_1 doesn’t vary w.r.t. x_2 on that edge of the domain. Then you do not need to include a point evaluation of the function in your variational problem. Some other notes:

  • If you’re satisfied with a regularization of the constraint (as caused by Exp from your example code), it would be simpler to just add a quadratic penalty term to the energy. Then you don’t need to worry about adding a Lagrange multiplier (or choosing a stable space for it).

  • If you do want to use the multiplier, there are some nicer solutions for solving problems with unknown fields that are only supported on the boundary. See @francesco-ballarin’s library multiphenics or @cdaversin’s mixed dimension branch.