I’m trying to implement a stochastic filter for a system based on a time-dependent PDE (heat equation). The PDE depends on an unknown potential that I eventually want to estimate. For this, I need the derivative of the solution at a time t with respect to the potential. How do I get this?
I know that what I am looking for is basically the solution to the tangent linear equation, but I don’t want to implement this by myself, because I assume that dolfin-adjoint is able to do this much better.
from fenics import * from fenics_adjoint import * T = 50.0 #time of evaluation N_t = 200 #number of time steps dt = T / N_t #time step size N_s = 64 #space mesh size # Set initial value u_0 = Expression('0', degree=1) # Set boundary value u_D = Expression('0', degree=1) # Basic setup mesh = UnitIntervalMesh(N_s) V = FunctionSpace(mesh, 'P', 2) def boundary(x, on_boundary): return on_boundary bc = DirichletBC(V, u_D, boundary) # Define the potential q = interpolate(Expression('0.1-0.05*sin(2*pi*(x-0.25))', degree=1), V) control = Control(q) # Want to differentiate with respect to q # Define the variational problem (using backwards Euler) u_next = Function(V) u = interpolate(u_0,V) v = TestFunction(V) r = Expression('x >= 0.215 && x <= 0.315 ? (4*sin(4*pi*t)+0.001*pow(t,2)) : 0', \ degree=1, t=0) F = ( (u_next - u)*v /dt + dot(q*grad(u_next), grad(v)) - r*v )*dx # Start solution for t in range(N_t): # Update current time t += dt # At time n, we want r = r((n+1)dt) r.t = t # Compute solution solve(F == 0, u_next, bc) # Update u u.assign(u_next) # Now, I want the derivative of u(T) with respect to q dudq = compute_gradient(u, control) # But obviously, this does not work.