[SciPy-Dev] GSoC' 15 Ideas: define general integrators; parallelize integration; higher-order ODEs

[Benny Malengier]

2015-03-17 11:12 GMT+01:00 Max Mertens <max.mail at dameweb.de>:

> @Geordie McBain, Benny Malengier: Thank you for your feedback, and for
> the online/paper references.
>
> @Benny:
> > This for stiff problems and implicit methods! The python code would
> > not care as long as normal arrays can be passed.
> Do you mean it is needed, or difficult instead to implement for
> stiff/implicit ODEs?
>

I mean it would be great if parallel backend is an option for the cases
where a sparse parallel linear solver is possible.

 > Most specific implementations have their own datastructures however.
> >
> > Code like fipy casts PDE problems to equations which can be solved via
> > pytrilinos, http://trilinos.org/packages/pytrilinos/ in parallel via
> > MPI. You get out numpy arrays, internally something else is used.
> >
> > The cvode solver of sundials for stiff problems has a fully parallel
> > implementation, and sundials has some RK methods in ARKode. I don't
> > think any of the python bindings expose that at the moment. I suppose
> > the mapping to parallel array is not straightforward. Having a look at
> > the parallel implementation there might also give ideas though for a
> > general interface.
> > Pytrilinos interfaces parts of sundials via Rythmos
> > (
> http://trilinos.org/docs/dev/packages/rythmos/doc/html/classRythmos_1_1ImplicitBDFStepper.html
> ),
> > but it's unclear to me if that is parallel or not. In any case, the
> > data structures used are not ndarray, but must be mapped
> > to Epetra.Vector. All far too high level abstraction to be usable in
> > scipy. Scipy should focus on simple interfaces for ode/dae, but the
> > existing examples of parallel implementation seem to indicate simple
> > is no longer possible then.
> >
> > Whatever the approach, scipy should not redo work present in high
> > level packages like pytrilinos, but instead offer the basis to start
> > from, so people can evolve to those packages if needed.
> If I understand you correctly, you suggest to not implement parallel
> solvers, as those exist in other libraries already, but rather provide a
> general interface to those?
>

No, what I mean is that if a parallel solver is offered, it must remain a
simple API, and not a very detailed API geared towards a specific area
(like molecular modelling). Scipy will not replace dedicated packages. The
methods present should be great in general. I'm not a core scipy developer
though, they might have other ideas on what is fitting.

Also, for integrate, scipy mostly is a wrapper layer around academic, well
tested, codes. I'm not sure own written solvers, or pure python solvers,
would be accepted. So you would need to select an existing well tested
parallel solver, and then wrap that.

Personally I like the approach of maple with their student package:
http://www.maplesoft.com/support/help/Maple/view.aspx?path=Student%2fNumericalAnalysis
The basic integrators are present there. It paints a clear picture: look,
you learn this stuff in numerical analysis, so here are the methods, but
for really doing ode or dae, use dsolve numberic
http://www.maplesoft.com/support/help/maple/view.aspx?path=dsolve%2fnumeric
instead (though that has the classical option too with the common classical
methods (
http://www.maplesoft.com/support/help/maple/view.aspx?path=dsolve%2fclassical
).

Some sort of 'student/classical' version for integrate in scipy is
something I always missed.


>
> What about the other ideas, a general interface to define RK (or other)
> integrators, and methods to automate higher-order to first-order ODE
> transformation?
>

I did not react to RK as I have no experience there. As I said, if there is
a good paper that allows an implementation, then yes, it probably is an
option.
But mostly, a wrapper to an existing well tested codebase for RK seems the
fastest and best approach for a student.
Current dopri5 and dop853 in scipy seem a very specific choice for RK, but
the fact they have step size control makes them in practise better than
whatever standard RK you could devise based on textbooks. If the aim is
just to have an implementation for classical methods, eg RK variants,
available in scipy, then adding those via a clear student/classical package
is something I would like, but perhaps not core scipy developers.

Benny


> Do you suggest something similar as a project for me?
>
> Regards,
> Max
>
> _______________________________________________
> SciPy-Dev mailing list
> SciPy-Dev at scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-dev
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/scipy-dev/attachments/20150317/5ebd19c9/attachment.html>