Hi,
I'm experimenting with the use of PETSc to replace a DAE solver in an
existing code that I use to solve stiff nonlinear problems. I expect to
use TSBDF in the final instance, and so am currently playing with it but
applied to a simpler linear problem - just to get some experience with
the SNES/KSP/PC controls before diving in to the hard problem.
Below is some output from TSAdapt for the simple linear problem, using
TSBDF and PCLU, so that the linear algebra solve in the newton loop is
direct:
TSAdapt basic bdf 0:2 step 0 accepted t=0 + 1.000e-03
dt=2.000e-03 wlte=2.51e-07 wltea= -1 wlter= -1
TSResidual...
TSJacobian... calculate
TSResidual...
TSAdapt basic bdf 0:2 step 1 accepted t=0.001 + 2.000e-03
dt=4.000e-03 wlte=2.83e-07 wltea= -1 wlter= -1
TSResidual...
TSJacobian... calculate
TSResidual...
TSAdapt basic bdf 0:2 step 2 accepted t=0.003 + 4.000e-03
dt=8.000e-03 wlte=1.22e-07 wltea= -1 wlter= -1
TSResidual...
TSJacobian... calculate
TSResidual...
I have added the "TSResidual..." and "TSJacobian..." echoes so that I
can see when PETSc is requesting residuals and Jacobians to be
computed. (This is the Jacobian routine specified via TSSetIJacobian.)
Regarding the above output, it appears that TS / SNES always requests a
new (I)Jacobian at each new timestep (after the first residual is
calculated). I can see mathematically why this would be the default
choice, but had hoped that it might be possible for out-of-date
Jacobians to be used until they become inefficient. My reason for
wanting this is that the Jacobian calculations for the intended
application are particularly expensive, but for small enough timesteps
out-of-date Jacobians may be good enough, for a few steps.
Is there any way of specifying that out-of-date (I)Jacobians can be
tolerated (at the expense of increased Newton iterations, or smaller
timesteps)? Alternatively would it make sense to include callbacks to
TS / SNES from the Jacobian evaluation function to determine whether
sufficiently few iterations have been used that it might be safe to
return the previously calculated Jacobian (if I store a copy)? If so,
is there any advice on how I should do this?
NB. I see that there is an option for TSRHSJacobianSetReuse(), but this
only applies to the RHS component of the DAE (the G(t,u) part, using the
terminology from the manual), but I am not using this as ultimately I
expect to be solving strongly nonlinear problems with no "slow" G(t,u) part.
Any advice would be greatly appreciated.
