<div dir="ltr"><p>Hi all,<br></p><p>I am considering using <code>petsc4py</code> instead of <code>scipy.integrate.odeint</code>
 (which is a wrapper for Fortran solvers) for a problem involving the 
solution of a system of ODEs. The problem has the potential to be stiff.
 Writing down its Jacobian is very hard.</p>

<p>So far, I have been able to produce reasonable speed gains by writing the RHS functions in "something like C" (using either <code>numba</code> or <code>Cython</code>). I'd like to get even more performance out, hence my consideration of PETSc. <br></p><p>Due to the large number of equations involved, it is already tedious to 
think about writing down a Jacobian. Even worse though, is that some of the functions governing a particular 
interaction do not have neat analytical forms (let alone whether or
 not their derivatives have neat analytical forms), so we might have a 
mess of piecewise functions needed to approximate them if we were to go 
about still trying to produce a Jacobian...</p><p>All the toy examples I see of PETSc time stepping problems have 
Jacobians defined, so I wonder if I would even get a speed gain going 
from switching to it, if perhaps one of the reasons why I have a high 
computational cost is due to not being able to provide a Jacobian 
function?</p><p>I described the sort of problem I am working with in more detail in this scicomp.stackexchange question, which is where most of this question is copied from, except it also comes with a toy version of the problem I am dealing with: <a href="http://scicomp.stackexchange.com/questions/21501/is-it-worth-switching-to-timesteppers-provided-by-petsc-if-i-cant-write-down-a">http://scicomp.stackexchange.com/questions/21501/is-it-worth-switching-to-timesteppers-provided-by-petsc-if-i-cant-write-down-a</a></p><p>All your advice would be most helpful :)<br></p><p>Kind regards,Brian<br></p></div>