Issue solving the Poisson equation on a 3D surface mesh

Errors
1

Hi, I’m trying to solve the Poisson equation on a 3D surface (here testing just a unit cube). Following this presentation for the moebius strip I think the code below should work. But I either get errors like “Warning: Found no facets matching domain for boundary condition.” or “ArityMismatch: Adding expressions with non-matching form arguments (‘v_0’,) vs (‘v_0’, ‘v_1’)”. The minimal non-working example is below:

mesh = BoundaryMesh(UnitCubeMesh(1, 1, 1), 'exterior')
# mesh = Mesh('unit_cube_shell.xml')

V = FunctionSpace(mesh, 'Lagrange', 1)
u = TrialFunction(V)
v = TestFunction(V)

a = inner(grad(u), grad(v))*dx
f = Constant(1)
L = f*v*dx

u = Function(V)
# bc = DirichletBC(V, Constant(0), "on_boundary")
bc = DirichletBC(V, Constant(0), \
   lambda x, on_boundary: on_boundary and near(x[2], 0))
solve(a == L, u, bc)
# solve(a - L == 0, u, bc)

Maybe Fenics has changed since that presentation from '13? Is there a way to do this?

2

Hi @patrick-kwok,

Your domain is a closed surface, and therefore all of the mesh facets are interior facets. Thus on_boundary is False for every facet in your mesh. This is why your DirichletBC finds no matching facets.

Modifying your boundary condition as follows takes care of the warning

bc = DirichletBC(V, Constant(0), lambda x, on_boundary: near(x[2], 0))
1 Like
3

Yes, of course, thank you now its so obvious.