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))
and then store sol somewhere for later access and/or plotting.
(in this case, how do i store the solutions correctly?)
Any smarter ideas, this feels really dumb.
EDIT: added the fix in the second answer
You are right about both points, i fixed the code and now at least runs.
But now if i try to do
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:
g = 9.81
L = 2
alpha = u
omega = u
du = omega
du = -g/L*sin(alpha)
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,tspan)
sol = solve(prob)
for u in us
sol = solve(remake(prob,u0 = u))
The result is as follows (not very interesting…):
Thanks a bunch!
I had tried before without the final
plot!() and it did not work since i forgot the
display inside the for loop!
https://docs.juliadiffeq.org/latest/features/ensemble/ is for that.
I had trouble finding the right terminology in the documentations, thanks for pointing me to it!