Time-dependent expressions and preassembled solves interact in a funny way

[Patrick Farrell]

When the script below is run the replay_dolfin() function fails, but if the time-stepping information is added then the replay works! Doing the same thing with a series of variational problem solves without the timestepping information is ok as well. Is there something special about the matrix vector solves that requires them to always need time-stepping information?

from dolfin import *
from dolfin_adjoint import *

mesh = UnitSquareMesh(4,4)
V = FunctionSpace(mesh, "CG", 1)

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

a = inner(grad(u), grad(v))*dx
u = Function(V, name = "solution")

bc = DirichletBC(V, 0, "on_boundary")

source = Expression("t*sin(x[0])*sin(x[1])", t = 0.0)
f = source*v*dx

t = 0.0
dt = 0.1

A = assemble(a)
#adj_start_timestep()
for i in range(2):
    t += dt
    source.t = t
    F = assemble(f)
    bc.apply(A)
    bc.apply(F)
    solve(A, u.vector(), F, "cg", "ilu")
    #adj_inc_timestep(t)

adj_html("forward.html", "forward")
assert(replay_dolfin())

[Patrick Farrell (pefarrell)]

Interestingly, uncommenting those lines doesn't change that the replay fails, for me.

[Patrick Farrell (pefarrell)]

Hi Gabriel,

It was a race condition, well spotted!

dolfin-adjoint needs to bundle state (such as values of Expression parameters) with each equation. To do that, it takes a hash. If two equations had the same hash, then they would both fetch the same values for expression parameters, even if they should be different. To that end, I included time.time() in the hash -- but your example was so small that time.time() would (usually) return the same value, the expression parameters would be incorrect, and the replay would fail. Hmm!

That's why adding the adj_inc_timestep calls made a difference for you -- it made it so that time.time() returned something different.

I've fixed this by hashing with random.random() instead, in r726.

Patrick

[Gabriel Balaban (gabrib)]

Hi Patrick!
Great fix! I suspected that what I was seeing in the code example was a "smell" indicating some other issue. It's great to hear that you found the cause.
  Cheers,
   Gabriel