Description
Background: We are trying to convert a number of places in the ASPECT project (https://github.com/geodynamics/aspect/) to use KINSOL instead of hand-rolled nonlinear solvers. ASPECT is a code for large-scale simulation of processes in the Earth based on deal.II.
In a test case for a really simple case extracted by my colleague @bobmyhill (reduced from geodynamics/aspect#5698), we are trying to solve for a root of a one-argument function f(x) for which we get the following sequence of callbacks when using KIN_NONE (=Newton):
Computing residual at x=-69.77069997, f(x)=-176.46700446
Recomputing J at x=-69.77069997
Solving for dx with rhs=176.46700446
Computing residual at x=-69.77070101, f(x)=-176.4670055
Computing residual at x=106.69630449, f(x)=203.06842599
Recomputing J at x=106.69630449
Solving for dx with rhs=-203.06842599
Computing residual at x=106.69630608, f(x)=203.06843156
Computing residual at x=48.676754204, f(x)=1.4210854715e-14
I've got a number of questions about this sort of log:
- I am confused by the fact that it calls the residual function twice in a row for nearly identical points (x=-69.77069997 and x=-69.77070101; and then again for x=106.69630449 and x=106.69630608). This almost looks like it's doing a finite difference approximation, but I know that's not it because the calls happen before/after solving for dx. What is the purpose for the call after solving the Newton system at a slightly perturbed point?
- We get exactly the same sequence of events and points whether we choose
KIN_NONE(=Newton) orKIN_PICARD. This is surprising.
If it is of any help, I've annotated every one of the residual calls above with the last KINSOL functions that invoked the call to the user-supplied residual:
Computing residual at x=-69.77069997, f(x)=-176.46700446
from KINSolInit ()
from KINSol ()
Recomputing J at x=-69.77069997
Solving for dx with rhs=176.46700446
Computing residual at x=-69.77070101, f(x)=-176.4670055
from kinLsDQJtimes ()
from kinLsATimes ()
from kinLsSolve ()
from KINPicardFcnEval ()
from KINPicardAA ()
from KINSol ()
Computing residual at x=106.69630449, f(x)=203.06842599
from KINPicardAA ()
from KINSol ()
Recomputing J at x=106.69630449
Solving for dx with rhs=-203.06842599
Computing residual at x=106.69630608, f(x)=203.06843156
from kinLsDQJtimes ()
from kinLsATimes ()
from kinLsSolve ()
from KINPicardFcnEval ()
from KINPicardAA ()
from KINSol ()
Computing residual at x=48.676754204, f(x)=1.4210854715e-14
from KINPicardAA ()
from KINSol ()
I assume (but don't know) that the functions named "Picard" are also doing the Newton method?