[dune-pdelab] Why the lambda_*** methods?

[Steffen Müthing]

Hi Oliver,

writing as someone who wasn’t around when PDELab was originally written,
but as one of the people who know the kind pretty well right now:

> Am 12.05.2017 um 22:14 schrieb Oliver Sander <oliver.sander at tu-dresden.de>:
> 
> Dear pdelab team,
> 
> I have recently tried to explain to somebody why the lambda_*** methods
> are separate in the LocalOperator interface, rather than having alpha_***
> assemble the entire residual.  As it turns out, I am not sure I know
> the reason myself.
> 
> Initially, I thought this was to make Newton solvers more efficient.
> A Newton solver needs to call alpha_*** at each iteration, but lambda_***
> only once.  However, the Newton implementation in dune-pdelab calls
> lambda_*** at each iteration.  In fact, there is no way to access the
> alpha_*** contribution throught the GridOperator interface at all!
> So that's not the reason.

I think this was actually the original intent (in fact, it was probably driven by the
desire to avoid loading the solution in the PDELab machinery, as the lambda_***
methods don’t have access to the solution. The assembler machine contains lots
of flags and hooks to do exactly that. But as you say, there is no useful way of
accessing this functionality (the only way would be to write a local operator that
sets all of the doAlpha*** to false and the doLambda*** to true, and off the top of
my head, I can’t think of a good reason to do so).

> 
> Then I learned that the finite difference code in NumericalJacobian uses
> only alpha_***, but disregards lambda_***.  That makes sense: in first-order
> FD formulas, the lambda_*** parts cancel out anyway.

Yes, that’s the other reason.

Finally, there’s one more thing (that almost nobody uses): If you move the constant
part to lambda_***, you can implement jacobian_apply_*** by just forwarding to
alpha_*** for linear operators. But I’ll be the first to say that that’s not a good reason
for the interface ;-)

> 
> But now I am wondering: Is that really all?  Because that reason does
> not really convince me:
> 
> - While having lambda_*** separate may speed up FD approximations of the
>  Jacobian, it seems likely that it will slow down other things.
>  For example, calling lambda_*** and alpha_*** together means that
>  shape functions have to be computed twice, there are two separate
>  quadrature loops, etc.

Yes, for all high-perfomance code that I write, all lambda_*** methods are off.
Performance-wise it’s horrific.

> 
> - Having separate lambda_*** methods makes LocalOperator implementations
>  much larger (statically).  Adding a normal source term to a Poisson assembler
>  is simply two extra lines in alpha_volume.  A separate lambda_volume
>  method is at least 20 lines.

Absolutely.

> 
> - Its doubtful whether FD approximations are desirable to begin with:
>  Algorithmic differentiation is faster and more precise.

I think we can all agree on that, but it is still a good thing to have the FD tools for
quick prototyping etc. For systems, writing the jacobian by hand can be a lot of work
 - which you will end up doing for stability/performance reasons, but only after you
know that your basic operator makes sense.

That said, the current state of FD is awful and Christian and I want to change it (as
talked about in Heidelberg this spring).

Apart from that, I also consider the current LocalOperator interface really ugly, but that
won’t change for the "old" PDELab (without dune-functions).

TL;DR: I would probably skip the lambda_*** in a book.

Steffen

> 
> Can I please get some insight from the pdelab team?  Am I missing something?
> 
> Thanks,
> Oliver
> 
> _______________________________________________
> dune-pdelab mailing list
> dune-pdelab at dune-project.org
> http://lists.dune-project.org/mailman/listinfo/dune-pdelab

-------------- next part --------------
A non-text attachment was scrubbed...
Name: signature.asc
Type: application/pgp-signature
Size: 488 bytes
Desc: Message signed with OpenPGP
URL: <https://lists.dune-project.org/pipermail/dune-pdelab/attachments/20170516/ee833525/attachment.sig>