Dirichlet boundary conditions and ALE.move()

The following code produces a satisfactory Taylor test:

from dolfin import *
from dolfin_adjoint import *
import numpy as np

mesh = UnitSquareMesh(10, 10)
V = FunctionSpace(mesh, "CG", 1)
VD = VectorFunctionSpace(mesh, "CG", 1)

W = Function(VD)
W.vector()[:] = .005
ALE.move(mesh, W)

u = TrialFunction(V)
v = TestFunction(V)

a = inner(grad(u), grad(v)) * dx(mesh)
l = v * dx(mesh)

u = Function(V)

dbc = DirichletBC(V, Constant(1.0), "on_boundary")

# works
# solve(a == l, u, dbc)

# works
A = assemble(a)
b = assemble(l)
dbc.apply(A, b)
solve(A, u.vector(), b)


J = assemble(u ** 2 * dx(mesh))
j = ReducedFunctional(J, [Control(W)])

H = Function(VD)
H.vector()[:] = (np.random.rand(VD.dim())*.5 - 1) * 1e-3

taylor_test(j, W, H)

Changing the Dirichlet boundary like below, yields a bad Taylor test:

dbc = DirichletBC(V, Expression("x[0]", degree=3), "on_boundary")

But with these new boundary conditions and imposing the displacements W, H to be zero on the boundary, we recover the convergence in the residuals:

mesh = UnitSquareMesh(10, 10)
V = FunctionSpace(mesh, "CG", 1)
VD = VectorFunctionSpace(mesh, "CG", 1)

W = Function(VD)
W.vector()[:] = .005
dbcv = DirichletBC(VD, Constant([0,0]), "on_boundary")
dbcv.apply(W.vector())
ALE.move(mesh, W)

u = TrialFunction(V)
v = TestFunction(V)

a = inner(grad(u), grad(v)) * dx(mesh)
l = v * dx(mesh)

u = Function(V)

dbc = DirichletBC(V, Expression("x[0]", degree=3), "on_boundary")

# works
# solve(a == l, u, dbc)

# works
A = assemble(a)
b = assemble(l)
dbc.apply(A, b)
solve(A, u.vector(), b)


J = assemble(u ** 2 * dx(mesh))
j = ReducedFunctional(J, [Control(W)])

H = Function(VD)
H.vector()[:] = (np.random.rand(VD.dim())*.5 - 1) * 1e-3
dbcv.apply(H.vector())
taylor_test(j, W, H)

From here I am tempted to conlcude that the code:

A = assemble(a)
b = assemble(l)
dbc.apply(A, b)
solve(A, u.vector(), b)

works for shape optimization, as long as:

  • the dirichlet conditions don’t depend on the mesh
  • or the boundary of the mesh isn’t moved

Could anybody perhaps further explain what’s going on and in case confirm/deny that my observation is correct? Thanks in advance!

Addendum

… and now imposing Neumann conditions on part of the boundary, in the last code snippet, yields the residuals: [1.9974521264767453, 1.9355106278012318, 1.895501764619647], which are close to 2 but decreasing, is this just numerical noise?

There is no way of differentiating Expressions. If you instead project spatial coordinates into the function space of the boundary condition and use this function in the boundary condition, the Taylor test should pass. (You could Also impose the DirichletBC weakly with Nitsches method).

Doing Taylor tests without mesh smoothing can affect the convergence rate, see section 5.1 of: https://arxiv.org/pdf/2001.10058.pdf

1 Like

Thanks a lot. The article is also very interesting: about section 5.1.1, are you aware by chance of any code out there demonstrating the implementation of a custom Riesz map in the context of shape optimization?