Velocity Reviews > need help with ODE solver from scipy

# need help with ODE solver from scipy

T.Crane
Guest
Posts: n/a

 04-24-2007
Hi,

OK, I'm trying to figure out how to use the ODE solver
(scipy.integrate.ode.ode). Here's what I'm doing (in iPython)

y0 = [0,1,1]
dt = 0.01
tEnd = 12
t0 = 0
Y = zeros([tEnd/dt, 3])

As an aside, I've used this assignment for Y in the past and it
worked. When I tried it this morning I got a TypeError and a message
saying I needed to use an integer. So, instead...

Y = zeros([int(tEnd/dt), 3])
T = zero([int(tEnd/dt)])
index = 0
def foo(t,y):
dy = zeros([3])
dy[0] = y[1]*y[2]
dy[1] = -y[0]*y[2]
dy[2] = -0.51*y[0]*y[1]
return dy

r = ode(foo).set_integrator('vode').set_initial_value( y0,t0)
while r.successful() and r.t < tEnd:
r.integrate(r.t+dt)
Y[index] = r.y
T[index] = r.t
index += 1

As a further aside, this system of three coupled linear ODE's is from
an example in MATLAB's documentation on their function ode45, so I
know what I 'should' get. The while loop and call to ode is adapted
from scipy's (very limited) documentation on scipy.integrate.ode.ode.

At the end of this while loop I've gotten a few different results
depending on what I have no idea. This morning another TypeError
exception was thrown, saying:

TypeError: Array can not be safely cast to required type.

Although, to be fair this is after the output from one iteration:

array([0.01, 1. , 1. ])

So, clearly this isn't working right. Does anyone have any experience
using this function or anything they can contribute?

thanks,
trevis

Robert Kern
Guest
Posts: n/a

 04-24-2007
T.Crane wrote:
> Hi,
>
> OK, I'm trying to figure out how to use the ODE solver
> (scipy.integrate.ode.ode).

You will get more help from the scipy-user list than you will here.

http://www.scipy.org/Mailing_Lists

> Here's what I'm doing (in iPython)
>
> y0 = [0,1,1]
> dt = 0.01
> tEnd = 12
> t0 = 0
> Y = zeros([tEnd/dt, 3])
>
> As an aside, I've used this assignment for Y in the past and it
> worked. When I tried it this morning I got a TypeError and a message
> saying I needed to use an integer.

Well, the current behavior is correct. When you give us a float instead of an
int for a dimension, we shouldn't guess what you want. Particularly, in this
case, truncating is incorrect. I ran your code and eventually wound up with an
IndexError because your iteration ran past the end of Y and T.

>
> Y = zeros([int(tEnd/dt), 3])
> T = zero([int(tEnd/dt)])

Also, a better way to do this is simply to leave these as lists and simply
append to them. That way, you don't have to manage indices.

> index = 0
> def foo(t,y):
> dy = zeros([3])
> dy[0] = y[1]*y[2]
> dy[1] = -y[0]*y[2]
> dy[2] = -0.51*y[0]*y[1]
> return dy
>
> r = ode(foo).set_integrator('vode').set_initial_value( y0,t0)
> while r.successful() and r.t < tEnd:
> r.integrate(r.t+dt)
> Y[index] = r.y
> T[index] = r.t
> index += 1
>
> As a further aside, this system of three coupled linear ODE's is from
> an example in MATLAB's documentation on their function ode45, so I
> know what I 'should' get. The while loop and call to ode is adapted
> from scipy's (very limited) documentation on scipy.integrate.ode.ode.
>
> At the end of this while loop I've gotten a few different results
> depending on what I have no idea. This morning another TypeError
> exception was thrown, saying:
>
> TypeError: Array can not be safely cast to required type.
>
> Although, to be fair this is after the output from one iteration:
>
> array([0.01, 1. , 1. ])
>
> So, clearly this isn't working right. Does anyone have any experience
> using this function or anything they can contribute?

I tried your example and it worked for me (after fixing the IndexError). Could
you come to scipy-user with the complete code that you ran (the one above isn't
the one you ran; there is at least one typo), and the full traceback of the error?

--
Robert Kern

"I have come to believe that the whole world is an enigma, a harmless enigma
that is made terrible by our own mad attempt to interpret it as though it had
an underlying truth."
-- Umberto Eco