Extend odeint to handle matrices and complex numbers
I wrote a wrapper of odeint called odeintw that is currently available on github: https://github.com/WarrenWeckesser/odeintw To quote from README.md in the repository: odeintw provides a wrapper of scipy.integrate.odeint that allows it to handle complex and matrix differential equations. That is, it can solve equations of the form dZ/dt = F(Z, t, param1, param2, ...) where t is real and Z is a real or complex array. What do folks think of adding this functionality directly to odeint in scipy? It would be completely backwards compatible, and would not introduce any overhead in the callbacks if the inputs are 1-d real arrays. Warren
Warren Weckesser writes:
I wrote a wrapper of odeint called odeintw that is currently available on github: https://github.com/WarrenWeckesser/odeintw
To quote from README.md in the repository:
odeintw provides a wrapper of scipy.integrate.odeint that allows it to handle complex and matrix differential equations. That is, it can solve equations of the form
dZ/dt = F(Z, t, param1, param2, ...)
where t is real and Z is a real or complex array.
What do folks think of adding this functionality directly to odeint in scipy? It would be completely backwards compatible, and would not introduce any overhead in the callbacks if the inputs are 1-d real arrays.
Warren
I cannot really comment on the code internals, but to have transparent nd-array-valued vector fields would be a major improvement (so don't stop at "matrix"!). I have encountered the need over and over again, and one can always flatten out the structures, but a lot of code clarity, in particular the closeness between mathematical representation and the code, is getting lost. And why stop at odeint? There is integrate.ode which currently has the same limitation. --Marcel
On Mon, Nov 16, 2015 at 4:46 PM, Marcel Oliver < m.oliver@jacobs-university.de> wrote:
Warren Weckesser writes:
I wrote a wrapper of odeint called odeintw that is currently available on github: https://github.com/WarrenWeckesser/odeintw
To quote from README.md in the repository:
odeintw provides a wrapper of scipy.integrate.odeint that allows it to handle complex and matrix differential equations. That is, it can solve equations of the form
dZ/dt = F(Z, t, param1, param2, ...)
where t is real and Z is a real or complex array.
What do folks think of adding this functionality directly to odeint in scipy? It would be completely backwards compatible, and would not introduce any overhead in the callbacks if the inputs are 1-d real arrays.
Warren
I cannot really comment on the code internals, but to have transparent nd-array-valued vector fields would be a major improvement (so don't stop at "matrix"!).
I shouldn't have said "matrix", because in fact it already handles n-dimensional arrays. The case of a 2-d matrix is just the most common request that I've seen. For example, in the following a system in the shape of a 2x2x2 array is solved: In [21]: from odeintw import odeintw In [22]: def func(x, t): ....: return -x ....: In [23]: x0 = np.arange(8.0).reshape(2,2,2) In [24]: t = np.arange(5.0) In [25]: sol = odeintw(func, x0, t) In [26]: sol.shape Out[26]: (5, 2, 2, 2)
I have encountered the need over and over again, and one can always flatten out the structures, but a lot of code clarity, in particular the closeness between mathematical representation and the code, is getting lost.
And why stop at odeint? There is integrate.ode which currently has the same limitation.
That's certainly possible. I referred to odeintw as a wrapper of odeint (and it is), but most of the code is actually about wrapping the user's callback functions (i.e. the differential equations and the Jacobian function, if given). Warren
--Marcel _______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org https://mail.scipy.org/mailman/listinfo/scipy-dev
Warren, It would be most useful to add that to odes, https://github.com/bmcage/odes which interfaces the modern (Krylov precond, ...) cvode code instead of the old codebases included in scipy. We have been working on a new api for odes which returns named tuples, see the examples here: https://github.com/bmcage/odes/tree/master/docs/ipython (chrome might be needed). This as a way to move forward the ode solvers in scipy which don't have a nice/good/liked API (many suggestions for rewrites have been suggested in the past). The solution as named tuple can be accessed as sol.values.t and roots eg sol.roots.t, ... One drawback to deprecate the old solvers in scipy has been the complex solver zvode. Allen Hindmarsh ( http://history.siam.org/oralhistories/hindmarsh.htm) did not include a complex solver in the new codes in sundials. Several people have indicated their need for (transparent) complex ODE solvers. Your approach might solve this. Personally, I advise against using old solvers as lsode or vode as available in scipy for actual academic (or industrial) use. Benny 2015-11-16 23:12 GMT+01:00 Warren Weckesser <warren.weckesser@gmail.com>:
On Mon, Nov 16, 2015 at 4:46 PM, Marcel Oliver < m.oliver@jacobs-university.de> wrote:
Warren Weckesser writes:
I wrote a wrapper of odeint called odeintw that is currently available on github: https://github.com/WarrenWeckesser/odeintw
To quote from README.md in the repository:
odeintw provides a wrapper of scipy.integrate.odeint that allows it to handle complex and matrix differential equations. That is, it can solve equations of the form
dZ/dt = F(Z, t, param1, param2, ...)
where t is real and Z is a real or complex array.
What do folks think of adding this functionality directly to odeint in scipy? It would be completely backwards compatible, and would not introduce any overhead in the callbacks if the inputs are 1-d real arrays.
Warren
I cannot really comment on the code internals, but to have transparent nd-array-valued vector fields would be a major improvement (so don't stop at "matrix"!).
I shouldn't have said "matrix", because in fact it already handles n-dimensional arrays. The case of a 2-d matrix is just the most common request that I've seen.
For example, in the following a system in the shape of a 2x2x2 array is solved:
In [21]: from odeintw import odeintw
In [22]: def func(x, t): ....: return -x ....:
In [23]: x0 = np.arange(8.0).reshape(2,2,2)
In [24]: t = np.arange(5.0)
In [25]: sol = odeintw(func, x0, t)
In [26]: sol.shape Out[26]: (5, 2, 2, 2)
I have encountered the need over and over again, and one can always flatten out the structures, but a lot of code clarity, in particular the closeness between mathematical representation and the code, is getting lost.
And why stop at odeint? There is integrate.ode which currently has the same limitation.
That's certainly possible. I referred to odeintw as a wrapper of odeint (and it is), but most of the code is actually about wrapping the user's callback functions (i.e. the differential equations and the Jacobian function, if given).
Warren
--Marcel _______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org https://mail.scipy.org/mailman/listinfo/scipy-dev
_______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org https://mail.scipy.org/mailman/listinfo/scipy-dev
On Mon, Nov 16, 2015 at 6:00 PM, Benny Malengier <benny.malengier@gmail.com> wrote:
Warren,
It would be most useful to add that to odes, https://github.com/bmcage/odes which interfaces the modern (Krylov precond, ...) cvode code instead of the old codebases included in scipy. We have been working on a new api for odes which returns named tuples, see the examples here: https://github.com/bmcage/odes/tree/master/docs/ipython (chrome might be needed). This as a way to move forward the ode solvers in scipy which don't have a nice/good/liked API (many suggestions for rewrites have been suggested in the past). The solution as named tuple can be accessed as sol.values.t and roots eg sol.roots.t, ...
One drawback to deprecate the old solvers in scipy has been the complex solver zvode. Allen Hindmarsh ( http://history.siam.org/oralhistories/hindmarsh.htm) did not include a complex solver in the new codes in sundials. Several people have indicated their need for (transparent) complex ODE solvers. Your approach might solve this.
Personally, I advise against using old solvers as lsode or vode as available in scipy for actual academic (or industrial) use.
Benny
Hi Benny, I'm familiar with `odes`, but I haven't looked at it in a while. I just cloned the repo, and I'll start poking around. I don't think it will be difficult to do a wrapper for `odes` similar to what `odeintw` does for scipy's `odeint`. Warren
2015-11-16 23:12 GMT+01:00 Warren Weckesser <warren.weckesser@gmail.com>:
On Mon, Nov 16, 2015 at 4:46 PM, Marcel Oliver < m.oliver@jacobs-university.de> wrote:
Warren Weckesser writes:
I wrote a wrapper of odeint called odeintw that is currently available on github: https://github.com/WarrenWeckesser/odeintw
To quote from README.md in the repository:
odeintw provides a wrapper of scipy.integrate.odeint that allows it to handle complex and matrix differential equations. That is, it can solve equations of the form
dZ/dt = F(Z, t, param1, param2, ...)
where t is real and Z is a real or complex array.
What do folks think of adding this functionality directly to odeint in scipy? It would be completely backwards compatible, and would not introduce any overhead in the callbacks if the inputs are 1-d real arrays.
Warren
I cannot really comment on the code internals, but to have transparent nd-array-valued vector fields would be a major improvement (so don't stop at "matrix"!).
I shouldn't have said "matrix", because in fact it already handles n-dimensional arrays. The case of a 2-d matrix is just the most common request that I've seen.
For example, in the following a system in the shape of a 2x2x2 array is solved:
In [21]: from odeintw import odeintw
In [22]: def func(x, t): ....: return -x ....:
In [23]: x0 = np.arange(8.0).reshape(2,2,2)
In [24]: t = np.arange(5.0)
In [25]: sol = odeintw(func, x0, t)
In [26]: sol.shape Out[26]: (5, 2, 2, 2)
I have encountered the need over and over again, and one can always flatten out the structures, but a lot of code clarity, in particular the closeness between mathematical representation and the code, is getting lost.
And why stop at odeint? There is integrate.ode which currently has the same limitation.
That's certainly possible. I referred to odeintw as a wrapper of odeint (and it is), but most of the code is actually about wrapping the user's callback functions (i.e. the differential equations and the Jacobian function, if given).
Warren
--Marcel _______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org https://mail.scipy.org/mailman/listinfo/scipy-dev
_______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org https://mail.scipy.org/mailman/listinfo/scipy-dev
_______________________________________________ SciPy-Dev mailing list SciPy-Dev@scipy.org https://mail.scipy.org/mailman/listinfo/scipy-dev
2015-11-17 0:19 GMT+01:00 Warren Weckesser <warren.weckesser@gmail.com>:
On Mon, Nov 16, 2015 at 6:00 PM, Benny Malengier < benny.malengier@gmail.com> wrote:
Warren,
It would be most useful to add that to odes, https://github.com/bmcage/odes which interfaces the modern (Krylov precond, ...) cvode code instead of the old codebases included in scipy. We have been working on a new api for odes which returns named tuples, see the examples here: https://github.com/bmcage/odes/tree/master/docs/ipython (chrome might be needed). This as a way to move forward the ode solvers in scipy which don't have a nice/good/liked API (many suggestions for rewrites have been suggested in the past). The solution as named tuple can be accessed as sol.values.t and roots eg sol.roots.t, ...
One drawback to deprecate the old solvers in scipy has been the complex solver zvode. Allen Hindmarsh ( http://history.siam.org/oralhistories/hindmarsh.htm) did not include a complex solver in the new codes in sundials. Several people have indicated their need for (transparent) complex ODE solvers. Your approach might solve this.
Personally, I advise against using old solvers as lsode or vode as available in scipy for actual academic (or industrial) use.
Benny
Hi Benny,
I'm familiar with `odes`, but I haven't looked at it in a while. I just cloned the repo, and I'll start poking around. I don't think it will be difficult to do a wrapper for `odes` similar to what `odeintw` does for scipy's `odeint`.
Warren
Great. Instead of a wrapper an option that activates complex/matrix is also a possible way. Possibly somewhat more transparent in actual use (eg set option to true if any of the given parameters are complex, instead of an if structure to select another integrator obejct (eg odew instead of ode). There are 3 pull requests being worked on, https://github.com/bmcage/odes/pulls, which when ready should finalize the new API. After that a release is planned. All inputs on improvements to the new API are welcome. Idea is to later use cvodes instead of cvode as soon as options indicate user wants also sensitivity output, after which also sensitivities would be present in the output named tuple. Handling complex/matrix with an option would fit with this approach also. Benny
participants (3)
-
Benny Malengier -
Marcel Oliver -
Warren Weckesser