[petsc-dev] Discussion about time-dependent optimization moved from PR

[Jed Brown]

Barry Smith <bsmith at mcs.anl.gov> writes:

>> On Oct 17, 2017, at 10:00 AM, Jed Brown <jed at jedbrown.org> wrote:
>> 
>> Barry Smith <bsmith at mcs.anl.gov> writes:
>> 
>>>> On Oct 17, 2017, at 9:20 AM, Jed Brown <jed at jedbrown.org> wrote:
>>>> 
>>>> Barry Smith <bsmith at mcs.anl.gov> writes:
>>>> 
>>>>>> On Oct 17, 2017, at 7:43 AM, Stefano Zampini <stefano.zampini at gmail.com> wrote:
>>>>>> 
>>>>>> 
>>>>>> 
>>>>>> 
>>>>>>  What happens if you call TSCreateAdjointsTS() on a TS obtained with TSCreateAdjointsTS()? Is the resulting TS useful?
>>>>>> 
>>>>>> 
>>>>>> I don't think so.
>>>>> 
>>>>>  Hmm, this is rather bothersome. So the new TS that comes out is not
>>>>>  really a full-featured TS (using your new extended TS object that
>>>>>  knows about adjoints)? 
>>>> 
>>>> Is it useful to precondition Newton with Newton preconditioned with
>>>> Newton?  
>>>> 
>>>>>  It is like a subclass of TS that can only do timestepping? This
>>>>>  goes back to my concern that TS shouldn't be overloaded with all
>>>>>  the adjoint sensitivity stuff, there should be a higher level
>>>>>  PetscSensi object that does sensitivities and it contains multiple
>>>>>  TS for the ODE/DAE integrator and one for "reverse" mode.
>>>> 
>>>> Let's not mix together Stefano's TSComputeObjectiveAndGradient (which
>>>> could easily go into a helper PetscSensi) with the control flow needed
>>>> to create/integrate an adjoint system.
>>> 
>>>   Jed,
>>> 
>>>    I am not mixing them!  I am not even considering the TSComputeObjectiveAndGradient() at all I am only considering
>>> the "the control flow needed to create/integrate an adjoint system." My PetscSensi is not for TSComputeObjectiveAndGradient level stuff.
>> 
>> Okay, interesting.
>> 
>>> 
>>>> In my opinion, it must be
>>>> possible to have a comparable/same interface to create a discrete
>>>> adjoint (which depends on the forward implementation) as to create a
>>>> continuous adjoint (with customizable space and time discretization).
>>> 
>>>   Absolutely and fundamentally. No argument at all about this.
>>> 
>>>> 
>>>> I think we should postpone discussion of higher level convenience
>>>> functions like TSComputeObjectiveAndGradient or some new organization
>>>> thereof until we agree about the basic representation of adjoints.
>>> 
>>>   Totally agree and that is what I have done. 
>>> 
>>>> My
>>>> preference is for TSCreateAdjointTS() which the user can customize to
>>>> use arbitrary time integration methods (thus becoming a
>>>> "continuous-in-time" adjoint) and to use different spatial
>>>> discretizations (for consistency with the adjoint PDE).
>>> 
>>>   I understand, I think TSCreateAdjointTS() is not the right model. I think that the correct model is, for example
>>> 
>>> PetscSensiCreate()
>>> PetscSensiSetObjective.... or whatever etc etc
>>> PetscSensiGetIntegrator(,&TS)
>>> PetscSensiGetAdjointIntegrator(,&TS).
>> 
>> How would PetscSensiGetAdjointIntegrator be implemented?  Since it needs
>> to be able to create a discrete adjoint (which depends on the forward TS
>> implementation), the control flow *must* go through TS.  What does that
>> interface look like?
>
>    I am not sure what you mean. 
>
>    PetscSensi will contain the information about objective function etc as well as a pointer to the trajectory history and (of course) a pointer to the original TS. All the  auxiliary vectors etc related to sensitivities. 
>
>    The adjoint TS that is returned has the calls need to get access to the forward solution that it needs etc from the PetscSensi object (which likely will call the forward TS for missing timestep values etc. But the calls go through the PetscSensi object from the adjoint TS, they don't go directly to the forward TS.)

To be concrete with TSRK, PetscSensi does not have access to the Butcher
table for the forward method because that only exists in
src/ts/impls/explicit/rk/rk.c.  Since the discrete adjoint depends on
this detail of the forward method, there MUST be some code somewhere
that depends on the forward Butcher table.  Since that data ONLY EXISTS
in rk.c, whatever it is that PetscSensi does, there must be some code in
rk.c to define the adjoint method.

What is that code and what is the call stack to reach it?

>    Both proposed models shove one heck of a lot of stuff into TS; I think it would be a better design and better user API if we don't do that; we reserve the TS for being able to integrate something (and that is about all). Maybe it is not possible but I don't see why it is not possible.
>
>    Since I don't know the details I asked for the people who do know the details (Hong and Stefano) to try to develop such an API (or clearly articular why it is impossible.)