Hi all, I am integrating a stiff system of ODEs/PDEs using PETSc TS (typically with BDF or other implicit time-stepping schemes), and I would like to exploit the fact that I can efficiently compute the action of the Jacobian on a vector (Jv) without assembling the full Jacobian matrix. Since for a large system it becomes expensive to assemble the Jacobian in each iteration. In scikits.odes (SUNDIALS/CVODE), there is a native API for passing only a J*v routine to the time integrator. In my experience, when I use only a Jacobian-vector product routine (without assembling the full matrix), the performance improves significantly for large systems. However, in PETSc TS, the workflow seems more matrix-centric, and I have only found the possibility to use MatShell for the Jacobian
Is there a way to do something similar in PETSc TS (for BDF or other implicit schemes)? Currently, I use the matrix-free Newton-Krylov method to approximate the Jacobian and have adjusted the tolerances to achieve convergence, as recommended by Barry. In that case, I obtain similar integration times with scikits.odes CVODE without using the Jacobian times vector. Best regards, Art