Consensus on the future of integrate.ode
Dear list, I think it would be interesting to start a discussion regarding the direction of integrate.ode. Reading the list, there seems to be many people interested seeing improvements. Is there a roadmap for development? If not can we come to a consensus about how best to contribute to this part of scipy? I think it is fair to say that other projects have overtaken scipy in terms integrating new solvers and general development in this area. Which is a shame. Personally, speaking as a user, I would like to see the use of more modern solvers. For example the code could be updated to use CVODE. This has many improvements/fixes over the Fortran VODE solver that is currently used. This would be a nice place to start, clearly that is still a nontrivial task. Also, is there a way that developers could be encourage to contribute additional solvers? For example by making an elegant low level API that makes this easier. Moreover, say that in order for a developer to connect a new solver to scipy all that had to be done was implement a particular, well defined, interface for their wrapper. Once implemented the new solver could be called with the exciting integrate.ode interface. Has anybody got experience is such design? Surely something like this has been attempted? The GNU Scientific Library interface springs to mind as a good example, http://www.gnu.org/software/gsl/manual/html_node/Ordinary-Differential-Equat... Finally, I have read some discussion here regarding changing the API and or unifying the odeint and ode. I don't want to get too distracted with this as my main point is outlined above. Personally the current API doesn't bother me too much, I find both easy to use. But if this is being considered for an update I think the crucial thing is to have the interface have no reliance on any particular solver. Abstraction is key. Happy to know your thoughts and if we can make a push together for some improvements in this area then all the better. Best regard, Dan
On 09/04/2013 06:14 PM, boyfarrell@gmail.com wrote:
Dear list,
I think it would be interesting to start a discussion regarding the direction of integrate.ode. Reading the list, there seems to be many people interested seeing improvements. Is there a roadmap for development? If not can we come to a consensus about how best to contribute to this part of scipy? I think it is fair to say that other projects have overtaken scipy in terms integrating new solvers and general development in this area. Which is a shame.
Personally, speaking as a user, I would like to see the use of more modern solvers. For example the code could be updated to use CVODE. This has many improvements/fixes over the Fortran VODE solver that is currently used. This would be a nice place to start, clearly that is still a nontrivial task.
Also, is there a way that developers could be encourage to contribute additional solvers? For example by making an elegant low level API that makes this easier. Moreover, say that in order for a developer to connect a new solver to scipy all that had to be done was implement a particular, well defined, interface for their wrapper. Once implemented the new solver could be called with the exciting integrate.ode interface. Has anybody got experience is such design? Surely something like this has been attempted? The GNU Scientific Library interface springs to mind as a good example, http://www.gnu.org/software/gsl/manual/html_node/Ordinary-Differential-Equat...
Finally, I have read some discussion here regarding changing the API and or unifying the odeint and ode. I don't want to get too distracted with this as my main point is outlined above. Personally the current API doesn't bother me too much, I find both easy to use. But if this is being considered for an update I think the crucial thing is to have the interface have no reliance on any particular solver. Abstraction is key.
Happy to know your thoughts and if we can make a push together for some improvements in this area then all the better.
Regarding the last part, I created a PR last week with some changes to the odeint function, but I didn't aim to unify odeint and ode or do a big rework of the API. I introduced a couple of convenience changes: odeint objetive function has different signature than ode (one is (y, t) and the other is (t, y)) and some minor (but backwards incompatible) changes, done to solve some annoying issues that odeint has had. The thing is that even this little changes require a proper deprecation cycle to not break everyone's code too soon, so this major uplifting in ode/odeint must be taken very carefully. The good thing about odeint is that it's a simple function to GetThingsDone™ when you just want to quickly integrate some equation. You just need a lambda, an initial condition and a vector of time values. Having to configure a solver with its options is not very straightforward, specially for interactive work, as Gaël Varoquaux had already described in a previous discussion some years ago: http://osdir.com/ml/python-scientific-devel/2009-02/msg00042.html (this is true for everything else IMHO, also for building splines and stuff like that). Perhaps the most important thing is that someone with interest in this patiently iterates with this idea, searches consensus in the mailing list, write a proposal somewhere if needed, thinks about the proper deprecation process with the maintainers, writes the code and, well, takes responsibility of this until the very end. I guess it can take some months. Juan Luis Cano
I commented my ideas in the pull request here: https://github.com/scipy/scipy/issues/2818#issuecomment-23646270 Main point: I think ode class should be redone to drop all the vode specific stuff and keep to the minimum, using like matlab set_options for all the rest, and then current odeint implemented as a simple wrapper around ode. Benny 2013/9/5 Juan Luis Cano <juanlu001@gmail.com>
Dear list,
I think it would be interesting to start a discussion regarding the
On 09/04/2013 06:14 PM, boyfarrell@gmail.com wrote: direction of integrate.ode. Reading the list, there seems to be many people interested seeing improvements. Is there a roadmap for development? If not can we come to a consensus about how best to contribute to this part of scipy? I think it is fair to say that other projects have overtaken scipy in terms integrating new solvers and general development in this area. Which is a shame.
Personally, speaking as a user, I would like to see the use of more
modern solvers. For example the code could be updated to use CVODE. This has many improvements/fixes over the Fortran VODE solver that is currently used. This would be a nice place to start, clearly that is still a nontrivial task.
Also, is there a way that developers could be encourage to contribute
additional solvers? For example by making an elegant low level API that makes this easier. Moreover, say that in order for a developer to connect a new solver to scipy all that had to be done was implement a particular, well defined, interface for their wrapper. Once implemented the new solver could be called with the exciting integrate.ode interface. Has anybody got experience is such design? Surely something like this has been attempted? The GNU Scientific Library interface springs to mind as a good example, http://www.gnu.org/software/gsl/manual/html_node/Ordinary-Differential-Equat...
Finally, I have read some discussion here regarding changing the API and
or unifying the odeint and ode. I don't want to get too distracted with this as my main point is outlined above. Personally the current API doesn't bother me too much, I find both easy to use. But if this is being considered for an update I think the crucial thing is to have the interface have no reliance on any particular solver. Abstraction is key.
Happy to know your thoughts and if we can make a push together for some
improvements in this area then all the better.
Regarding the last part, I created a PR last week with some changes to the odeint function, but I didn't aim to unify odeint and ode or do a big rework of the API. I introduced a couple of convenience changes: odeint objetive function has different signature than ode (one is (y, t) and the other is (t, y)) and some minor (but backwards incompatible) changes, done to solve some annoying issues that odeint has had. The thing is that even this little changes require a proper deprecation cycle to not break everyone's code too soon, so this major uplifting in ode/odeint must be taken very carefully.
The good thing about odeint is that it's a simple function to GetThingsDone™ when you just want to quickly integrate some equation. You just need a lambda, an initial condition and a vector of time values. Having to configure a solver with its options is not very straightforward, specially for interactive work, as Gaël Varoquaux had already described in a previous discussion some years ago:
http://osdir.com/ml/python-scientific-devel/2009-02/msg00042.html
(this is true for everything else IMHO, also for building splines and stuff like that).
Perhaps the most important thing is that someone with interest in this patiently iterates with this idea, searches consensus in the mailing list, write a proposal somewhere if needed, thinks about the proper deprecation process with the maintainers, writes the code and, well, takes responsibility of this until the very end. I guess it can take some months.
Juan Luis Cano _______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org http://mail.scipy.org/mailman/listinfo/scipy-dev
Could odeint be useful in this context? http://headmyshoulder.github.io/odeint-v2/ "Odeint is a modern C++ library for numerically solving Ordinary Differential Equations. It is developed in a generic way using Template Metaprogramming which leads to extraordinary high flexibility at top performance." Best, Arnd
Please, please don't implement anything in the SciPy ODE solvers with Boost. odeint unless it's a separate interface that requires extreme generality. Boost.odeint uses Boost.uBLAS for its linear algebra (see http://headmyshoulder.github.io/odeint-v2/doc/boost_numeric_odeint/concepts/...), and Boost.uBLAS is slow because it tries to be general (everything is templated), and it has debugging features added by default. A lot of users (myself included) don't need that generality, and would prefer the extra performance that you'd get from a code that can use a faster dense linear algebra library (like an optimized LAPACK, Eigen, Elemental, maybe MTL4). Geoff On Thu, Sep 5, 2013 at 5:25 AM, Arnd Baecker <arnd.baecker@web.de> wrote:
Could odeint be useful in this context?
http://headmyshoulder.github.io/odeint-v2/
"Odeint is a modern C++ library for numerically solving Ordinary Differential Equations. It is developed in a generic way using Template Metaprogramming which leads to extraordinary high flexibility at top performance."
Best, Arnd _______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org http://mail.scipy.org/mailman/listinfo/scipy-dev
-- Geoffrey Oxberry, Ph.D., E.I.T.
Let's say that boost.ublas is not optimized, not because it is templated, but because of something else. Eigen is very efficient, as is nt2, and there are both pure templates! Cheers, Le 5 sept. 2013 22:50, "Geoff Oxberry" <goxberry@mit.edu> a écrit :
Please, please don't implement anything in the SciPy ODE solvers with Boost .odeint unless it's a separate interface that requires extreme generality.
Boost.odeint uses Boost.uBLAS for its linear algebra (see http://headmyshoulder.github.io/odeint-v2/doc/boost_numeric_odeint/concepts/...), and Boost.uBLAS is slow because it tries to be general (everything is templated), and it has debugging features added by default.
A lot of users (myself included) don't need that generality, and would prefer the extra performance that you'd get from a code that can use a faster dense linear algebra library (like an optimized LAPACK, Eigen, Elemental, maybe MTL4).
Geoff
On Thu, Sep 5, 2013 at 5:25 AM, Arnd Baecker <arnd.baecker@web.de> wrote:
Could odeint be useful in this context?
http://headmyshoulder.github.io/odeint-v2/
"Odeint is a modern C++ library for numerically solving Ordinary Differential Equations. It is developed in a generic way using Template Metaprogramming which leads to extraordinary high flexibility at top performance."
Best, Arnd _______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org http://mail.scipy.org/mailman/listinfo/scipy-dev
-- Geoffrey Oxberry, Ph.D., E.I.T.
_______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org http://mail.scipy.org/mailman/listinfo/scipy-dev
Hello all,
Perhaps the most important thing is that someone with interest in this patiently iterates with this idea, searches consensus in the mailing list, write a proposal somewhere if needed, thinks about the proper deprecation process with the maintainers, writes the code and, well, takes responsibility of this until the very end. I guess it can take some months.
I was hoping that we have enough people that we could avoid the Rambo developer phenomenon. Besides this seems to be exactly what we have right now. What I would like to do is find is a group of people that are interested in seeing improvements and working together to fix things. Remark 1. Maybe we should set a few simple goals? For example, would people be interested in replacing vode with the cvode solver from sundials? This would be a nice improvement. What would you like to see? Remark 2. I don't understand the problem well enough to design an general new interface that solvers can be connected to. I'm not even sure is this will help. I think someone with more experience writing numerical code could architect something like this. Best wishes, Dan
On 09/07/2013 09:40 AM, boyfarrell@gmail.com wrote:
Hello all,
Perhaps the most important thing is that someone with interest in this patiently iterates with this idea, searches consensus in the mailing list, write a proposal somewhere if needed, thinks about the proper deprecation process with the maintainers, writes the code and, well, takes responsibility of this until the very end. I guess it can take some months. I was hoping that we have enough people that we could avoid the Rambo developer phenomenon. Besides this seems to be exactly what we have right now. What I would like to do is find is a group of people that are interested in seeing improvements and working together to fix things.
Remark 1. Maybe we should set a few simple goals? For example, would people be interested in replacing vode with the cvode solver from sundials? This would be a nice improvement. What would you like to see?
Remark 2. I don't understand the problem well enough to design an general new interface that solvers can be connected to. I'm not even sure is this will help. I think someone with more experience writing numerical code could architect something like this.
(I am sorry: long email ahead) I have not been in this list for long, but as far as I can tell, the problem is: there are two *unrelated* interfaces to integrate differential equations. The simple one is `odeint`, which uses a manually C-wrapped version of lsoda.f and does the iteration on the C side, and the complex one (also in the mathematical sense) is the `ode` class, which uses f2py generated interfaces to vode and other packages and does the iteration in the Python side. These two have several incompatibilities, the most silly one being the different signature of the objective function: `odeint` expects `func(y, t)` and `ode` expects `func(t, y)`. So let's say we want to solve the problem of the interfaces being unrelated. Then you would make `odeint` a simplified wrapper to the `ode` class, providing some nice defaults to rapidly integrate an ode in a one-liner. Then we step on the `ode` territory: the routines it wraps are outdated and other packages, and they probably don't see bug fixes anymore. So it is desirable to replace these with newly developed codes, like the mentioned sundials. Benny Malengier already created the odes scikit, which if I understand correctly wraps sundials. Plus, related to this replacement, it's been suggested also that we should make the `ode` class more solver-agnostic, providing an more unified interface for options and a MATLAB-style way of setting extra parameters, because as Benny already pointed out in the previously linked comment, "the scipy ode class is also too much linked to the old vode solver". In that comment he also points out which he considers are the main pieces for a full DE solution. I think this is a decent summary of what we are discussing (others please correct my mistakes and omissions). Now my opinion is: from an API point of view, I think all these are good things, both making odeint a wrapper to ode and simplifying and reworking the parameters handling in the ode class. From an implementation point of view, my only concern is: now `odeint` iterates in a compiled language (C in current SciPy, Fortran 90 in my pull request), whereas `ode` iterates on Python. There are C/Fortran <-> Python callbacks anyway, but while working on my `odeint` rewrite I noticed that the Fortran iteration was ~30 % faster. If after rewriting `ode` we keep iterating on Python and make `odeint` a wrapper around it, we might see performance losses. Other than that, I am not really familiarized with the state-of-the-art differential equations solving and such low-level numerical code (just reading lsoda.f was an enormous pain). The reason I called for a Rambo developer here is that, at the end of the day, both `odeint` and `ode` work decently, people use them and probably most really don't care about this "problems". That's why we need enough motivated people to carry the change. I ain't be that Rambo, but definitely will help. I hope this clarified something, apart from cluttering everybody's inbox :) Regards Juan Luis Cano
Hello Juan Luis, Yes, I saw your pull request and this, in part, motivated me to start this discussion. Looks like you put a lot of work into your odeint modifications! I completely agree with you, that ode should be fixed first, then odeint should be a nice interface wrapped around it. integrate.ode becomes a package aimed at people who understand the complexities of solving ODE problems and integrate.odeint is for people who just want an answer. Also, thank you for doing a nice summary, that was really helpful. Remark 1: Backend. I don't think we need to reinvent the wheel, Benny's work looks excellent so I think it make sense to base a new integrate.ode off scikits.odes. It has a nice array of modern solvers, in both C and Fortran. Remark 2 Frontend. We can define a new integrate.ode interface that wraps scikits.odes. It seems that you already have some ideas in that area, mentioning the MATLAB style. I looked at the MATLAB docs just now and spend 20 minutes thinking how this might look if it were more pythonic. You can take a look here, fork your own version to make changes or add comments, https://gist.github.com/danieljfarrell/6482713 Best wishes, Dan
2013/9/8 boyfarrell@gmail.com <boyfarrell@gmail.com>
Hello Juan Luis,
Yes, I saw your pull request and this, in part, motivated me to start this discussion. Looks like you put a lot of work into your odeint modifications!
I completely agree with you, that ode should be fixed first, then odeint should be a nice interface wrapped around it. integrate.ode becomes a package aimed at people who understand the complexities of solving ODE problems and integrate.odeint is for people who just want an answer.
Also, thank you for doing a nice summary, that was really helpful.
Remark 1: Backend.
I don't think we need to reinvent the wheel, Benny's work looks excellent so I think it make sense to base a new integrate.ode off scikits.odes. It has a nice array of modern solvers, in both C and Fortran.
Thanks for the thumbs up. However, odes is only one of 3 implementations, and then the one that exposes least of the C solvers. Odes was based originally on the previous pysundials which was no longer maintained and was a non-cython wrapper. The other implementations are geared towards a specific problem domain though. One is: https://github.com/casadi/casadi/wiki (LGPL tough) The other: http://www.jmodelica.org/assimulo (annoying copyright transfer contribution license though to a company,http://www.jmodelica.org/page/14) The non sundials solvers in odes are also not as state of the art as some that are added in above two interfaces. The problem with above packages is their license and the fact that they are packages with parsing language, ..., see eg: https://github.com/casadi/casadi/blob/master/examples/python/dae_single_shoo.... Note that you can pass an array of times to sundials at which you want output, and then all iteration is done inside sundials, so the point of time step outside of sundials is not needed (though in many problems you'll want to do stepping controlled in a python loop). Benny
Remark 2 Frontend.
We can define a new integrate.ode interface that wraps scikits.odes. It seems that you already have some ideas in that area, mentioning the MATLAB style. I looked at the MATLAB docs just now and spend 20 minutes thinking how this might look if it were more pythonic. You can take a look here, fork your own version to make changes or add comments,
https://gist.github.com/danieljfarrell/6482713
Best wishes,
Dan _______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org http://mail.scipy.org/mailman/listinfo/scipy-dev
There's also petsc4py (https://code.google.com/p/petsc4py/), which is a set of Python bindings for the scientific library PETSc. PETSc implements a lot of solvers, including an interface to CVODE in SUNDIALS. However, it generally does things using sparse linear algebra, so if you want to solve a system of ODEs using a direct linear solver, you essentially have to use LU decomposition as a preconditioner, and pass the inverted linear system to CVODE (which will then solve it trivially in one iteration of GMRES, TFQMR, or BiCGStab). There's no interface for IDA, although PETSc has a number of DAE solvers implemented itself. The license is BSD, and is well-maintained. I'd say the main drawbacks of the petsc4py library are that: - the library itself is not well-documented, but PETSc is, and the library generally follows the PETSc interface, plus has a number of Python examples - PETSc's a huge library, and it takes a while to get a handle on the API PETSc and scipy are trying to do different, but related things. PETSc and petsc4py are trying to provide a software platform for developing computational science applications that can run on anything from desktops with a single processor to a Blue Gene/Q with tens or hundreds of thousands of processors. SciPy is more geared towards rapid prototyping and computational exploration on a single processor, so it helps to have a more intuitive API, preferably something like MATLAB since that's so prevalent in education, and a lot of people who use SciPy come from a MATLAB background. Even though the licenses are compatible, I'm not sure you'd want to incorporate elements of petsc4py into SciPy (although it would be incredibly cool to be able to use any of their extremely long list of ODE or DAE solvers). Assimulo does not require you to write your ODE system in a parsing language, though in the Python implementation of JModelica, you could certainly do that. I used it for my thesis without too much difficulty, and sent the authors a few feature requests. Assimulo makes assumptions about the arguments being passed to CVODE, so you can't pass arbitrary data to CVODE via the Assimulo interface without some kludging; they assume you're going to pass an array of doubles, essentially. (Or at least, this was the case 15 months ago; I have no idea if they've changed the interface.) On Mon, Sep 9, 2013 at 3:37 AM, Benny Malengier <benny.malengier@gmail.com>wrote:
2013/9/8 boyfarrell@gmail.com <boyfarrell@gmail.com>
Hello Juan Luis,
Yes, I saw your pull request and this, in part, motivated me to start this discussion. Looks like you put a lot of work into your odeint modifications!
I completely agree with you, that ode should be fixed first, then odeint should be a nice interface wrapped around it. integrate.ode becomes a package aimed at people who understand the complexities of solving ODE problems and integrate.odeint is for people who just want an answer.
Also, thank you for doing a nice summary, that was really helpful.
Remark 1: Backend.
I don't think we need to reinvent the wheel, Benny's work looks excellent so I think it make sense to base a new integrate.ode off scikits.odes. It has a nice array of modern solvers, in both C and Fortran.
Thanks for the thumbs up. However, odes is only one of 3 implementations, and then the one that exposes least of the C solvers. Odes was based originally on the previous pysundials which was no longer maintained and was a non-cython wrapper. The other implementations are geared towards a specific problem domain though.
One is: https://github.com/casadi/casadi/wiki (LGPL tough)
The other: http://www.jmodelica.org/assimulo (annoying copyright transfer contribution license though to a company,http://www.jmodelica.org/page/14)
The non sundials solvers in odes are also not as state of the art as some that are added in above two interfaces. The problem with above packages is their license and the fact that they are packages with parsing language, ..., see eg: https://github.com/casadi/casadi/blob/master/examples/python/dae_single_shoo....
Note that you can pass an array of times to sundials at which you want output, and then all iteration is done inside sundials, so the point of time step outside of sundials is not needed (though in many problems you'll want to do stepping controlled in a python loop).
Benny
Remark 2 Frontend.
We can define a new integrate.ode interface that wraps scikits.odes. It seems that you already have some ideas in that area, mentioning the MATLAB style. I looked at the MATLAB docs just now and spend 20 minutes thinking how this might look if it were more pythonic. You can take a look here, fork your own version to make changes or add comments,
https://gist.github.com/danieljfarrell/6482713
Best wishes,
Dan _______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org http://mail.scipy.org/mailman/listinfo/scipy-dev
_______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org http://mail.scipy.org/mailman/listinfo/scipy-dev
-- Geoffrey Oxberry, Ph.D., E.I.T.
On 09/09/2013 10:27 AM, Geoff Oxberry wrote:
PETSc and scipy are trying to do different, but related things. PETSc and petsc4py are trying to provide a software platform for developing computational science applications that can run on anything from desktops with a single processor to a Blue Gene/Q with tens or hundreds of thousands of processors.
SciPy is more geared towards rapid prototyping and computational exploration on a single processor, so it helps to have a more intuitive API, preferably something like MATLAB since that's so prevalent in education, and a lot of people who use SciPy come from a MATLAB background. Even though the licenses are compatible, I'm not sure you'd want to incorporate elements of petsc4py into SciPy (although it would be incredibly cool to be able to use any of their extremely long list of ODE or DAE solvers).
I think this is a good point that has been made. Maybe we could aim for a less complex package that fits our (more modest) needs.
The non sundials solvers in odes are also not as state of the art as some that are added in above two interfaces. The problem with above packages is their license and the fact that they are packages with parsing language, ..., see eg: https://github.com/casadi/casadi/blob/master/examples/python/dae_single_shoo... .
What if we just stay with your odes? It's still a big improvement, perhaps a good tradeoff between state of the art algorithms, simplicity and licensing compatibility.
Note that you can pass an array of times to sundials at which you want output, and then all iteration is done inside sundials, so the point of time step outside of sundials is not needed (though in many problems you'll want to do stepping controlled in a python loop).
Benny
Remark 2 Frontend.
We can define a new integrate.ode interface that wraps scikits.odes. It seems that you already have some ideas in that area, mentioning the MATLAB style. I looked at the MATLAB docs just now and spend 20 minutes thinking how this might look if it were more pythonic. You can take a look here, fork your own version to make changes or add comments,
https://gist.github.com/danieljfarrell/6482713
Best wishes,
Dan _______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org <mailto:SciPy-Dev@scipy.org> http://mail.scipy.org/mailman/listinfo/scipy-dev
_______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org <mailto:SciPy-Dev@scipy.org> http://mail.scipy.org/mailman/listinfo/scipy-dev
-- Geoffrey Oxberry, Ph.D., E.I.T.
_______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org http://mail.scipy.org/mailman/listinfo/scipy-dev
2013/9/9 Juan Luis Cano <juanlu001@gmail.com>
On 09/09/2013 10:27 AM, Geoff Oxberry wrote:
PETSc and scipy are trying to do different, but related things. PETSc and petsc4py are trying to provide a software platform for developing computational science applications that can run on anything from desktops with a single processor to a Blue Gene/Q with tens or hundreds of thousands of processors.
SciPy is more geared towards rapid prototyping and computational exploration on a single processor, so it helps to have a more intuitive API, preferably something like MATLAB since that's so prevalent in education, and a lot of people who use SciPy come from a MATLAB background. Even though the licenses are compatible, I'm not sure you'd want to incorporate elements of petsc4py into SciPy (although it would be incredibly cool to be able to use any of their extremely long list of ODE or DAE solvers).
I think this is a good point that has been made. Maybe we could aim for a less complex package that fits our (more modest) needs.
The non sundials solvers in odes are also not as state of the art as
some that are added in above two interfaces. The problem with above packages is their license and the fact that they are packages with parsing language, ..., see eg: https://github.com/casadi/casadi/blob/master/examples/python/dae_single_shoo....
What if we just stay with your odes? It's still a big improvement, perhaps a good tradeoff between state of the art algorithms, simplicity and licensing compatibility.
Yes, for current scipy it would be a logical evolution. For my own work though I might need more complex stuff present in sundials and not exposed yet like in PETsc (preconditioning, krylov, parallel). But then, these things can be added as needed, and as PETsc has a nice license, we are allowed to look at it without fear of contamination. Benny
2013/9/9 Benny Malengier <benny.malengier@gmail.com>
2013/9/8 boyfarrell@gmail.com <boyfarrell@gmail.com>
Hello Juan Luis,
Yes, I saw your pull request and this, in part, motivated me to start this discussion. Looks like you put a lot of work into your odeint modifications!
I completely agree with you, that ode should be fixed first, then odeint should be a nice interface wrapped around it. integrate.ode becomes a package aimed at people who understand the complexities of solving ODE problems and integrate.odeint is for people who just want an answer.
Also, thank you for doing a nice summary, that was really helpful.
Remark 1: Backend.
I don't think we need to reinvent the wheel, Benny's work looks excellent so I think it make sense to base a new integrate.ode off scikits.odes. It has a nice array of modern solvers, in both C and Fortran.
Thanks for the thumbs up. However, odes is only one of 3 implementations, and then the one that exposes least of the C solvers. Odes was based originally on the previous pysundials which was no longer maintained and was a non-cython wrapper. The other implementations are geared towards a specific problem domain though.
One is: https://github.com/casadi/casadi/wiki (LGPL tough)
The other: http://www.jmodelica.org/assimulo (annoying copyright transfer contribution license though to a company,http://www.jmodelica.org/page/14)
The non sundials solvers in odes are also not as state of the art as some that are added in above two interfaces. The problem with above packages is their license and the fact that they are packages with parsing language, ..., see eg: https://github.com/casadi/casadi/blob/master/examples/python/dae_single_shoo....
Note that you can pass an array of times to sundials at which you want output, and then all iteration is done inside sundials, so the point of time step outside of sundials is not needed (though in many problems you'll want to do stepping controlled in a python loop).
What would really be interesting is find some way to use the parallel implementation of sundals in python. Some mapping of a numpy array to the parallel sundial vector ( https://github.com/bmcage/odes/blob/master/scikits/odes/sundials/nvector/nve...) would then be needed. Has there ever been work to contruct a parallel aware array? Would be nice to find funding for that somewhere Benny
Benny
Remark 2 Frontend.
We can define a new integrate.ode interface that wraps scikits.odes. It seems that you already have some ideas in that area, mentioning the MATLAB style. I looked at the MATLAB docs just now and spend 20 minutes thinking how this might look if it were more pythonic. You can take a look here, fork your own version to make changes or add comments,
https://gist.github.com/danieljfarrell/6482713
Best wishes,
Dan _______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org http://mail.scipy.org/mailman/listinfo/scipy-dev
On Mon, Sep 9, 2013 at 5:53 AM, Benny Malengier <benny.malengier@gmail.com>wrote:
2013/9/9 Benny Malengier <benny.malengier@gmail.com>
2013/9/8 boyfarrell@gmail.com <boyfarrell@gmail.com>
Hello Juan Luis,
Yes, I saw your pull request and this, in part, motivated me to start this discussion. Looks like you put a lot of work into your odeintmodifications!
I completely agree with you, that ode should be fixed first, then odeintshould be a nice interface wrapped around it. integrate.ode becomes a package aimed at people who understand the complexities of solving ODE problems and integrate.odeint is for people who just want an answer.
Also, thank you for doing a nice summary, that was really helpful.
Remark 1: Backend.
I don't think we need to reinvent the wheel, Benny's work looks excellent so I think it make sense to base a new integrate.ode off scikits.odes. It has a nice array of modern solvers, in both C and Fortran.
Thanks for the thumbs up. However, odes is only one of 3 implementations, and then the one that exposes least of the C solvers. Odes was based originally on the previous pysundials which was no longer maintained and was a non-cython wrapper. The other implementations are geared towards a specific problem domain though.
One is: https://github.com/casadi/casadi/wiki (LGPL tough)
The other: http://www.jmodelica.org/assimulo (annoying copyright transfer contribution license though to a company, http://www.jmodelica.org/page/14)
The non sundials solvers in odes are also not as state of the art as some that are added in above two interfaces. The problem with above packages is their license and the fact that they are packages with parsing language, ..., see eg: https://github.com/casadi/casadi/blob/master/examples/python/dae_single_shoo....
Note that you can pass an array of times to sundials at which you want output, and then all iteration is done inside sundials, so the point of time step outside of sundials is not needed (though in many problems you'll want to do stepping controlled in a python loop).
What would really be interesting is find some way to use the parallel implementation of sundals in python. Some mapping of a numpy array to the parallel sundial vector ( https://github.com/bmcage/odes/blob/master/scikits/odes/sundials/nvector/nve...) would then be needed. Has there ever been work to contruct a parallel aware array?
petsc4py does this already, using data structures from PETSc. You can instantiate some of these data structures with NumPy arrays. I'm not sure if such work has already been done in SciPy, nor am I sure if SciPy is "the right place" for it. I don't know a ton about parallel computing, but my limited experience with it backs up the common warning people give about parallelizing code -- you really have to parallelize your data structures in order for most algorithms to be effective. My worry is that people will write SciPy code as they have done for the previous releases -- that is, in serial -- and then expect parallel ODE functionality to just work, when continuing to use serial data structures won't let you take advantage of parallelism. The other thing about wrapping the parallel implementation of SUNDIALS I'd be concerned about is the preconditioning of iterative linear solvers; in most cases, these solvers have to be preconditioned in order to converge quickly. For that, do you expect users to provide their own preconditioners, like SUNDIALS does? Something that is useful about solvers like PETSc/petsc4py or Trilinos/PyTrilinos is that they provide a number of methods that can also be used to precondition iterative linear solvers (variants of incomplete LU/Cholesky factorization, block Jacobi, algebraic multigrid, and others), which makes it easier for users to select from a number of "black box" methods in the event that they cannot easily construct a preconditioner themselves. These sorts of things would probably be better implemented in other parts of SciPy, since they could be used for other purposes besides solving ODEs (namely, for any application involving the solution of a linear system using iterative methods). Geoff
Would be nice to find funding for that somewhere
Benny
Benny
Remark 2 Frontend.
We can define a new integrate.ode interface that wraps scikits.odes. It seems that you already have some ideas in that area, mentioning the MATLAB style. I looked at the MATLAB docs just now and spend 20 minutes thinking how this might look if it were more pythonic. You can take a look here, fork your own version to make changes or add comments,
https://gist.github.com/danieljfarrell/6482713
Best wishes,
Dan _______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org http://mail.scipy.org/mailman/listinfo/scipy-dev
_______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org http://mail.scipy.org/mailman/listinfo/scipy-dev
-- Geoffrey Oxberry, Ph.D., E.I.T.
Hello all,
Even though the licenses are compatible, I'm not sure you'd want to incorporate elements of petsc4py into SciPy (although it would be incredibly cool to be able to use any of their extremely long list of ODE or DAE solvers).
PETSc seems powerful and I love the idea of having many solver than can take advantage of multiple cores. That is certainly very attractive. I agree with your downsides Geoff; it looks terrifying! This is compounded by the fact that I can't even tell how to solve a simple problem because the docs are poor (well not as friendly as scipy at least). Here is what I think we should do: 1) Update interface to integrate.ode. - see my suggestion here, https://gist.github.com/danieljfarrell/6482713 - comments welcomed 2) Make integrate.odeint a functional wrapper around integrate.ode. - We should use the concept of 'event' from the MATLAB/sundials API here, it would be a nice addition. 3) Update the integrate.ode solvers - replace VODE with CVODE - add CODES - add a DAE solver(s). Point 3 is the most unclear. What technology do we use as the foundation? We have mentioned (in not particular order), * petsc4py, https://code.google.com/p/petsc4py/ * scikit.odes, http://cage.ugent.be/~bm/progs.html * assimulo, http://www.jmodelica.org/assimulo * casadi, https://github.com/casadi/casadi/wiki I really don't know where we should start as I have used none of the above. Generally the higher level we start at the more quickly we can get things done. Maybe Geoff and Benny (and others) can comment. Best wishes, Dan
On Sep 9, 2013, at 9:08 AM, boyfarrell@gmail.com wrote:
PETSc seems powerful and I love the idea of having many solver than can take advantage of multiple cores. That is certainly very attractive. I agree with your downsides Geoff; it looks terrifying! This is compounded by the fact that I can't even tell how to solve a simple problem because the docs are poor (well not as friendly as scipy at least).
The petsc4py docs are effectively nonexistent, although the examples are not bad. When I was working on a FiPy-PETSc interface (not quite there, yet) I spent a lot of time in the petsc4py code to figure out what it was actually doing and what it was calling in PETSc-proper (not that I know that interface, either, but it's at least documented).
participants (7)
-
Arnd Baecker -
Benny Malengier -
boyfarrell＠gmail.com -
Geoff Oxberry -
Guyer, Jonathan E. Dr. -
Juan Luis Cano -
Matthieu Brucher