scipy/scipy — pull request #2818 by Juanlu001 (2013-09-01 20:38:18+00:00)
Title: Odeint rewrite

The `odeint` loop has been rewritten in standards-compliant Fortran 95 using
the interface generated by f2py. The code is simpler and briefier now, more
robust since further checks were added and some issues regarding the old
odeint have been addressed (see below). No performance regressions were
detected - I even saw insignificant, but existent, speedups.

There are two relevant backwards-incompatible changes. The first one: `func`
has now the signature `func(t0, y)` instead of `func(y, t0)` and same goes
for `dfunc`. Others have requested this change in the past, as the `ode` class
and also other software such as MATLAB follow this convention.

The second one: now `odeint` returns a `Result` object with the solution and
information about the status of the solver, so the new API works like this:

>>> res = odeint(func, y0, t)
>>> res.y # Solution
[...]
>>> res.success # boolean
True
>>> res.infodict # Same as old infodict
[...]

The parameter `full_output` has been therefore eliminated. This was inspired
by the `Result` object of the root finding functions in `scipy.optimize` and
has been implemented using a `namedtuple`. This can be extended to other parts
of SciPy.

Several optional parameters to `odeint` changed its name, inspired by the `ode`
class: `dfunc` (was Dfun), `first_step` (was h0), `max_step` (was hmax), ...

I didn't bother to implement the extra `args` to the objective function,
inspired by the convention in MATLAB to use partials or lambdas. I think it
simplifies the code on the SciPy side and adds minimal overhead on the user
side, plus giving more flexibility. This is a change worth extending to other
parts of SciPy too:

>>> def func(t, y, a):
[...]
>>> res = odeint(lambda t, y: func(t, y, 1.0), y0, t) # a = 1.0
>>> from functools import partial
>>> func1 = partial(func, a=1)
>>> res = odeint(func1, y0, t)

Some other parameters have disappeared: `printmessg` and `col_deriv` exist no
more. Talking about the latter, f2py creates copies if the input arrays are
C-contiguous. If performance or memory usage are an issue, using `order='F'`
inside the objective function and the jacobian prevents the creation of copies.

These issues were fixed with this rewrite:

* print exception instead of raising (#2570) - fixed, regression test added.

These issues were addressed but I'm not sure if this is the best solution
and no regression tests were added:

* exception on failure (#1801) - warnings are shown every time, to supress
the print messages modifications to lsoda.f are necessary.

* silently discards imaginary parts (#1567) - a ComplexWarning is raised
(hence not silently) but should `complex_ode` be explicitly recommended in
these cases?

These issues were not addressed:

* ignores Ctrl-C (#2515) - not sure if the problem appears only on Windows.

Finally, these changes were reflected in the docstring and the tutorial and
a note was added pointing out that the Jacobian is only used when the
problem is stiff.