Solve and plot Lorenz equations for two different initial conditions and two values of rho in julia

My goal is to solve lorenz equations and plot them as it shows in the figure. I have two different initial conditions [x0, 1, 0] and x0= 0 then x0 =1* 10^-5 the two values of rho are ρ= 14 and ρ=28. I wrote this code in julia but I’m getting errors says “Got exception outside of a @test UndefVarError: u not defined lorenz_solver(x₀::Float64, ρ::Int64)” any idea on how to solve this issue.

image

using DifferentialEquations
using Plots
using LaTeXStrings

function lorenz!(t,p,u,du)
    x, y ,z = u
    sigma, ρ, beta = p
    du[1] = sigma*(u[2]-u[1])
    du[2] = u[1]*(ρ-u[3]) - u[2]
    du[3] = u[1]*u[2] - beta*u[3]
end

function lorenz_solver(x₀, ρ)
    u[0] = [x₀, 1, 0]
    tspan = (0.0, 100.0)
    p = [10, ρ, 8/3] 
    x₀, ρ = params
    prob = ODEProblem(lorenz!, u0, tspan, params)
    sol = solve(prob) 
    return sol
    #return the ODESolution
end

function lorenz_plot()
    params = [0, 14]
    sol1 = solve(prob)

    params = [1 * 10^-5 , 14]
    sol2 = solve(prob)

    p1 = plot(t, [sol1 sol2], ylabel=L"x1, x2", title="ρ = 14")
    p2 = plot(abs(sol1 - sol2), xlabel=L"t", ylabel=L"|x1 - x2|")
    p3 = plot(sol1, vars = (1,2,3), xlabel=L"x1", ylabel=L"z1")

    params = [0, 28]
    sol3= solve(prob)
   
    params = [1 * 10^-5, 28]
    sol4 = solve(prob)

    p4 = plot(t, [sol3 sol4], title="ρ = 28")
    p5 = plot(abs(sol3 - sol4), xlabel=L"t")
    p6 = plot(sol3, vars = (1,2,3), xlabel=L"x1")
    link=:all
    plot(p1, p2, p3, p4, p5, p6, layout(3, 2))
    #return the final plot
end

export lorenz_solver, lorenz_plot

end

Here you might have wanted to type u0 = [x₀, 1, 0].

I just changed it and got an error says " Got exception outside of a @test
BoundsError: attempt to access 2-element Vector{Float64} at index [3]"

You have 2 different things as parameters. 1 for the general problem you want to solve (x₀, ρ) and 1 for the ODE problem of Lorenz equation. For the ODE problem, you should put p as an argument not params.

@tomaklutfu
I’m not sure if I’m getting since I’m new with Julia.
I added p as an argument here:

function lorenz!(t,p,u,du)

then I defined sigma and beta to solve the problem while keeping ρ undefined until the plots section since it’s going to take 2 different values

The problem is in the construction of the ODEProblem, especially in the above quote. Confusing yourself with the parameter of the general question in your hand, you pass a 2-valued (x₀ and ρ) variable named params as an argument to it. In fact, the parameter that should be passed is 3-valued p ([10, ρ, 8/3]) instead.

You should read this section of the tutorials carefully. The order of the arguments that you are using when defining the lorenz! function is not what the DifferentialEquations package expects.

If you check out the tutorial, you will see

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

That is the order of arguments you shoud use. Now, if you want to start simulating the system for different parameters, you can pass them via the p vector, but always defining the function as lorenz!(du,u,p,t) and not defining the arguments in a different order.

In fact, a little below on the tutorial you actually have what you want, the lorentz system parametrized.

2 Likes

I want to add to the last answer that it is necessary to keep track of what type of variables you are passing. Everything should be, for example, Float32. This can greatly affect performance