[SciPy-user] The other solvers in ODEPACK

[Anne Archibald]

2008/6/28 James <tomo.bbe at gmail.com>:

> Recently I have been using ODEPACK (direct from Netlib) in a pure Fortran
> project and have noticed dramatically different execution times between two
> of the solvers (LSODA and LSODES).
>
>
> I have also noticed the source for LSODES is in scipy's repository - so is
> it possible to access the other solvers through scipy at the moment ?

No, it's not, as Robert said. I have written a (partial) wrapper for
one of them (LSODAR) before deciding it wasn't very useful. I do have
a few comments:

* The FORTRAN interface is a colossal pain to use in a python program,
so you'd want to come up with some sort of python wrapper. There are
two very different examples in scipy already, odeint() and the ode
class. the latter of which was intended to support multiple backends
(but only one was implemented). Neither interface did what I wanted,
though, which was part of my motivation for trying to wrap things
myself.

* If your concern is straight runtime, you should be aware that for
ODE solvers implemented in an ODEPACK-like way, where the right-hand
side is an arbitrary python function, there is a very substantial
overhead associated with calling back into python for each RHS
evaluation. I realize that different algorithms may conceivably make a
bigger difference than this, but do be aware that solveing ODEs in
python in this way is quite a slow process. For this reason I would be
much more interesting in wrapping an ODE solver that was more powerful
- that is, allowed solution of more problems - than one that was more
efficient.

* If you want a more full-featured ODE solver, you should look at
PyDSTool. As its name suggests, it is oriented towards the study of
dynamical systems, but it provides many useful tools for simply
solving ODEs. High on the list is a system for providing symbolic RHSs
which can be compiled to C and (I think) automatically differentiated.
It also includes some more modern ODE solvers (though I seem to recall
none of them support automatic stiff/non-stiff switching). There was
some discussion of integrating some of this into scipy.

Anne