Using FacetNormal (and other things only defined on the boundary) in DirichletBCs currently not possible

[florianwechsung]

The following fails

from firedrake import *
mesh = UnitSquareMesh(10, 10)
V = VectorFunctionSpace(mesh, "CG", 1)
n = FacetNormal(mesh)
u = TrialFunction(V)
v = TestFunction(V)
F = inner(grad(u), grad(v)) * dx + inner(Constant((1, 1)), v) * dx
bcs = [DirichletBC(V, n, [1, 2, 3, 4])]
solution = Function(V)
solve(lhs(F) == rhs(F), solution, bcs=bcs)

[wence-]

I think one could do this by extending ufl interpolation to support intepolating onto facets, rather than just cells. Proposed api: interpolate(expr, V, entity=ds) (or interpolate(restricted_expr, V, entity=dS) for interior facets). The way the UFL interpolation works is that we ask FIAT for the location in the reference cell of the point evaluation nodes, then compile a "point kernel" in TSFC that evaluations that expression at those points. Currently this is all hard-coded to work on cells. Instead we would need to ask FIAT for the location in the reference facets of the entities in the closure of the facet, and send that through.

[dorugeber]

Minor note: I was initially worried that DirichletBC(V, n, [1, 2, 3, 4]) doesn't behave well at corners. However, I guess it's no different to having bcs = [DirichletBC(V, 1.0, 1), DirichletBC(V, 1.5, 2)].

[wence-]

correct