Hi, I am solving an optimization problem on the unit square, with the control being a scalar function V \in \mathbb{R}. For this problem, V needs to be constant in the y direction: V= V(x). (This can be interpreted as the thickness of a sheet of material, which varies only with x.) It has become apparent I need to discard the redundant degrees of freedom. To enforce this, I tried solving a simple PDE (prior to the main PDE) of the form \int \frac{\partial{u}}{\partial{y}}v=0, with boundary conditions specified only on the bottom portion of the domain y=0. (This seemed simpler than trying to enforce a variational equality.) I cannot get the taylor test to pass. Iād be grateful for any suggestions on how to get this to work or alternative ways to enforce this type of constraint, thanks!
from dolfin import *
from dolfin_adjoint import *
mesh = UnitSquareMesh(50, 50)
V = FunctionSpace(mesh, "CG", 2)
bd_thickness = Function(V, name="boundary_thickness")
bd_thickness.assign(interpolate(Expression('x[0]', degree=1), V))
bd = CompiledSubDomain("near(x[1],0) && on_boundary")
bc = DirichletBC(V, bd_thickness, bd)
u_ = TrialFunction(V)
v = TestFunction(V)
j_hat = Constant(('0', '1'), name="j_hat")
a = inner(inner(grad(u_), j_hat), v)*dx
L = v*Constant(0, name="zero")*dx
u = Function(V, name="u")
solve(a == L, u, bc)
plot(u)
This seems to be correct: u varies with x only, and is described only by the values of bd_thickness along the bottom boundary y=0.
J = assemble(u*dx)
dJdnu = compute_gradient(J, Control(bd_thickness))
plot(dJdnu)
The gradient appears constant along y.
perturbation = interpolate(Expression("0.01*sin(x[0]*pi)", pi=pi, degree=2), V)
Jhat = ReducedFunctional(J, Control(bd_thickness))
conv_rate = taylor_test(Jhat, bd_thickness, perturbation)
tape = get_working_tape()
tape.visualise()
The taylor test fails. However, the tape looks correct.