Brian,
Could send an example of your "rhs" function; not a totally trivial example
Barry
> On Dec 7, 2015, at 11:21 AM, Brian Merchant <[email protected]> wrote:
>
> Hi all,
>
> I am considering using petsc4py instead of scipy.integrate.odeint (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.
>
> So far, I have been able to produce reasonable speed gains by writing the RHS
> functions in "something like C" (using either numba or Cython). I'd like to
> get even more performance out, hence my consideration of PETSc.
>
> 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...
>
> 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?
>
> 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:
> http://scicomp.stackexchange.com/questions/21501/is-it-worth-switching-to-timesteppers-provided-by-petsc-if-i-cant-write-down-a
>
> All your advice would be most helpful :)
>
> Kind regards,Brian
>