Accessing Tape

Asked by lvleph

I am trying to access the tape so that I can create error indicators. However I am having trouble, because in the tape there are iterations and timesteps. How can I access just the last iteration of a timestep? In the following code snippet I compute the adjoint and then look for the corresponding adjoint solution. Once I have the adjoint solution I try to use it and the tape value to create my error indicators, but as I stated some of the vars in compute_adjoint are for an iteration. This makes it difficult for me to just step through and assume what I have corresponds to a timestep rather than an iteration for a timestep.

        i = int(math.ceil(T/k)) - 1 #last time step

        adjoint = compute_adjoint(J,forget=False) #adjoint solution
        for (phi, var) in adjoint: #step through the adjoint solution
            if var.name == 'w':
                # Compute error indicators ei
                wtape = DolfinAdjointVariable(w).tape_value(timestep=i)
                wtape_ = DolfinAdjointVariable(w_).tape_value(timestep=i)
                LR1 = k*self.weak_residual(W, wtape, wtape_, phi, ei_mode=True)
                ei.vector()[:] += assemble(LR1,annotate=False).array() #error indicators
                i -= 1 #decrease the timestep

Question information

Language:
English Edit question
Status:
Solved
For:
dolfin-adjoint Edit question
Assignee:
No assignee Edit question
Solved by:
lvleph
Solved:
Last query:
Last reply:
Revision history for this message
lvleph (erichlf) said :
#1

I had to look into the dolfin-adjoin and libadjoint source to figure this out. I ended up doing it in two steps.

        timestep = None
        wtape = []
        phi = []

        adjoint = compute_adjoint(J,forget=False) #adjoint solution
        for (adj, var) in adjoint:
            if var.name == 'w' and timestep != var.timestep:
                timestep = var.timestep
                # Compute error indicators ei
                wtape.append(DolfinAdjointVariable(w).tape_value(timestep=timestep))
                phi.append(adj)

        for i in range(0, len(wtape)-1):
            LR1 = k*self.weak_residual(W, wtape[i], wtape[i+1], phi[i], ei_mode=True)
            ei.vector()[:] += assemble(LR1,annotate=False).array()

Revision history for this message
Patrick Farrell (pefarrell) said :
#2

Hi Erich,

Glad to hear you figured it out; it wasn't obvious!

Best regards,

Patrick

Revision history for this message
lvleph (erichlf) said :
#3

I think this is my biggest complaint about DOLFIN-Adjoint, and I think I seen others make the same comment; it is a bit too black box. Maybe the documentation could be better, or maybe better examples? Most of the examples I found were for testing and certainly none of them attempted to access the tape or build error indicators. I certainly don't mind contributed to an example, once I have it completely worked out.

Revision history for this message
Patrick Farrell (pefarrell) said :
#4

Yep, it's something we're quite aware of -- we added a few more examples a few weeks a go, with more on the way. In general I agree with you that there's too much magic, and there will be significant redesigns (of time-dependent functionals, and how the controls are handled) over the next few months to make things a bit clearer and more explicit.

It would be great if you would contribute to an example (especially on error estimation); I'd be delighted to see it.