# Solving same DiffEq repeatedly but with different initial conditions

I have a set of coupled ODEs that I’d like to solve, which is simple enough using DifferentialEquations.jl. My problem is that I need to solve the same ODE repeatedly for many different initial conditions parametrized by an array `x`.

Currently I just do this in a loop. I’m sure there’s a smarter way but I haven’t quite found it in the documentation yet. This is very inefficient of course since I need to keep setting up the problem again every time.

Could someone give me some pointers?

Here is a minimal non-working example of what I’m doing:

``````using Elliptic.Jacobi

function dn_ab!(du,u,p,t)
du[1] .= +im*p[1].*u[1] .+ im.*u[2].*dn.(t,real(p[2]))
du[2] .= -im*p[1].*u[2] .+ im.*u[1].*dn.(t,real(p[2]))
end

function ab_dn_t0(x)
A = +χ .+ 0.5.*Ω*λ.*(x - xs) .- π/4
B = -χ .+ 0.5.*Ω*λ.*(x - xs).- π/4
return [2*im*sin(A); 2*cos(B)]
end

# How do you get rid of this whole loop?
tspan = (0.0,abs(minimum(c.box.t)))
for i in eachindex(x)
u0 = ab_dn_t0(x[i])
prob = ODEProblem(dn_ab!, u0, tspan, [λ, m])
sol = solve(prob, Tsit5(), saveat=t) # t is an array
# do some processing on sol and save the result
end
# do some more processing on the saved result for all x

``````

You can use the integrator interface with `reinit`.

Or you can use the ensemble interface to auto-parallelize it.

This includes automated GPU parallelism with DiffEqGPU.jl

4 Likes

Thanks Chris, that resolves it for my use case.

That’s wonderful, I’ll definitely move to this at some point if my system gets large enough.