Can you make a PR with a test to cover this feature? We desperately need coverage as a CI integration.
Stefano Zampini <stefano.zamp...@gmail.com> writes: > Jed, > > 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 > adjoints. > 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 > thing) > > 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 > %g\n",ts->rhsjacobian.shift,ts->rhsjacobian.scale); > if (ts->rhsjacobian.shift != 0) { > ierr = MatShift(A,-ts->rhsjacobian.shift);CHKERRQ(ierr); > } > > -- > Stefano