[SciPy-User] odeint python and mathematica
Nasser Mohieddin Abukhdeir
nasser@udel....
Wed Feb 10 12:10:53 CST 2010
It seems as though there are two different integration methods being
used, so you really need to get to the bottom of that. There is one
other important difference between SciPy and Mathematica in this case,
Mathematica has just-in-time compilation functionality:
http://www.wolfram.com/technology/guide/TransparentAutoCompilation/
while SciPy does not, although it is using compiled FORTRAN modules.
Your case is also an extreme one, a system of two ODEs, where I
typically deal with systems of 10^3+. The Python overhead is much less
of a problem for me since all of the matrix operations are done in
C-compiled code. Look into SciPy's Weave functionality, I have not had
much success with it and am working on using Cython to speed-up my code
and so I can easily interface with external C-code. It might seem like
alot of work now, but switching from Matlab to Python (and not C) has
benefited me greatly. I think it is well worth the effort, access to
VODE through scipy.integrate.ode, excellent 3D plotting using Mayavi,
and wrappers exist for the SUNDIALS ODE package if you ever need a
matrix-free implicit ODE solver (see PySundials).
Nasser
On Wed, 2010-02-10 at 12:17 -0500, Mark P wrote:
> In the Python routine:
> - number of time step evaluations 50,000
> - number of function evaluations 120,000
>
>
> In the Mathematica routine:
> - time steps 231,770
> - function evaluations 500,000
>
>
> I'm not sure how to work out the relative and absolute error
> tolerances, but I did set the 'atol' in Python to a number just before
> the answer became inaccurate and it did speed it up 2.5x but still 4x
> slower than Mathematica.
>
>
> Thank you,
>
>
> Mark
>
>
>
>
>
> On 2010-02-10, at 10:58 AM, Nasser Mohieddin Abukhdeir wrote:
>
>
>
> > Hello Mark:
> > I am quite new to Python and Scipy, but I have some experience
> > with SciPy's ODE functionality and it has been comparable to
> > compiled code, depending on my implementation of the integrand.
> > There are a few important parameters for the integrator that you
> > should take into account in order to assure a fair comparison. Based
> > upon odeint's default values, you should confirm that the relative
> > and absolute error tolerances are the same for both computations.
> > Something else that is important is the integration scheme, it looks
> > like both methods use LSODA, so it should be an even comparison. Can
> > you rerun both simulations and check the number of time steps and
> > function evaluations from each implementation?
> >
> >
> > Nasser Mohieddin Abukhdeir
> > Postdoctoral Fellow (Vlachos Research Group)
> > University of Delaware - Department of Chemical Engineering
> > (Personal Web Page) http://udel.edu/~nasser/
> > (Group Web Page) http://www.dion.che.udel.edu/
> >
> > On Wed, 2010-02-10 at 09:08 -0600, scipy-user-request@scipy.org
> > wrote:
> >
> >
> > > Message: 5
> > > Date: Wed, 10 Feb 2010 10:08:29 -0500
> > > From: Mark P <m.perrin@me.com>
> > > Subject: [SciPy-User] odeint (python vs mathematica)
> > > To: Scipy <scipy-user@scipy.org>
> > > Message-ID: <18B95518-6CA6-4C5F-BB57-3BBECC863C0B@me.com>
> > > Content-Type: text/plain; charset="us-ascii"
> > >
> > > Hello all,
> > >
> > > (my first post to the list)
> > >
> > > I have been performing modelling of ionic currents in Mathematica for the past three years as a student. I have recently passed into the hallowed halls of post-grad, and am keen to be done with Mathematica and use a free and open solution for my projects, hence Python and numpy/scipy.
> > >
> > > I've made some progress with the transition; the cookbooks are especially helpful. There is one in particular I have been looking at - the rabbit/fox model ( Lotka-Volterra Tutorial), that is helping me to understand odeint.
> > >
> > > =====
> > > from numpy import *
> > > import pylab as p
> > > # Definition of parameters
> > > a = 1.
> > > b = 0.1
> > > c = 1.5
> > > d = 0.75
> > > def dX_dt(X, t=0):
> > > """ Return the growth rate of fox and rabbit populations. """
> > > return array([ a*X[0] - b*X[0]*X[1] ,
> > > -c*X[1] + d*b*X[0]*X[1] ])
> > >
> > > from scipy import integrate
> > >
> > > t = linspace(0, 10000, 10000) # time
> > > X0 = array([10, 5]) # initials conditions: 10 rabbits and 5 foxes
> > > X = integrate.odeint(dX_dt, X0, t)
> > > ===
> > >
> > > A simple form of this code in Mathematica is:
> > >
> > > ===
> > > a = 1
> > > b = 0.1
> > > c = 1.5
> > > d = 0.75
> > >
> > > eqn1 := r'[t] == a * r[t] - b * r[t] * f[t]
> > > eqn2 := f'[t] == -c * f[t] + d * b * r[t] * f[t]
> > >
> > > solution = NDSolve[{eqn1, eqn2, r[0] == 10, f[0] == 5}, {r[t], f[t]}, {t, 0, 10000}, MaxSteps -> Infinity]
> > > ===
> > >
> > > My question is about speed. When computed in Mathematica 7 the code is roughly about 10x faster than odeint (I've increased t from 15 in the cookbook to 10000 to allow the measurement of speed to be made). But I know from the web that numpy/scipy can be of an equivalent or faster speed. I was wondering whether the solution can be sped up; and would that require the use of C or Fortran.
> > >
> > > I see the inherent advantage of using Python with numpy/scipy and I do not expect the speeds of a compiled language but I was hoping for equivalent speeds to other interpreted solutions.
> > >
> > > Thank you,
> > >
> > > Mark Perrin
> > > Electrophysiology Research Fellow
> > > Ottawa
> > > -------------- next part --------------
> > > An HTML attachment was scrubbed...
> > > URL: http://mail.scipy.org/pipermail/scipy-user/attachments/20100210/0e2124a7/attachment.html
> > >
> > > ------------------------------
> > >
> > > _______________________________________________
> > > SciPy-User mailing list
> > > SciPy-User@scipy.org
> > > http://mail.scipy.org/mailman/listinfo/scipy-user
> > >
> > >
> > > End of SciPy-User Digest, Vol 78, Issue 19
> > > ******************************************
> >
> >
> >
> >
> >
> >
>
>
Nasser Mohieddin Abukhdeir
Postdoctoral Fellow (Vlachos Research Group)
University of Delaware - Department of Chemical Engineering
(Personal Web Page) http://udel.edu/~nasser/
(Group Web Page) http://www.dion.che.udel.edu/
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mail.scipy.org/pipermail/scipy-user/attachments/20100210/b9f429d1/attachment.html
More information about the SciPy-User
mailing list