DifferentialEquations - vector or one by one on rhs

Specific Domains Modelling & Simulations
1

I wonder what the difference is between

function lorenz!(du, u, p, t)
    du = [10.0 * (u[2] - u[1]), u[1] * (28.0 - u[3]) - u[2], u[1] * u[2] - (8 / 3) * u[3]]
end

and

function lorenz!(du, u, p, t)
    du[1] = 10.0 * (u[2] - u[1])
    du[2] = u[1] * (28.0 - u[3]) - u[2]
    du[3] = u[1] * u[2] - (8 / 3) * u[3]
end

The first one doesn’t work when I do:

using DifferentialEquations
u0 = [1.0; 0.0; 0.0]
tspan = (0.0, 100.0)
prob = ODEProblem(lorenz!, u0, tspan)
sol = solve(prob)

and just assigns u0:

image
image837×352 4.01 KB

I have a massive rhs and whereas I know I can assign elements to du in a loop, I am just curious why the first approach fails. I thought that was because du is a matrix, but it isn’t.

2

the first one is creating a new vector, not mutating a vector. In the in-place form, the return does not matter, it’s about updating the memory. So using .= is fine:

function lorenz!(du, u, p, t)
    du .= [10.0 * (u[2] - u[1]), u[1] * (28.0 - u[3]) - u[2], u[1] * u[2] - (8 / 3) * u[3]]
end

but is slow of course.

2 Likes
3

See also the Julia manual on assignment vs. mutation.