1. Потік

Розглянемо автономну систему звичайних диференціальних рівнянь:

(1)     dx_j / dt = f_j (x₁ , x₂ , …, x_n)

при дійсних x_j , функціях f_j : Rⁿ → R ,   j = 1, 2, …, n.

Словом «автономність» позначають незалежність у явному вигляді правої частини рівнянь (1) від аргумента t. Систему рівнянь (1) ще записують у такому вигляді:

(2)     \(x' = f(x),\)

де \(x\) — вектор-стовпчик з координатами \(x_1,~x_2,~...,~x_n,\)
\(f(x)\) — вектор-стовпчик з такими координатами: \(f_1(x_1,~x_2,~...,~x_n)\), \(f_2(x_1,~x_2,~...,~x_n)\), …, \(f_n(x_1,~x_2,~...,~x_n)\).

Природно розглядати такі праві частини рівняння (2), при яких роз’язок рівняння існує та єдиний хоча б в околі вибраних початкових значень x. Таким умовам задовольняють, наприклад, неперервно диференційовні за координатами x функції f — див. традиційні підручники з теорії диференціальних рівнянь для студентів вищих навчальних закладів (В.И.Арнольд. Обыкновенные дифференциальные уравнения. М: Наука, 1971, 240 с. та інші) або адаптований логічно послідовний виклад теоретичних основ шкільного курсу у посібнику: О.Б.Рудик. Початки алгебри, аналізу, аналітичної геометрії і теорії ймовірностей. Тернопіль: «Навчальна книга–Богдан», 2005, 415 с. При розв'язуванні нелінійних рівнянь малоймовірно отримати розв'язок у термінах елементарних функцій. Але можна використати методи наближеного чисельного інтегрування. Наприклад, модифікацію Мерсона методу Рунге-Кутти. Та часто обмежуються такою різнецевою для рівняння (2): \[x(t+h)=x(t)+h\cdot f(x(t))\] для достатньо малого приросту часу \(h\).

Означення 1. Запровадимо такі поняття для рівняння (2).

1. Потоком, породженим полем напрямків f називають множину {f ^t} відображень

   f ^t : Rⁿ → Rⁿ,

   що початковій величині x ставлять у відповідність розв’язок рівняння (2) за час t:

   (f ^tx)' = f (f ^tx)   при   f ⁰x = x.

2. Означимо дію потоку {f ^t} на С(Rⁿ) — просторі неперервних на Rⁿ функцій — таким чином:

   (f ^tu) (x) = u (f ^tx).

3. Означимо дію потоку {f ^t} на просторі обмежених мір µ, заданих на Rⁿ, таким чином, щоб при довіль­ній не­пе­рер­вній обмеже­ній функ­ції u ∈ С(Rⁿ) справджувалися такі рівності:

   ∫ u(x) d(f ^tµ)(x) = ∫ (f ^tu)(x) dµ(x) = ∫ u(f ^tx) dµ(x).

Зауваження 1. Якщо рівняння (2) моделює певний природний чи суспільний процес, то з конкретно-наукових міркувань бажано, щоб траекторії-розв’язки для додатного часу t не прямували до нескінченості, а належали до певної замкненої обмеженої множини.

Означення 2. Потоки {f ^t} і {g^t} називають еквівалентними (спряженими, подібними), якщо існує взаємно однозначне відображення

що відображає потік {f ^t} у потік {g^t} при заміні координат h. Інакше кажучи,

при довільних t ∈ R, x ∈ Rⁿ. При цьому потоки називають:

• лінійно еквівалентними, якщо відображення h є лінійним оператором;

• аналітично еквівалентними, якщо відображення h є аналітичним, тобто h можна подати степеневим рядом приростів координат x;

• диференціально еквівалентними, якщо відображення h є дифероморфізмом, тобто диференційовним відображенням;

• топологічно еквівалентним, якщо відображення h є взаємно неперервним.

Зауваження 2. З лінійної еквівалентності потоків випливає їхня аналітична еквівалентність, з аналітичної — диференціальна, з диференціальної — топологічна.

2. Стійкість розв'язку

Надалі через |…| в R позначимо модуль (абсолютну величину), в Rⁿ — довільну норму (довжину вектора). Наприклад,

при p ≥ 1. При p = 2 це збігається зі звичайною евклідовою нормою.

Означення 3. Розв’язок f ^tx рівняння (2) називають:

• стійким за Ляпуновим, якщо вибором достатньо малого околу x початкових значень y можна зробити як завгодно малим відхилення |f ^ty – f ^tx| рівномірно для всіх невід’ємних t:

  ∀ε > 0 ∃δ > 0 ∀y ∈ Rⁿ : |y – x| < δ ⇒ ∀t ≥ 0 | f ^ty – f ^tx| < ε ;

• асимптотично стійким, якщо:

  • він є стійким за Ляпуновим;

  • існує окіл x початкових значень y, при яких відстань між розв’язками з початковими значеннями y та x збігається до нуля:

    ∃Δ > 0 ∀y ∈ Rⁿ : |y – x| < Δ ⇒ | f ^ty – f ^tx| → 0 при t → + ∞ ;

• асимптотично стійким у цілому, якщо:

  • він є стійким за Ляпуновим;

  • для довільного початкового значення y відстань між розв’язками з початковими значеннями y та x збігається до нуля:

    ∀y ∈ Rⁿ : | f ^ty – f ^tx| → 0 при t → + ∞ ;

• експотенціально стійким, якщо

  ∃ (Δ > 0, M > 0, γ > 0) ∀y ∈ Rⁿ : |y – x| < Δ ⇒ ∀t ≥ 0 | f ^ty – f ^tx| ≤ Me ^{– γ t} ;

• експотенціально стійким у цілому, якщо

  ∃ (M > 0, γ > 0) ∀y ∈ Rⁿ ∀t ≥ 0 | f ^ty – f ^tx| ≤ Me ^{– γ t} .

Зауваження 3. З експотенціальної стійкості випливає асимптотична стійкість, з якої, у свою чергу, випливає стійкість. Зворотні імплікації не справджуються.

Зауваження 4. Поняття стійкості можна узагальнити на випадок неавтономних систем, формально вказавши залежність від змінної t та вибраного початкового значення t₀ цієї змінної.

Зауваження 5. Нехай квадратна матриця A за домомогою матриці С зводиться до нормальної форми Жордана B:

A = СBС ^–1.

Тоді розв’язок рівняння (2) при f (x) = Ax має таке подання:

e^Atx = Сe^BtС ^–1x.

1. Якщо дійсні частини усіх власних чисел A від’ємні, то розв’язок x = 0 є експотенціально стійким.

2. Якщо дійсні частини усіх власних чисел A недодатні, а деякі з них (дійсних частин) дорівнють нулю при розмірності 1 відповідних підпросторів Rⁿ, то розв’язок x = 0 є стійким, але не є асимптотично стійким.

3. Якщо дійсні частини усіх власних чисел A недодатні, а деякі з них (дійсних частин) дорівнють нулю при розмірності відповідних підпросторів Rⁿ більшій за 1, то розв’язок x = 0 не є асимптотично стійким.

4. Якщо дійсна частина хоча б одного з власних чисел A додатна, то розв’язок x = 0 не є стійким.

Означення 4. Точку x простору Rⁿ називають стаціонарною (нерухомою) точкою потоку {f ^t}, якщо f ^tx = x при довільних дійсному t та x з Rⁿ. Розглянемо власні числа матриці часткових похідних ∂f/∂x = ||∂f_j/∂x_k|| у стаціонарній точці x за умови, що відповідні похідні існують:

• якщо всі власні числа ∂f/∂x дійсні й розташовані по один бік від нуля на числовій прямій, то стаціонарну точку називають вузлом;

• якщо всі власні числа ∂f/∂x дійсні, але розташовані по різні боки від нуля на числовій осі, то стаціонарну точку називають сідлом;

• якщо серед власних чисел ∂f/∂x є такі, що мають відмінні від нуля уявні частини, то стаціонарну точку називають фокусом.

Зауваження 6. Необхідною і достатною умовою нерухомості точки x відносно потоку {f ^t} є справдження рівності: f(x) = 0.

Нерухому точку математичної моделі еволюції складної системи традиційно пов’язують зі станом рівноваги (стійкої чи нестійкої).

Зауваження 7. Анрі Пуанкаре у своїй дисертаційній роботі довів, що при:

• аналітичній залежністі поля напрямків f(x) від координат x;

• відсутності резонансу власних чисел матриці ∂f(x)/∂x при нерухомій точці x потоку, тобто при справдженні усіх нерівностей такого вигляду λ_j ≠ (λ | m), де:
  λ₁, λ₂, …, λ_n — власні числа матриці ∂f(x)/∂x — координати вектора λ;
  m₁, m₂, …, m_n — невід’ємні цілі числа, сума яких не менша, ніж 2, — координати вектора m;
  (λ | m) = λ₁m₁ + λ₂m₂ + ... + λ_nm_n — скалярний добуток векторів λ і m

формальною заміною змінних рівняння (2) в околі нерухомої точки x зводиться до лінеаризованої задачі:

(3)     \(y' = {\large ∂ f (x)\over\large ∂ x}\cdot y\)

в околі початку координат (про доведення див. В.И.Арнольд. Дополнительные главы теории обыкновенных дифференциальных уравнений. М.: Наука, 1978, 304 с.).

При локальній еквівалентності рівнянь (2) і (3):

• якщо дійсні частини всіх власних чисел матриці ∂f(x)/∂x від’ємні, то стаціонарна точка x є стійкою;

• якщо дійсна частина хоча б одного власного числа матриці ∂f(x)/∂x додатна, то стаціонарна точка x є не стійкою;

В усіх інших випадках для встановлення стійкості потрібно досліджувати вищі похідні поля напрямків f.

3. Дотичне відображення

Означення 5. Позначимо через x(t) вектор-стовпчик з координатами x₁(t), x₂(t), ... , x_n(t), що є розв’язком рівняння (2).

1. Дотичним відображенням потоку {f ^t} на відрізку [t₀; t₁], називають лінійний опера­тор, що подають такою матрицею часткових похідних: \[{∂x(t_1)\over ∂x(t_0)} = \left|\left| ∂x_j(t_1)\over ∂ x_k(t_0)\right|\right|.\]

2. Характеристичним показником Ляпунова χ напряму вектора v при дотичному відображенні ∂x(t)/∂x(0) потоку {f ^t} називають таку верхню часткову границю:

   \[χ(v, x(0)) = \varlimsup\limits_{t\to +\infty} {1\over t} \ln\left|{∂x(t) \over ∂x(0)} v\right|.\]

Зауваження 8. Згідно з правилом диференціювання складеної функції маємо: \[{∂x(t+dt)\over ∂x(t_0)} = {∂x(t+dt)\over ∂x(t)} \cdot {∂x(t)\over ∂x(t0)}.\] Тут знаком · позначено операцію множення матриць. Враховуючи співвідношення: \[x(t + dt) = x(t) + f (x(t)) dt + o(dt),\] отримаємо таке рівняння для матриці дотичного відображення: \[{d \over dt} {∂x(t)\over ∂x(t0)} = {∂f\over ∂x} (x(t)) · {∂x(t)\over ∂x(t0)}.\] Згідно з правилом диференціювання складеної функції маємо: \[ {d f(x(t))\over dt} = {∂f\over ∂x}(x(t)) · f(x(t)).\] Таким чином, напрям потоку зберігається при дотичному відображенні.

Якщо розв'язки рівняння (2) з плином часу \(t \to + \infty\) наближаються до деякої нерухомої точки потоку або до циклу (замкненої кривої, яку поток відображає саму у себе), поведінка розв'язку стає добре прогнозованою. Та існують динамічні системи, у яких такі довгострокові прогнози неможливі. Наприклад, у більярдах з опуклими всередину краями траєкторії руху куль без тертя "розбігаються і перемішуються". У цих випадках використовують ергодичну теорію — розділ теорії динамічних систем для випадків, коли середні значення функції вздовж траекторій збігаються до інтегралу цієї функції за певним розподілом (мірою) майже всюди за цією мірою. У цьому випадку фазовий простір (множину, на якій діє потік чи відображення) розбивають на ергодичні компоненти — підмножини фазового простору, кожна з яких є носієм міри, інваріантної відносно дії потоку, і яку неможливо подати об'єднанням хоча б двох множин з такою властивістю. Ознайомитися з основами ергодичної теорії можна за таким виданням:

• Корнфельд И.П., Синай Я.Г., Фомин С.В. Эргодическая теория. — М.: Наука. Главная редакция физико-математической литературы, 1980. 384 стр. c илл.

після ознайомлення з традиційним курсом функціонального аналізу. Для існування інваріантних мір та ергодичних компонент (і для потоків, і для відображень) є достатня умова: наявність обмеженої замкненої множини, в яку потрапляють усі траєкторії, але жодна не виходить з цієї множини. Це узгоджено з висловлюванням Зауваження 1, сформульованого раніше.

Характеристичні показники Ляпунова явно описують асимптотичну швидкість збігання чи розбігання траєкторій у певних напрямках. Ґрунтовна публікація:

• Оселедец В.И. Мультипликативная эргодическая теорема. Харак­терис­тические показатели Ляпунова // Труды Московского матем. общ-ва, 1968, том 19, c. 197–231.

містить алгоритм обчислення усіх характеристичних показників Ляпунова. Цей алгоритм легко програмно втілити довільною мовою програмування. Він має такий вигляд.

1. Вибираємо точку x(0) ергодичної компоненти динамічної системи. В ході чисельного експерименту це здійснюють, обчислюючи розв’язок рівняння (2) на досить великому інтервалі для довільного початкового значення. Таким чином можна знайти точку лише в околі атрактора — множини, до яких збігаються розв'язки (2) (які "притягують траєкторії") .

2. Вибираємо X(0) — довільну невироджену квадратну матрицю розміру n×n. Наприклад, матрицю, діагональні елементи якої дорівнюють 1, а решта — 0.

3. Вибираємо T — довільне додатне дійсне число.

4. Для натурального k = 1, 2, ... робимо такe.

   • Розглянемо розширення у дотичне розшарування рівняння (2) на проміжку [(k – 1)T; kT). Інакше кажучи, на цьому проміжку розв’яжемо таку систему рівнянь: \[x'(t) = f(x);\qquad X'(t) = {∂f\over ∂x} (x(t)) \cdot X(t).\] і знайдемо \[X(kT-0) = {∂x(kT)\over ∂x((k – 1)T)} · X((k – 1)T).\]

   • Знаходимо матрицю С ∈ SO(n, R), при якій матриця D = С ^– 1 · X(kT – 0) має верхньо-трикутний вигляд. Для матриць з дійсними елементами умова С ∈ SO(n, R) означає справдження таких рівностей:

     CC ^t = C ^tC = diag(1, 1, …, 1).

     Інакше кажучи, рядки (стовпчики) матриці С утворюють ортонормовану базу Rⁿ: скалярний добуток кожного такого вектора з собою дорівнює 1, а з іншим таким вектором — 0.
     
     Опишемо детально одну з можливих процедур визначення матриці С, скориставшись такими позначеннями:
     с_j — вектор Rⁿ, що є рядком j матриці С;
     x_j — вектор Rⁿ, що є рядком j матриці X(kT – 0);
     d_jl — елемент матриці D, розташований у рядку j та стовчику l, d_jl = 0 при j > l;
     (u|v) — скалярний добуток (сума добутків відповідних координат) векторів u і v.
     
     Співвідношення C·D = X(kT – 0) запишимо за допомогою запроваджених позначень. При l = 1, 2, …, n:
     

     l
     ∑	djl cj = xl .
     j = 1

     При l = 1 маємо:

     • d₁₁ = (x₁ | x₁)^1/2 — при обчисленні наближень характеристичних показників Ляпунова на цю величину є посилання як на d₁₁(kT );

     • с₁ = (d₁₁)^{– 1}x₁.

     При l = 2, 3, …, n маємо:

     • d_ml = (c_m | x_l) при m = 1, 2, …, l – 1, бо коефіцієнти розкладу вектора за ортонор­мованою базою збігаються зі скалярними добутками вектора і векторів бази;

     • u_l = x_l – d_1l c₁ – d_2l c₂ – ... – d_{(l – 1) l} c_l – 1 — складова x_l, перпендикулярна до c₁, c₂, c_l – 1;

     • d_ll = (u_l | u_l)^1/2 — при обчисленні наближень характеристичних показників Ляпунова на цю величину є посилання як на d_ll(kT );

     • с_l = (d_ll)^{– 1}u_l.

   • Обчислимо наближення характеристичних показників Ляпунова

     	k		k
     χl(kT ) = (kT ) – 1	∑	ln dl( jT ) = (kT ) – 1 ln (	П	dll( jT )).
     	j = 1		j = 1

   • Покладемо X(kT) = С.

Пункт 4 виконують до тих пір, поки певна кількість цифр після десяткової коми у десяткових записах наближень характеристичних показників Ляпунова не стануть сталими. При виконанні алгоритму послідовно обчислюють верхньо-трикутні матриці — добутки по три множники, рахуючи справа наліво, — такого виразу: \[\cdots~\cdot X^{-1}(2T) \cdot {∂x(2T)\over ∂x(T)} \cdot X(T) \cdot X^{-1} (T) \cdot {∂x(T)\over ∂x(0)} \cdot X(0).\] Зауваження 9. Мультиплікативну ергодичну теорему доведено для так званих коциклів, частковим випадком яких є дотичне перетворення потоку або неперервно диференційованого відображення. В останньому випадку немає потреби розв'язувати диференціальне рівняння для обчислення матриці дотичного відображення — достатньо явно вичначити матрицю похідної відображення.

Зауваження 10. Нехай потік {f ^t} з неперервно диференційовним f має інваріантну замкнену обмежену множину.

1. Для траекторій всередині цієї множини:

   • або траекторія наближається до нерухомої точки і характеристичний показник Ляпунова напряму потоку збігається з дійсною частиною відповідного власного значення матриці ∂f/∂x, обчисленій у цій нерухомій точці;

   • або величина | f | відділена від нуля та нескінченності і характеристичний показник Ляпунова на­пряму по­то­ку дорів­нює нулю.

2. Якщо на ергодичній компоненті фазового простору характеристичні показники Ляпунова усіх векторів v для дотичного відображення ∂x(t)/∂x(0) потоку {f ^t} від’ємні, можливо, за виключенням поля напрямів та сум його з векторами, чиї показники від’ємні, то розв’язок x(t) рівняння (5) є експотенціально стійким.

3. Якщо для ергодичної компоненті фазового простору існує вектор, чий характеристичний показник Ляпунова для дотичного відображення ∂x(t)/∂x(0) потоку {f ^t} додатний, то розв’язок x(t) рівняння (2) не є стійким.

При чисельному інтегруванні при наявності потоків з додатними характеристичними показниками спочатку різницева схема спотворює напрям потоку, і лише потім відбувається збіг наближень показників Ляпунова. Інакше кажучи, порахувати можна характеристичні показники лише різнецевої схеми, а не потоку. Це було зауважено в роботі:

• Rudyk A.B. The Lyapunov Characteristic Exponents of Dissipative Systems with a Strange Attractor. К.: 1983. — 17 с. Препpинт / АН Украины, Ин-т теор. физики; 83–156E

при обчисленні характеристичних показників Ляпунова системи трьох звичайних диференціальних рівнянь, побудованої як приклад до такого твердження: потоки, чиї поля напрямків є зваженою сумою полів потоків з різними стійкими циклами, можуть демонструвати складну поведінку: каскад подвоєння періоду стійкого циклу, наявність атрактора з дробовою розмірністю тощо залежно від коефіцієнтів зваженої суми.