Визначальна риса хаотичної системи — чутлива залежність від початкових умов: дві траєкторії, що починаються як завгодно близько, розходяться експоненційно. Показники Ляпунова кількісно характеризують цю швидкість розходження. Додатний максимальний показник Ляпунова (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 показує спад до сильно від'ємних значень.

← Діаграми біфуркацій Теорема Такенса →