Taylor test for time dependent problem
Issue #81
new
Hello, I am trying to conduct a Taylor test for a simple time dependent problem. The following code appears to work:
from dolfin import *
from dolfin_adjoint import *
n = 30
mesh = UnitSquareMesh(n, n)
V = VectorFunctionSpace(mesh, "CG", 2)
ic = project(Expression(("sin(2*pi*x[0])", "cos(2*pi*x[1])"), degree=2), V)
def main(nu):
u = ic.copy(deepcopy=True)
u_next = Function(V)
v = TestFunction(V)
timestep = Constant(0.005)
F = (inner((u_next - u)/timestep, v)
+ inner(grad(u_next)*u_next, v)
+ nu*inner(grad(u_next), grad(v)))*dx
bc = DirichletBC(V, (0.0, 0.0), "on_boundary")
t = 0.0
Jnu = 0.0
end = 0.1
times = [t]
adj_start_timestep(time=t)
while (t <= end):
solve(F == 0, u_next, bc)
u.assign(u_next)
Jnu += assemble(inner(u,u)*dx)
t += float(timestep)
times.append(t)
adj_inc_timestep(time=t, finished=t>end)
return u,times, Jnu
if __name__ == "__main__":
nu = Constant(0.0001)
u,times, Jnu = main(nu)
Jfull = sum([inner(u,u)*dx*dt[t] for t in times])
J = Functional(Jfull)
#J = Functional(inner(u, u)*dx*dt)
dJdnu = compute_gradient(J, Control(nu))
print "Gradient is: ", float(dJdnu)
parameters["adjoint"]["stop_annotating"] = True # stop registering equations
def Jhat(nu): # the functional as a pure function of nu
u,times, Jnu = main(nu)
return Jnu
conv_rate = taylor_test(Jhat, Control(nu), Jnu, dJdnu)
However, if I change the functional to the following
J = Functional(inner(u, u)*dx*dt)
and change Jnu to this:
Jnu += assemble(timestep*inner(u,u)*dx)
I do not get the expected rate of convergence. This is a pretty simple integration scheme, but I thought it would work if the timestep was small. How should I be calculating Jnu?
Thank You!
Jacob
Hi Jacob,
The mistake is that the
dt
integral implements a trapezoidle rule, whic you implemented a right Riemann sum. In other words, you need to update yourJnu
calculation to get the second order convergence again.Best wishes,
Simon