1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
|
using DifferentialEquations
using Plots
# Simulate Lorenz 63 system and investigate sensitivity to initial conditions
function tendency!(du, u, p, t)
x,y,z = u
σ,ρ,β = p
du[1] = dx = σ*(y-x)
du[2] = dy = x*(ρ-z) - y
du[3] = dz = x*y - β*z
end
p = [10.0, 28.0, 8/3] # parameters of the Lorentz 63 system
tspan = (0.0, 1000.0)
u0 = [1.0, 0.0, 0.0]
δu0 = [0.0, 0.001, 0.0] # small deviation in initial condition
prob = ODEProblem(tendency!, u0, tspan, p)
sol1 = solve(prob, Tsit5(), reltol=1e-8, abstol=1e-8) # control simulation
prob = ODEProblem(tendency!, u0 + δu0, tspan, p)
sol2 = solve(prob, Tsit5(), reltol=1e-8, abstol=1e-8) # perturbed initial condition
plot(sol1, idxs = (1, 2, 3)) # 3d plot
# plot(sol1, idxs = (1, 2))
plot([sol1[1, :], sol2[1, :]], [sol1[2, :], sol2[2, :]])
|
