[aspect-devel] Newton solver

Wolfgang Bangerth bangerth at tamu.edu
Fri Sep 25 15:09:05 PDT 2015

On 09/25/2015 11:47 AM, Menno Fraters wrote:
> Thanks for your detailed response. I really like that you structured it in a
> way in which it can be implemented in small pieces.

Thanks for the kind words -- I don't have all that much time to program myself 
any more, but I think I'm still pretty good with this sort of conceptual work 
of designing code and laying out the steps to get there :-)

> My first question is if you could elaborate a bit on why the first option has
> your preference?

The problem with matrix free methods is always the issue with preconditioners. 
If you don't have access to a matrix, you don't have access to matrix entries, 
and therefore don't have any easy way to define preconditioners that are based 
on matrix entries. This is true for simple preconditioners like Jacobi or SSOR 
as well as more complicated ones like algebraic multigrid (AMG). The only 
option you have is to use geometric multigrid methods because there you can 
assemble the level matrices as you need them on each level, just like you 
would for the system matrix.

A second problem is that if you needed to assemble (cell-local) matrix 
contributions many times over because you want to work matrix-free (where you 
need to do assembly on every cell for every GMRES iteration and every inner CG 
iteration), then you really need to be fast with assembly. deal.II has some 
facilities that could be used for this, but the current scheme would almost 
certainly be too slow for it. It would also make it a far bigger problem than 
it currently is that some of the material models need to do a *lot* of work at 
every quadrature point to compute coefficient values. This is not so much of a 
problem in the current scheme where we do the assembly only once per timestep, 
but it would become a major obstacle if you had to do it in every iteration.

> My second question is if it would be possible the make evaluating the
> derivative optional, because otherwise if I would implement this, I would need
> to implement the derivative for all current material models of which some work
> with quite complex rheologies and cut of the viscosity at a upper and lower
> threshold, which can make calculating the derivative very hard.

Correct: you'd make it optional at first and only allow using the Newton 
solver for models that implement computing derivatives.

We could come up with ways to compute derivatives automatically using 
automatic differentiation (see step-33 for an example). There clearly are also 
cases where exact derivatives don't exist and in that case one would have to 
run Newton methods with inexact derivatives -- quite a common thing to do, in 


Wolfgang Bangerth               email:            bangerth at math.tamu.edu
                                 www: http://www.math.tamu.edu/~bangerth/

More information about the Aspect-devel mailing list