Визначальна риса хаотичної системи — чутлива залежність від початкових умов: дві траєкторії, що починаються як завгодно близько, розходяться експоненційно. Показники Ляпунова кількісно характеризують цю швидкість розходження. Додатний максимальний показник Ляпунова (MLE) — середня експоненційна швидкість, з якою дві початково близькі траєкторії віддаляються одна від одної — це еталонне означення хаосу, що перетворює інтуїтивне уявлення про «ефект метелика» на точне, обчислюване число.

Коротко: Показники Ляпунова вимірюють, як швидко розходяться близькі траєкторії системи; додатне значення означає хаос. Стаття розглядає повний спектр показників, формулу розмірності Каплана–Йорка для фрактальних атракторів, три методи обчислення (пряма формула, Вольф–Свінні–Вастано, Грама–Шмідта), метод вкладення Такенса для реальних даних та зв'язок із межею передбачуваності, наприклад у прогнозі погоди.

1. Означення та геометричний зміст

Розгляньмо гладку динамічну систему ẋ = f(x) (потік) або xₙ₊₁ = F(xₙ) (відображення). Починаючи з x₀, збуримо початкову умову крихітним вектором δ₀ з |δ₀| = ε. Через час t (або n ітерацій) збурення еволюціонує до δ(t) через лінеаризовану (дотичну) динаміку:

δ(t) = DFᵗ(x₀) · δ₀,    DFᵗ ≡ ∂Fᵗ/∂x|_{x₀}

Показник Ляпунова для напрямку δ₀ дорівнює:

λ(x₀, δ₀) = lim_{t→∞} (1/t) · ln |δ(t)| / |δ₀|

За мультиплікативною ергодичною теоремою Оселедця (1968) для майже всіх x₀ (відносно ергодичної інваріантної міри) ця границя збігається до одного з щонайбільше d значень (де d — розмірність простору станів), незалежно від початкового напрямку δ₀ для множини напрямків повної міри. Ці значення λ₁ ≥ λ₂ ≥ … ≥ λ_d утворюють спектр Ляпунова.

Знак показника Геометричний зміст Типова система
λ₁ > 0 Експоненційне розходження — хаос Лоренц, логістичне r > r∞
λ₁ = 0 Нейтральний / граничний — біфуркація або квазіперіодичність Граничний цикл, тори КАМ
λ₁ < 0 Збіжність до атрактора Стійка нерухома точка, стійкий граничний цикл
Σλᵢ < 0 Стиснення фазового об'єму (дисипативна) Усі фізичні атрактори
Σλᵢ = 0 Збереження об'єму (гамільтонова) Маятник, задача N тіл

2. Розмірність Каплана–Йорка

Фрактальна (інформаційна) розмірність дивного атрактора пов'язана зі спектром Ляпунова через гіпотезу Каплана–Йорка (1979), чисельно підтверджену для багатьох систем:

D_KY = j + Σᵢ₌₁ʲ λᵢ / |λⱼ₊₁|

де j — найбільший індекс, такий що Σᵢ₌₁ʲ λᵢ ≥ 0. Для системи Лоренца з σ=10, ρ=28, β=8/3 спектр приблизно дорівнює λ₁ ≈ 0.906, λ₂ = 0, λ₃ ≈ −14.57, що дає D_KY ≈ 2 + 0.906/14.57 ≈ 2.062. Знаменитий атрактор Лоренца справді є фрактальною множиною нецілочислової розмірності.

Обмеження на спектр Ляпунова

• Для автономних неперервних систем один показник завжди дорівнює 0 вздовж напрямку потоку.
• Для гамільтонових систем показники утворюють симетричні пари: λᵢ = −λ_{d+1−i}.
• Для дисипативної системи: Σλᵢ = ∫ ∇·f dμ (усереднена за часом дивергенція векторного поля).
• Для системи Лоренца: Σλᵢ = −(σ + 1 + β) = −13.67.

3. Обчислення MLE з одновимірного відображення

Для одновимірного відображення xₙ₊₁ = F(xₙ) MLE зводиться до:

λ = lim_{N→∞} (1/N) · Σₙ₌₀^{N−1} ln |F′(xₙ)|

Для логістичного відображення F(x) = rx(1−x):

F′(x) = r(1 − 2x), λ(r) = lim_{N→∞} (1/N) Σ ln |r(1 − 2xₙ)|
Параметр r Режим MLE λ
r = 2.0 Стійка нерухома точка ln(r−2) → −∞ для нерух. точки
r = 3.2 Цикл періоду 2 ≈ −0.345
r = 3.5 Цикл періоду 4 ≈ −0.163
r ≈ 3.5699 Початок хаосу (r∞) = 0
r = 3.8 Хаос (внутрішня криза) ≈ +0.433
r = 4.0 Повний хаос (спряженість із наметовим відображенням) = ln 2 ≈ 0.693

Точки біфуркації логістичного відображення — це точно нулі λ(r), що підтверджує: λ = 0 є водночас необхідною й достатньою умовою для біфуркацій у просторі параметрів цього сімейства.

4. Алгоритм Вольфа–Свінні–Вастано

Для неперервного потоку (Лоренц, Реслер тощо) або експериментального часового ряду пряме диференціювання недоступне. Алгоритм WSV (Wolf, Swift, Swinney & Vastano, 1985) оцінює максимальний показник Ляпунова з єдиної траєкторії:

Алгоритм WSV — псевдокод

1. Проінтегрувати орбіту, щоб знайти опорну траєкторію x(t)
2. Обрати близьку точку x̃(0) з |δ₀| = |x̃(0) − x(0)| = ε (малою)
3. Проінтегрувати обидві траєкторії протягом часу еволюції T_ev
4. Виміряти нове розділення d₁ = |x̃(T_ev) − x(T_ev)|
5. Накопичити: L += ln(d₁ / ε)
6. Перемасштабувати: замінити x̃ → x(T_ev) + ε · δ₁ / d₁  (ренормалізація, збереження напрямку)
7. Повторити з кроку 3; після N циклів:
   λ₁ = (1 / N·T_ev) · L

Ключові параметри: T_ev має бути достатньо коротким, щоб уникнути дотиків орбіт, але достатньо довгим, щоб зменшити статистичний шум; ε має бути малим відносно розміру атрактора, але великим відносно шуму чисел з рухомою комою. Перевірені типові значення: ε ~ 10⁻⁸, T_ev ~ 0.5–2.0 одиниці часу.

5. Метод повного спектра Грама–Шмідта

Щоб обчислити всі d показників одночасно, інтегруйте дотичні (варіаційні) рівняння паралельно з траєкторією. Дотичні вектори мають періодично ортонормуватися за Грамом–Шмідтом, щоб запобігти їхньому злипанню в максимальний напрямок.

dδᵢ/dt = Df(x(t)) · δᵢ для i = 1, …, d

Після кожного кроку ренормалізації в момент τ записуйте логарифм розтягу кожного базисного вектора до ортонормування. Накопичена сума / загальний час → λᵢ. Це основа алгоритму Бенеттіна та ін. (1980) і його чисельної реалізації у більшості пакетів наукових обчислень.

QR-інтерпретація Грама–Шмідта

На кожному кроці d дотичних векторів утворюють матрицю Φ. Розклад Φ = QR з ортогональною Q; діагональні елементи rᵢᵢ матриці R дають розтяг уздовж кожного напрямку Оселедця. Після N кроків:
λᵢ = (1 / N·τ) · Σₖ ln |rᵢᵢ⁽ᵏ⁾|

6. Показники Ляпунова з експериментальних часових рядів

Коли доступний лише скалярний ряд вимірювань {s₁, s₂, …, sₙ} (наприклад, ЕКГ, ЕЕГ, ціни акцій), фазовий простір спершу слід відновити за допомогою вкладення із затримкою за Такенсом:

x(t) = [s(t), s(t + τ), s(t + 2τ), …, s(t + (m−1)τ)] ∈ ℝᵐ

За теоремою Такенса відновлений атрактор дифеоморфний оригінальному, якщо m ≥ 2d_A + 1 (де d_A — розмірність атрактора). Rosenstein та ін. (1993) і Kantz (1994) пропонують чисельно стійкі алгоритми для MLE на вкладених даних — метод Розенштейна відстежує розходження найближчих сусідів для всіх початкових умов одночасно:

d_j(i) = 〈ln Δⱼ(i)〉_j,   нахил d_j(i) відносно i·Δt ≈ λ₁

Потрібні параметри: розмірність вкладення m (з методу хибних найближчих сусідів, FNN), затримка τ (з нуля автокореляції або першого мінімуму взаємної інформації), вікно Тейлера W (для виключення часових сусідів).

7. Зв'язок з ентропією та передбачуваністю

Ентропія Колмогорова–Сіная (КС) h_KS кількісно характеризує швидкість продукування інформації динамічною системою. Формула Пезіна прямо пов'язує її з додатними показниками Ляпунова:

h_KS = Σ_{λᵢ > 0} λᵢ (тотожність Пезіна)

Горизонт передбачуваності: початкова невизначеність ε зростає до невизначеності порядку одиниці через час t_pred ≈ −ln(ε) / λ₁. Для системи Лоренца з λ₁ ≈ 0.906 та ε = 10⁻⁵ горизонт передбачуваності становить ≈ 13 одиниць часу — що пояснює, чому чисельний прогноз погоди стає ненадійним за межами ~2 тижнів.

Система MLE λ₁ Передбачуваність (ε=10⁻⁵)
Лоренц (σ=10, ρ=28, β=8/3) ≈ 0.906 ≈ 12.7 одиниці часу
Реслер (a=0.2, b=0.2, c=5.7) ≈ 0.071 ≈ 162 одиниці часу
Логістичне r = 4.0 = ln 2 ≈ 0.693 ≈ 16.6 ітерацій
Подвійний маятник (енергія ≈ 10E_min) ≈ 3.5 рад⁻¹ ≈ 3.2 радіан-часу
Сонячна система (зовнішні планети) ≈ 1/(5 млн рр.) ≈ 58 млн рр. (1/λ₁)

8. Інтерактив: MLE для логістичного відображення в реальному часі

Перетягуйте повзунок r, щоб спостерігати оновлення показника Ляпунова в реальному часі. Верхня панель показує діаграму орбіти (останні 200 ітерацій після прогріву на 400); нижня панель будує λ(r), обчислений на льоту за формулою суми логарифмів похідних. Додатне λ — червоне (хаос), нульові перетини позначають біфуркації, а від'ємне λ — зелене (періодичність).

λ ≈ ?

Режим повного сканування проходить r від 2.5 до 4.0 й малює λ(r) як криву — точки біфуркації з'являються як різкі перетини нуля, а вікно періоду 3 при r≈3.828 показує спад до сильно від'ємних значень.