[firedrake] Problem with Jacobian

[Buesing, Henrik]

> -----Ursprüngliche Nachricht-----
> Von: firedrake-bounces at imperial.ac.uk [mailto:firedrake-
> bounces at imperial.ac.uk] Im Auftrag von Stephan Kramer
> Gesendet: 05 November 2015 11:27
> An: firedrake at imperial.ac.uk
> Betreff: Re: [firedrake] Problem with Jacobian
> 
> On 05/11/15 10:08, Lawrence Mitchell wrote:
> >
> >> On 5 Nov 2015, at 09:49, Buesing, Henrik <HBuesing at eonerc.rwth-
> aachen.de> wrote:
> >>
> >> Dear all,
> >>
> >> I'm having a variable Sw, which I calculate pointwise in a routine
> calc_Sw (see attachment). This variable depends on my primary unknown h:
> Sw = (h-hn)/(hw-hn) (hn,hw known values).
> >>
> >> But now it seems like automatic differentiation for this routine does
> not work. I'm getting zero entries for the Jacobian, whereas d(Sw)/dh =
> 1.0 should hold.
> >
> > Yes, this is because the AD doesn't know about the relationship
> between Sw and dh.  UFL has a facility for this, but I notice we don't
> expose it in firedrake (however, it is straightforward to do):
> >
> >
> > We want
> >
> > derivative(F, u)
> >
> > But F contains a coefficient, S, whose derivative wrt u is 1.0
> (however, they are not symbolically related in a way UFL understands).
> So we build a mapping from this coefficient to its derivative wrt u:
> >
> > coefficient_derivatives = {S: 1.0}
> >
> > and then pass this additional information to the derivative call.
> >
> > derivative(F, u, coefficient_derivatives=coefficient_derivatives)
> >
> > Firedrake uses the UFL derivative function, but does not expose this
> extra argument in the interface.  It is straightforward to alter the
> definition in firedrake/ufl_expr.py to take this extra argument and pass
> it through.  If this works for you, do you want to propose a patch that
> adds this functionality?
> >
> > Cheers,
> >
> > Lawrence
> >
> 
> Are you sure the problem isn't that he needs bounds for his solve?
> AFAICS dS/dh=0 for h<hw or h>hn. Otherwise shouldn't rewriting this with
> conditionals have worked? 
[Buesing, Henrik] I failed to rewrite my function with conditionals. Just because I didn't know how to properly write down the conditional() stuff and how to set Sw. 
If I could rewrite this function, I assume this would work!

> Henrik, I presume you are looking for a
> solution hw<h<hn ?
[Buesing, Henrik] Yes. 
> 
> Cheers
> Stephan
> 
> 
> 
> _______________________________________________
> firedrake mailing list
> firedrake at imperial.ac.uk
> https://mailman.ic.ac.uk/mailman/listinfo/firedrake