need help with ODE solver from scipy
T.Crane
trevis.crane at gmail.com
Tue Apr 24 11:53:15 EDT 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
More information about the Python-list
mailing list