# how to do a Backward integration using ode_solver?

Let's say that i want to integrate this from ti=0 to tf=10

T=ode_solver()
f = lambda t,y:[y[1],y[0]^2-y[1]]
T.function=f
T.y_0=[0,.1]
T.t_span=[0,10]
T.ode_solve(num_points=100)
H= [ [T.solution[i][0] , T.solution[i][1][0] ,T.solution[i][1][1] ] for i in range(len(T.solution))  ]
for i in range(5): H[i]

[0, 0, 0.100000000000000]
[0.10000000000000001, 0.0095163351611087642, 0.090486759458914337]
[0.20000000000000001, 0.018128063040304006, 0.081894953217975602]
[0.30000000000000004, 0.025923510483407335, 0.0741488103498774]
[0.40000000000000002, 0.032983605743555654, 0.06717622415610694]


How i do to integrate backwards, from ti=10 to tf=0? This doesn't works, because the solution that it gives me is a constant = ti:

T=ode_solver()
f = lambda t,y:[y[1],y[0]^2-y[1]]
T.function=f
T.y_0=[10,.05]
T.t_span=[10,0]
T.ode_solve(num_points=100)
H= [ [T.solution[i][0] , T.solution[i][1][0] ,T.solution[i][1][1] ] for i in range(len(T.solution))  ]
for i in range(5): H[i]

[10, 10, 0.0500000000000000]
[10.0, 10.0, 0.050000000000000003]
[9.9000000000000004, 10.0, 0.050000000000000003]
[9.8000000000000007, 10.0, 0.050000000000000003]
[9.7000000000000011, 10.0, 0.050000000000000003]


edit retag close merge delete

Sort by » oldest newest most voted

Unfortunately, ode_solve presently does not allow for backwards integration. (For reference, in the file "devel/sage/sage/gsl/ode.pyx", the line

while (t < t_end):

controls the loop that evolves the system) The documentation of GSL suggests that the latest versions of the library should be able to integrate backwards, but naive changes to the python interface routine didn't solve the problem for me.

work-around: you could rewrite your system in s=-t and integrate from s=-10 to 0 instead.

An alternative is scipy.integrate.odeint. I haven't found a thorough comparison of the GSL and SciPy (ODEPACK) ODE solvers, though, and I haven't checked if ODEPACK allows for backwards integration.

more

Many thanks! finally i decided to use the scipy odeint, since supports the backwards integration and is the fastest ode solver that i found

fa= lambda y,t: [y[1],y[0]^2-y[1]]

t=scipy.linspace(10, 0, 1/0.01)
len(t)
time l=integrate.odeint(fa,[0,.1],-t)

more