Oh, I missed something in all that: not only do we have yet to
automate \partial R / \partial p, but we leave it up to the user to
evaluate the quantity of interest derivative \partial Q / \partial u.
Trying to automate that is uglier because, unlike R, in general we
can't necessarily express Q as a sum of Q_e on each element. That
means finite differencing or automatic differencing is probably out of
the question, and that means there's probably nothing we can do to
make the users' task any simpler. I'd be happy to be contradicted,
though.

Roy
On Wed, 6 May 2009, Roy Stogner wrote:
> Copying this to the devel list since it's something a few of us might
> want to talk about. In particular I'd be happy to hear any better
> ideas for automatically computing \partial R / \partial p. A
> finitedifferenced residual would work for FEMSystem and
> NonlinearImplicitSystem, but LinearImplicitSystem doesn't even let you
> assemble the residual without getting a big fat matrix along with it.
> 
> Roy
>
>  Forwarded message 
> Date: Wed, 6 May 2009 11:42:45 0500 (CDT)
> From: Roy Stogner <roystgnr@...>
> ReplyTo: Roy Stogner <roy@...>
> To: Andrea Hawkins <andjhawkins@...>
> Subject: Re: libMesh Adjoint Capabilities
>
>
> On Wed, 6 May 2009, Andrea Hawkins wrote:
>
>> Dr. Oden told me that in the PECOS meetings earlier this week you
>> had mentioned something about libMesh having some remarkable adjoint
>> problem capabilities.
>
> Brand new development, but starting to get there.
>
>> From what he understood, if you had the forward problem R(u,v,p) =
>> 0, libMesh was somehow able to formulate \partial R/ \partial u and
>> \partial R / \partial p to set up the adjoint problem and the
>> sensitivities. Is this correct?
>
> That's the whole goal, but it's only about 25% correct, so far.
>
> The good news:
>
> In the FEMSystem framework you often don't even to write all of
> R(u,v,p) yourself  you can just give the weighted residual parts from
> du/dt=F(u,v) and G(u,v)=0, and if you're using conservation variables
> and unstabilized Galerkin, the library can do theta method Euler or
> trapezoidal time integration to give you R. Nonlinear or stabilized
> mass matrices are a little trickier but still doable. This is what
> we've started testing; Vikram claims that FEMSystem::adjoint_solve()
> gave him a verified solution for his quantity of interest.
>
> FEMSystem will also automatically do elementwise finitedifferenced
> Jacobians to give you \partial R / \partial u; that's been pretty well
> tested for years.
>
> For any ImplicitSystem based code, libMesh will do the Jacobian
> transpose automatically  take a System for which you've already got a
> good matrix/rhs assembly routine, attach a quantityofinterest
> assembly routine which replaces the rhs , and call
> System::adjoint_solve(), and it should give you the right output.
>
> The bad news:
>
> If you want to solve a discretizationoftheadjoint instead of
> adjointofthediscretization problem, libMesh can't help; best to try
> to make those things equivalent in the first place.
>
> The adjoint solution is done on your current mesh  if you need a
> finer adjoint solution for goaloriented error estimation then you
> still have to do the refinement yourself after your solve() and before
> your adjoint_solve(). That may not change soon; there's no telling
> which users will want to do p refinement vs. h refinement, use lower
> tolerances on a fine mesh, etc.
>
> We don't do adjoints for transient or decoupled problems yet. I
> probably won't add the former this year or the latter ever.
>
> Adjoint solve semantics may change soon, and frequently  we'll
> probably soon pick an "official" additional vector to store adjoint
> solutions in instead of overwriting System::solution, for an immediate
> example.
>
> This stuff is bleeding edge. For instance there's separate
> adjoint_solve() implementations for LinearImplicitSystem,
> NonlinearImplicitSystem, and FEMSystem, but only the latter has been
> tested.
>
> The worst news:
>
> Right now users have to assemble \partial R / \partial p by hand; I
> haven't even yet decided what the best way is to give the library
> access to p. My best guess at a first cut would be a
> finitedifferenced residual based on automatically adjusted
> EquationSystems::parameters values.
> 
> Roy
>
> 
> The NEW KODAK i700 Series Scanners deliver under ANY circumstances! Your
> production scanning environment may not be a perfect world  but thanks to
> Kodak, there's a perfect scanner to get the job done! With the NEW KODAK i700
> Series Scanner you'll get full speed at 300 dpi even with all image
> processing features enabled. http://p.sf.net/sfu/kodakcom
> _______________________________________________
> Libmeshdevel mailing list
> Libmeshdevel@...
> https://lists.sourceforge.net/lists/listinfo/libmeshdevel
>
