[SciPy-user] PyDSTool and SciPy (integrators)

[Neil Martinsen-Burrell]

Rob Clewley <rob.clewley <at> gmail.com> writes:

> >  I think the right return format depends what people want from the
> >  solution. I see two major kinds of problem one might use an ODE solver
> >  for:
> >
> >  * Integrate out to some stopping point and give the values of
> >  dependent and independent variable there.
> >  * Integrate out to some stopping point and give back a representation
> >  of the solution.
> >
> >  (The first case can in principle be obtained by evaluating the second
> >  at the endpoint, but at the cost of replacing O(1) space by O(n)
> >  space.)
> 
> I personally do not find the former very useful, as I almost always
> want to know something about the intermediate points. My understanding
> is the dopri and radau integrators inherently return the whole
> solution, so maybe an additional calling option supporting the former
> usage could just throw everything out but the last point before the
> whole solution gets copied into a numpy object.

I use both of these cases in different situations.  The first is useful in
contexts where the analysis is in terms of maps (e.g. data assimilation), but
the underling dynamics is continuous in time.  I like the ability to wrap up an
integration into a function that returns the solution after a specified time
interval.  I tend to use the second in cases where the dynamics of the solution
*is* the problem under investigation.  I think that a DS package should support
both of these modes of use: integration over a time interval as a simple
function and generation of solution objects.  It is certainly possible to
implement one in terms of the other and probably worth a convenience function to
that effect.

> >  For the second it is clear that some kind of object to represent the
> >  solution is a good idea. What exactly it should contain is open for
> >  discussion; I envision keeping estimated function values and
> >  derivatives at each step, along with convergence information, and
> >  providing an interface to view it as either an array of values or a
> >  function akin to our interpolating spline objects. Convergence
> >  information and auxiliary function values may also be useful. Even for
> >  the first, simple case, I think returning an object is a more
> >  manageable way to provide convergence information.
> 
> >  Are there other styles of using an ODE solver that people want?
> >  Perhaps a "running" mode, where the solver keeps track of whatever
> >  state information it has and is asked to advance by a certain amount
> >  (say, because more experimental data has become available)?
> 
> There are two issues that need separating here:
> 
> 1) the short term refactoring, reuse and adaptation of existing
> PyDSTool code for scipy
> 2) the long term design goals for supporting ODEs in scipy
> 
> Much as I want to keep (2) in mind, I think getting some short-term
> progress on improving scipy's ODE support based on existing code would
> be good. I foresee it taking a long time to get ODE support truly
> "right" for as many users as possible, and in the meantime I think
> some our code and other short-term solutions would be an appropriate
> stop gap.

Re (1): It looks to me as if the Trajectory object would be useful to include in
SciPy, particularly as the current functionality can be implemented on top of
it.  We should continue to provide the current interface to ODE solvers (all of
them?), while also providing a more capable object-oriented interface that we
will encourage people to use for new code.  I also think that adding a variety
of integrators is useful, particularly if those integrators are accessible
through a uniform interface (similar to what OpenOpt has done for optimizers). 
It seems that Trajectories, Coordinates and Models provide the necessary
abstractions to unify all of the existing solvers.  In addition, these objects
should be able to abstract whether or not these are constructed through
compiling generated C code, or by defining python obects.

-Neil