# derivatives with respect to boundary conditions

Is dolphin_adjoint capable of taking derivatives with respect to the constants that appear in Dirichlet boundary conditions? The code below, which produces the output:

Solving linear variational problem.
Solving linear variational problem.
Solving linear variational problem.
Solving linear variational problem.
------------------------------------------
f = 0.333333 (exact value is 0.333333)
df/da = None (exact value is 0.5)
df/db = None (exact value is 0.5)
df/dc = -0.0833333333333 (exact value is -0.083333)
------------------------------------------

is an example problem that illustrates my attempt to calculate such derivatives.

# The solution to the ODE u" = c with u(0) = a, u(1) = b is
# u(x) = (c/2) x(x-1) + b x + a (1-x).
# Integrating u from 0 to 1 gives -c/12 + b/2 + a/2.
# So the partial derivative of this integral with respect to c is -1/12,
# the partial derivative with respect to b is 1/2, and the partial
# derivative with respect to a is 1/2.
# This script uses dolfin-adjoint to try to find these derivatives.
# Written by Glen D. Granzow on September 19, 2012.

from dolfin import *

if __name__ == '__main__':

a = Constant(0.0)
b = Constant(1.0)
c = Constant(2.0)

mesh = UnitInterval(10)
functionSpace = FunctionSpace(mesh, 'Lagrange', 2)
u = TrialFunction(functionSpace)
v = TestFunction(functionSpace)
solution = Function(functionSpace)

bc = [DirichletBC(functionSpace, a, 'on_boundary && near(x[0],0)'),
DirichletBC(functionSpace, b, 'on_boundary && near(x[0],1)')]

solve(-u.dx(0)*v.dx(0)*dx == c*v*dx, solution, bc)

f = Functional(solution*dx*dt[FINISH_TIME])

print '------------------------------------------'
print 'f = %f (exact value is %f)' % (assemble(solution*dx), -float(c)/12 + float(b)/2 + float(a)/2)
print 'df/da = %s (exact value is 0.5)' % dfda
print 'df/db = %s (exact value is 0.5)' % dfdb
print 'df/dc = %s (exact value is %f)' % (dfdc, -1.0/12.0)
print '------------------------------------------'

The problem is that in the output "df/da = None" and "df/db = None" instead of the anticipated numerical values.

## Question information

Language:
English Edit question
Status:
Solved
For:
Assignee:
No assignee Edit question
Solved by:
Glen D. Granzow
Solved:
Last query:
 Revision history for this message Patrick Farrell (pefarrell) said on 2012-11-20: #1

Hi Glen,

That currently isn't possible, because I need to have a small amount of dolfin functionality exposed that I don't currently know how to do. However, differentiating with respect to boundary conditions is very important, and it's high on my todo list -- hopefully I will have it done in January.

In the meantime, Simon Funke has mentioned in a previous question how to do it if your BC value is weakly enforced (I don't know if that's feasible for you): https://answers.launchpad.net/dolfin-adjoint/+question/213089

Thanks for your interest, and I'll let you know when that feature is available.

Patrick

 Revision history for this message Glen D. Granzow (ggranzow) said on 2012-11-20: #2

Thank you for your response Patrick,

I'm not sure if what I have done is what you meant by weakly enforcing my boundary condition but the code below produces the desired output:

Solving linear variational problem.
Solving linear variational problem.
Solving linear variational problem.
Solving linear variational problem.
------------------------------------------
f = 0.333333 (exact value is 0.333333)
df/da = 0.5 (exact value is 0.5)
df/db = 0.5 (exact value is 0.5)
df/dc = -0.0833333333333 (exact value is -0.083333)
------------------------------------------

for my example problem.

# The solution to the ODE u" = c with u(0) = a, u(1) = b is
# u(x) = (c/2) x(x-1) + b x + a (1-x).
# Integrating u from 0 to 1 gives -c/12 + b/2 + a/2.
# So the partial derivative of this integral with respect to c is -1/12,
# the partial derivative with respect to b is 1/2, and the partial
# derivative with respect to a is 1/2.
# This script uses dolfin-adjoint to find these derivatives.
# Written by Glen D. Granzow on November 20, 2012.

from dolfin import *

if __name__ == '__main__':

a = Constant(0.0)
b = Constant(1.0)
c = Constant(2.0)

mesh = UnitInterval(10)
functionSpace = FunctionSpace(mesh, 'Lagrange', 2)
u = TrialFunction(functionSpace)
v = TestFunction(functionSpace)
solution = Function(functionSpace)

boundaries = MeshFunction("uint", mesh, mesh.topology().dim()-1)
boundaries.set_all(99)

class LeftBoundary(SubDomain):
def inside(self, x, on_boundary):
return on_boundary and near(x[0],0)

class RightBoundary(SubDomain):
def inside(self, x, on_boundary):
return on_boundary and near(x[0],1)

LeftBoundary().mark(boundaries, 0)
RightBoundary().mark(boundaries, 1)
ds = ds[boundaries]

solve(-u.dx(0)*v.dx(0)*dx + (u+u.dx(0))*v*ds(1) + (u-u.dx(0))*v*ds(0) == c*v*dx + a*v*ds(0) + b*v*ds(1), solution)

f = Functional(solution*dx*dt[FINISH_TIME])

print '------------------------------------------'
print 'f = %f (exact value is %f)' % (assemble(solution*dx), -float(c)/12 + float(b)/2 + float(a)/2)
print 'df/da = %s (exact value is 0.5)' % dfda
print 'df/db = %s (exact value is 0.5)' % dfdb
print 'df/dc = %s (exact value is %f)' % (dfdc, -1.0/12.0)
print '------------------------------------------'

 Revision history for this message Simon Funke (simon-funke) said on 2012-11-21: #3

Hi Glen,

Basically, dolfin-adjoint can currently differentiate only with respect to
Parameters that appear in the UFL form which is passed to the solve
routine.
This is the case in your last code example, that is why it works for you now.

Looking at your first code example that you provided, I saw that the parameters a, b and c
don't show up in the weak form in your equations. Hence differentiating with respect to these