[petsc-dev] TS RHS Jacobian caching (shift/scaling) in PETSc

Stefano Zampini stefano.zampini at gmail.com
Thu Apr 25 01:58:41 CDT 2019


It appears that the testsuite does not currently cover a nasty code we have
to cache RHS Jacobian computation when using the implicit interface. If I
run the testsuite with the below modification, no tests fail. (my
continuous adjoint codes instead show some reusage)

This caching mechanism is causing headaches to me and Hong when supporting
It would be nice if we can come up with a proper test (try to cover as much
user cases as possible) and provide a definitive fix for this (even
erroring if we don't support something, instead of silently doing the wrong

I'm pretty sure that when Jacobian computations are performed by TS itself,
using the proper matrices A,B obtained from TSGetIJacobian(ts,&A,B,...) the
caching is correct ( I have spent hours on this)

When different matrices are used, TSComputeIJacobian() without an IJacobian
callback  is buggy: note that  TSComputeIJacobian() is a public method.

diff --git a/src/ts/interface/ts.c b/src/ts/interface/ts.c
index 4b003c8f6e..02211b62ea 100644
--- a/src/ts/interface/ts.c
+++ b/src/ts/interface/ts.c
@@ -575,6 +575,7 @@ PetscErrorCode  TSComputeRHSJacobian(TS ts,PetscReal
t,Vec U,Mat A,Mat B)
   if (ts->rhsjacobian.time == t && (ts->problem_type == TS_LINEAR ||
(ts->rhsjacobian.Xid == Uid && ts->rhsjacobian.Xstate == Ustate)) &&
(rhsfunction != TSComputeRHSFunctionLinear)) {
     /* restore back RHS Jacobian matrices if they have been shifted/scaled
     if (A == ts->Arhs) {
+      PetscPrintf(PetscObjectComm((PetscObject)ts),"REUSE %g
       if (ts->rhsjacobian.shift != 0) {
         ierr = MatShift(A,-ts->rhsjacobian.shift);CHKERRQ(ierr);

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20190425/89abafec/attachment.html>

More information about the petsc-dev mailing list