I have a simple system of which i would like to build the phase portrait.
I thought of creating a mesh of interesting points and solve the problem for each of them.
This is what i naively tried until now but it does not work.
xs = -0.5:0.1:1.5
ys = -0.5:0.1:3
u0s = [[x,y] for x in xs, y in ys]
tspan = (0.0, 200.0)
λ = 0
prob = ODEProblem(f, u0s, tspan, λ)
sol = solve(prob, progress=true)
I could not find anything specific in the documentation on how to approach this problem, i know about remake but i am not sure if it is useful or a good idea to use it with so many new initialisations.
Could a for loop work where every time i store the solution for every initial point?
something on the likes of
for u0 in u0s
sol = solve(remake(prob; u0=u0))
end
and then store sol somewhere for later access and/or plotting.
(in this case, how do i store the solutions correctly?)
You are right about both points, i fixed the code and now at least runs.
But now if i try to do
plot!(sol)
with the new computed solution it does not get updated.
Anyhow, i was mainly wondering if the ODEProblem interface accepted some kind of “multiple initial point” kind of thing, the current approach feels rather hacky.
The following works – I’ve used a boring pendulum model, since you don’t give your model:
function pendulum(du,u,p,t)
g = 9.81
L = 2
alpha = u[1]
omega = u[2]
du[1] = omega
du[2] = -g/L*sin(alpha)
end
#
as = -2pi:pi/6:2pi
ws = -0.5:0.1:0.5
us = [[a,w] for a in as, w in ws]
tspan = (0,5.)
#
prob = ODEProblem(pendulum,us[1],tspan)
sol = solve(prob)
plot(sol,vars=(1,2),label="")
#
for u in us
sol = solve(remake(prob,u0 = u))
plot!(sol,vars=(1,2),label="")
end
#
plot!()
