Система уравнений Лоренца возникла в результате изучения физики атмосферы.

Он назван в честь физика Эдварда Лоренца, который в своем стремлении смоделировать тенденции в климатологии и метеорологии наткнулся на область хаотических динамических систем.



По сути, изучение хаотических систем направлено на выявление основных тенденций и корреляций в необъяснимых иначе явлениях.

Область обширна и заслуживает того, чтобы посвятить всю жизнь исследованиям. И есть те, кто в академических кругах делают именно это.

Однако, что, возможно, более ценно, так это использование этих уроков, блестящих в теории хаоса, для раскрытия идей в других областях.

Мы будем делать это сегодня с помощью инструментария дифференциальных уравнений в Джулии.

Следующие примеры кода изменены из приведенной ниже документации.



Мне особенно нравится Джулия за такие задачи, за безумную простоту и удобство, которые она обеспечивает, когда дело доходит до научных вычислений.

Это только вопрос времени, когда он станет королем языков научных вычислений, обогнав как R, так и Python.

Я буду исследовать и объяснять, как работает код Джулии, строка за строкой, как он работает и как мы работаем со следующими обыкновенными дифференциальными уравнениями (ОДУ).

Меня особенно интересуют дифференциальные уравнения, поскольку они возникают практически во всех областях физики, и значительное количество исследований в области информатики посвящено разработке алгоритмов для их решения.

Система Лоренца

На первый взгляд они кажутся достаточно простыми.

Всего 3 линейных уравнения, верно? Правильно?

Как мы увидим позже, они описывают богатую систему соотношений, невидимых человеческому глазу.

Джулия Код

#In the Julia command line# 
using Pkg
Pkg.add("DifferentialEquations")
Pkg.add("Plots")
using DifferentialEquations
using Plots

Построение функции.

Обратите внимание на научную нотацию. Мы можем просто ввести символы непосредственно как переменные.

#Writing a parameterized function
function parameterized_lorenz!(du,u,p,t)
    
    x, y, z = u      #variables are part of vector array u 
    σ, ρ, β, = p     #coefficients are part of vector array p
    
    du[1] = dx = σ*(y-x)
    du[2] = dy = x*(ρ-z) - y
    du[3] = dz = x*y - β*z
    
end
.
.
.
< parameterized_lorenz! (generic function with 1 method) >

Подача функции на соответствующие входы.

#Initial conditions. x=1.0, y=0.0, z=0.0
u0 = [1.0, 0.0, 0.0]
#Timespan of the simulation. 100 in this case. 
tspan = (0.0, 100.0)
#Coefficients of the functions. 
p = [10.0, 28.0, 8/3]
#Feeding the inputs to the solver 
prob = ODEProblem(parameterized_lorenz!, u0, tspan, p)
sol = solve(prob)

Числовой вывод решателя:

3 элемента в численном решении следующие:

  1. «автоматическая интерполяция переключения ордеров»

Различные дифференциальные уравнения требуют различных компромиссов с точки зрения детализации и скорости. Используемый алгоритм может быть дополнительно сконфигурирован для различных комбинаций. В этом случае мы используем алгоритм по умолчанию.

2. "t: 1242 элемента вектора".

Массив векторов решений возвращает массив с детализацией 1242 временных шагов.

3. «u: вектор из 1242 элементов»

Точно так же массив векторов решений возвращает решения для каждой переменной в каждый момент времени.

retcode: Success
Interpolation: automatic order switching interpolation
t: 1242-element Vector{Float64}:
   0.0
   0.00020190266322792406
   0.0022209292955071643
   0.008092024514959373
   0.017621941854931468
   0.029910970730013992
   0.046119886045152166
   0.06633761003486688
   0.09145306990996667
   0.12186146157177435
   0.15775265249378573
   0.19917912194536747
   0.23735551696369733
   ⋮
  99.14124015678328
  99.23005160643454
  99.3178426998792
  99.39089107520101
  99.46975368622516
  99.54941086677553
  99.63372391029684
  99.70543691911548
  99.79383354268911
  99.90535455139218
  99.99916037878022
 100.0
u: 1242-element Vector{Vector{Float64}}:
 [1.0, 7.0, 0.0]
 [1.0157322915791744, 7.0034699810052485, 0.001424457841035686]
 [1.1713386826104133, 7.042273078544965, 0.016891096442499136]
 [1.60812010145366, 7.196207031607893, 0.07444157353741501]
 [2.2799815664984635, 7.5686982713074045, 0.2083531444865257]
 [3.110193937388093, 8.256085485026166, 0.4611897477586723]
 [4.20373407541577, 9.492552123156415, 0.9606655860252622]
 [5.661209508205023, 11.525389422787784, 1.9400058922165015]
 [7.749623385496089, 14.731104165543684, 3.9773408171241535]
 [10.806830843846502, 19.246723699443496, 8.351046992399898]
 [14.925214594925762, 23.671518939598823, 17.361338489508825]
 [18.585923260717127, 22.195831734843853, 31.716044613009842]
 [18.071904208857134, 12.17331991676706, 41.081896247131]
 ⋮
 [-9.183457977837184, -12.682011018250817, 17.64472825931256]
 [-12.321421331276596, -13.355987197450428, 26.506979959195224]
 [-10.49509552307489, -6.614740940422851, 31.26114857062963]
 [-6.412316590547746, -2.347421430959372, 28.77038858949141]
 [-3.2434635259271523, -1.1729141982820335, 24.469982740650746]
 [-1.942226612879031, -1.347551644896042, 20.53334362551526]
 [-1.7799839534618225, -1.9964184236190687, 17.090940494285892]
 [-2.2255447014614846, -2.969365928483296, 14.77883507206152]
 [-3.5474340484603615, -5.177797852868663, 12.944117866721324]
 [-7.152802707828483, -10.573895913693827, 14.093251538190133]
 [-11.671892462551858, -14.878852040827208, 22.03887865658922]
 [-11.706741167586157, -14.885100008376064, 22.141663040463467]

Если бы мы построили график нашего решения, это то, что мы увидим

plot1 = plot(sol, 
             vars = (1,2,3), 
             xlabel="x",  
             ylabel="y", 
             zlabel="z",
             label="x against y against z")
plot2 = plot(sol, label = ["x" "y" "z"])
plot(plot1, plot2, layout = (1, 2))

Эстетически интересный синий график — это знаменитый аттрактор Лоренца, также известный как график бабочки.

Здесь мы видим первые проблески природы хаотических систем.

Если мы посмотрим, как они меняются со временем, то очень трудно определить характер корреляций друг с другом, и все, что мы сможем увидеть, это то, что они хаотично меняются со временем.

Однако, если бы мы преобразовали график в график фазового пространства, мы начали бы получать более четкое представление о характере корреляций, если они вообще существуют.

Чрезвычайная чувствительность к начальным условиям

В следующем разделе мы продемонстрируем ключевую характеристику хаотических систем. Как изменение одной переменной приводит к гигантским изменениям в дальнейшем.

Мы будем варьировать только один параметр, rho.

И посмотрите, как резко изменится поведение.

u0 = [1.0, 7.0, 0.0]
tspan = (0.0, 100.0)
p = [13.0, 24.0, 7/3]
prob = ODEProblem(parameterized_lorenz!, u0, tspan, p)
sol = solve(prob);
plot5 = plot(sol, vars = (1,2,3), xlabel="x",  ylabel="y", zlabel="z",label="x against y against z",title="ρ = 24")
plot6 = plot(sol, label = ["x" "y" "z"], title="ρ = 24")
plot(plot5, plot6, layout = (1, 2))

Хорошо, это было не так уж плохо, верно? Сюжеты не выглядят, сильно отличаются от оригинала.

Теперь давайте изменим ро с 24 на 23

u0 = [1.0, 7.0, 0.0]
tspan = (0.0, 100.0)
p = [13.0, 23.0, 7/3]
prob = ODEProblem(parameterized_lorenz!, u0, tspan, p)
sol = solve(prob);
plot3 = plot(sol, vars = (1,2,3), xlabel="x",  ylabel="y", zlabel="z",label="x against y against z", title="ρ = 23")
plot4 = plot(sol, label = ["x" "y" "z"], title="ρ = 23")
plot(plot3, plot4, layout = (1, 2))

И бум, поведение резко меняется.

Система меняется от широко колеблющихся повсюду к поведению, которое схлопывается в единое состояние.

Давайте просто рассмотрим их рядом для простоты сравнения.

plot(plot3, plot4, plot5, plot6, layout=(2,2))

Всего одно небольшое изменение, и система демонстрирует другие характеристики. От диких, неуправляемых колебаний к стабильной и упорядоченной системе, которая переходит в устойчивое состояние.

Вывод

И в этом сила хаотической динамики.

Эти переменные могут не соответствовать действительным физическим условиям, существующим в реальности. Они вполне могут соответствовать сложному взаимодействию факторов, которые происходят только один раз за голубую луну.

И поиск ключевого фактора часто является святым Граалем в моделировании динамических систем.