ODE.jl
ODE.jl copied to clipboard
ode45 hangs on backwards integration
Hello,
I came across this issue (recreated here using example code):
julia> using ODE
julia> function f(t, y)
# Extract the components of the y vector
(x, v) = y
# Our system of differential equations
x_prime = v
v_prime = -x
# Return the derivatives as a vector
[x_prime; v_prime]
end;
julia> # Initial condtions -- x(0) and x'(0)
const start = [0.0; 0.1]
WARNING: redefining constant start
2-element Array{Float64,1}:
0.0
0.1
julia> time = 0.:-.1:-4pi;
julia> t, y = ode45(f, start, time);
This hangs for a long time (I've waited up to a minute before killing it).
I belive I have the same problem.
ode45
seems to go into an infitine loop under some conditions.
For example:
function f(t,v,p)
(x,y) = v
dx = p[1]*x - x*y
dy = p[2]*x*y - y
return [dx;dy]
end
p = [4.82,-2.43]
init_vals = [0.5,0.5]
times = 1.0:1.0:10.0
g = (t,v) -> func(t,v,params)
(t,d) = ode45(g,init_vals,time_points,points=:specified)
I also waited longer than a minute before killing it.
Maybe this has to do with the 1.01 here: https://github.com/JuliaLang/ODE.jl/blob/a0ef97d6249e62a089370787850ae697e3d94962/src/runge_kutta.jl#L341 . Which I think is the wrong thing to do.
Could you check if setting it to 1 works?
Unfortunately, it does not. I changed it and still had to kill the process after a minute. Doesn't look like the method is even going into that if clause
This equation is stiff. Try a stiff solver (e.g. ode23s
or DASSL.jl
) instead, both can solve your problem very quickly.
EDIT: as a plan for a distant future we should add a stiffness detection system (sometimes it's difficult to see if the equation is stiff) and maybe switch between stiff and non stiff solvers.
Also, equation by @galbarel is unrelated to the original issue (there is no backward time integration).
For a temporary remedy to the original issue you might want to take a look at #49 where backward time integration works (at least in the case @crbinz described).
Thanks @pwl for noticing the stiff problem. The other problem seems to be something else, maybe related to backwards integration but not that as, note that e.g. ode78
or ODE.ode23
or ODE.ode45_fe
work fine.
Also x-ref to the original julia-users thread https://groups.google.com/d/msg/julia-users/zB2mhbX3KeU/Tqn_hxWACQAJ .