Удобным практическим критерием стохастичности в заданной области фазового пространства служит численное определение показателей Ляпунова. Максимальный показатель $\sigma_{1}$ широко использовался в качестве критерия стохастичности в ряде работ (см., например, $[19,41-43,73,133,371])$.
Вычисление $\sigma_{1}$. Ранее было показано (п. 5.2б), что $\sigma\left(\boldsymbol{x}, \boldsymbol{w}_{0}\right)=$ $=\sigma_{1}(\boldsymbol{x})$ для почти всех касательных векторов $\boldsymbol{w}_{0}$. Поэтому для вычисления $\sigma_{1}$ начальный вектор $w_{0}$ можно выбирать произвольно. Интегрируя совместно уравнения (5.2.6) и (5.2.7), находим величину
\[
d(t)=w(t) \mid,
\]

где для удобства принято $d_{0}=d(0)=1$. Трудности возникают в том случае, когда $|w|$ растет экспоненциально, что приводит к переполнению и другим численным ошибкам. Чтобы избавиться от этих неприятностей, выберем определенный интервал времени $\tau$, как показано ${ }^{2}$ ) на рис. 5.6, и будем перенормировывать
1) Согласно данным работы [73], $\mu_{r} \sim \exp (-C \sqrt{K})$, где $C \sim$ const, за исключением специальных значений $K$, близких к (4.1.14).- Прим. ред.
2) Фактически на рис. 5.6 показан другой метод вычисления $\sigma_{1}$ – по двум близким траекториям системы (5.2.6), а не по одной траектории и линеаризованному уравнению (5.2.7), как в тексте. Использование линеаризованного уравнения имеет то преимущество, что его решение не зависит от модуля |w|. Однако в некоторых специальных случаях удобнее интегрировать две траектории (см., например, [470]).- Прим. ред.

$|w|$ на единицу в конце каждого такого интервала ${ }^{1}$ ). Таким способом будем последовательно вычислять величины
\[
\begin{aligned}
d_{k} & =\left|w_{k-1}(\tau)\right|, \\
w_{k}(0) & =\frac{w_{k-1}(\tau)}{d_{k}},
\end{aligned}
\]

где $w_{k}(\tau)$ получается интегрированием (5.2.6), (5.2.7) на интервале $\tau$ с начальными значениями $\boldsymbol{x}(k \tau)$, $\boldsymbol{w}_{k}(0)$. Если ввести величину

Рис. 5.6. Схема вычисления максимального показателя Ляпунова (по данным работы [19]).
$y=x+w, \tau-$ конечный интервал временн.

то из (5.2.8) получаем
\[
\sigma_{1}=\lim _{n \rightarrow \infty} \sigma_{n} .
\]

Для регулярного движения $\sigma_{1}=0$, в то время как на стохастической компоненте $\sigma_{1}>0$ и не зависит от $\boldsymbol{x}$. Описанный метод применим как для непрерывных траекторий, так и для отображений.

Следуя работе Бенеттина и др. [19], проиллюстрируем использование этого метода на примере модели Хенона-Хейлеса (п. 1.4а). На рис. 5.7 приведена зависимость $\sigma_{n}$ от $n \tau$ для шести траекторий с энергией $E=0,125$ (см. рис. 1.13). Три траектории ( $1-3$ ) находятся, очевидно, в области устойчивого движения, тогда как три другие (4-6) – в широкой стохастической области, которая хорошо видна на рис. 1.13. Как и следовало ожидать, для первых траекторий с увеличением $n$ величина $\sigma_{n}$, очевидно, стремится к
1) То есть сокращать длину вектора $w\left(|w| \rightarrow d_{0}=1\right)$, не изменяя его направления (см. [73]).- Прим. ред.

нулю $^{1}$ ), а для вторых, по-видимсму,-к одному и тому же значению $\left.{ }^{2}\right) \sigma_{1}$.

Используя этот метод, Бенеттин с соавторами нашли зависимость $\sigma_{1}$ от энергии $E$ в областях стохастического и регулярного движения (рис. 5.8). В последнем случае $\sigma_{n} \rightarrow 0$ (черточки на оси абсцисс). Сплошная кривая соответствует подгонке значений $\sigma_{1}>0$ к экспоненциальной зависимости. Используя соотно-

Рис. 5.7. Зависимость показателя $\sigma_{n}$ от времени для траекторий на регулярной ( $1-3$ ) и стохастической (4-6) компонентах (по данным работы [19]).

шение (5.2.27) и зависимость $1-\mu_{s}$ от $E$, полученную Хеноном и Хейлесом [188], Бенеттин с соавторами вычислили КС-энтропию $h(E)$ (пунктирная кривая на рис. 5.8). Они высказывают предположение, что $h(E)>0$ для всех $E>0$, что согласуется с нашими общими представлениями о свойствах движения, хотя и не видно из рис. 5.8 из-за различных приближений при вычислении $\sigma_{1}$. Это указывает на некоторый присущий методу недостаток. Возможно, например, что одна из траекторий на рис. 5.7, которая кажется регулярной, на самом деле лежит в узком стохастическом слое. Такую траекторию одинаково трудно обнаружить как непосредственно, так и при помощи настоящего метода, поскольку выход на очень низкое плато $\sigma_{1}$ может оказаться за пределами фактической
1) Грубо говоря, как $1 / n \tau$ [см. (5.3.4)].- Прим. ред.
2) Вообще говоря, это не всегда так (см., например, [73], рис. 3). В данном примере близкие значения $\sigma_{n}$ указывают на принадлежность всех трех траекторий к одной и той же стохастической компоненте.- Прим. ред.

длительности счета ${ }^{1}$ ). Мы вернемся к этим вопросам в гл. 6 при обсуждении более сложного движения в многомерных системах.
Вычисление всех показателей Ляпунова ${ }^{2}$ ). Следуя Бенеттину и др. $\left.[20]^{3}\right)$, покажем, каким образом вычисляется полный набор показателей Ляпунова в $M$-мерном фазовом пространстве. Ясно, что любая попытка определить $\sigma_{2}, \sigma_{3}$ и т. д., выбирая касательный неудачу из-за неустойчивости этих направлений, так что любые ошибки повернут в конце концов w (t) вдоль $\boldsymbol{e}_{1}$. Вместо этого выберем начальный базис из $p$ ортонормированных касательных векторов и численно определим $p$-мерный объем $V_{p}(t)$, заданный этими векторами. Отсюда можно найти показатель Ляпунова $\sigma_{1}^{(p)}$ порядка $p$ (5.2.14). Проделав эту процедуру для $p=1,2, \ldots, M$, из $(5.2 .15)$ определим все показатели $\sigma_{1}$, $\sigma_{2}, \ldots, \sigma_{M}$. Здесь, однако, возникает следующая трудность. В процессе движения углы между касательными векторами, вообще говоря, экспоненциально уменьшаются и численные ошибки резко возрастают. Поэтому в дополнение к перенормировке длины риодически ортогонализовать их.
Рис. 5.8. Зависимость максимального показателя $\sigma_{1}$ от энергии на стохастической (крестики) и регулярной (черточки на оси абсщисс) компонентах (по данным работы [19]).
Сплошная кривая – экспоненциальная интерполяция для $\sigma_{1}>0$; пунктирная кривая – КС-энтропия, пересчитанная по $\sigma_{1}$.

При этом новые векторы !должны лежать в том же подпространстве, что и старые.

Оказывается, что такую процедуру одновременного вычисления всех $p$-мерных объемов можно свести к расчету эволюции $M$ векторов и их специальной ортогонализации по методу Грама-
1) Заметим, что хотя действительно таким методом нельзя установить регулярность движения (интегрируемость) в некоторой области фазового пространства, однако вполне можно обнаружить хаотическую (случайную) компоненту движения, т. е. неинтегрируемость. В некотором смысле этот численный метод дополняет аналитический метод обратной задачи рассеяния (см. примечание редактора на с. 56), который позволяет, наоборот, доказать только интегрируемость.- Прим. ред.
2) Заметим, что для критерия стохастичности нужен по существу только максимальный показатель $\sigma_{1}$ ввиду неравенства $h \geqslant \sigma_{1}>0$. Иначе говоря, точное значение $\mathrm{KC}$-энтропии $h$, которое требует знания всех показателей, несущественно.- Прим. ред.
) См. также [502].-Прим. ред.

Шмидта. Вспомним, что $w_{k-1}(\tau)$ – касательный вектор $w_{k-1}(0)$ через время $\tau$. Вычислим сначала для каждого временного интервала $\tau$ следующие величины:
\[
\begin{aligned}
d_{k}^{(1)} & =\left|w_{k-1}^{(1)}(\tau)\right|, \\
w_{k}^{(1)}(0) & =\frac{w_{k-1}^{(1)}(\tau)}{d_{k}^{(1)}} .
\end{aligned}
\]

९ис. 5.9. Все три положительных показателя Ляпунова для модели Контопулоса (см. текст) (промежуточная область) (по данным работы [20]).

Затем для $j=2$, . ., $M$ найдем последовательно величины
\[
\begin{aligned}
\boldsymbol{u}_{k-1}^{(i)}(\tau) & =\boldsymbol{w}_{k-1}^{(j)}(\tau)-\sum_{i=1}^{j-1}\left(\boldsymbol{w}_{k}^{(i)}(0) \cdot \boldsymbol{w}_{k-1}^{(j)}(\tau)\right) \boldsymbol{w}_{k}^{(i)}(0) \\
d_{k}^{(j)} & =\left|\boldsymbol{u}_{k-1}^{(j)}(\tau)\right| \\
\boldsymbol{w}_{k}^{(j)}(0) & =\frac{\boldsymbol{u}_{k-1}^{(j)}(\tau)}{d_{k}^{(j)}} .
\end{aligned}
\]

Тогда в течение ( $k-1$ )-го интервала времени $\tau$ объем $V_{p}$ возрастает в $d_{k}^{(1)} d_{k}^{(2)} \ldots d_{k}^{(p)}$ раз. Отсюда, согласно (5.2.14), получаем
\[
\sigma_{1}^{(p)}=\lim _{n \rightarrow \infty} \frac{1}{n \tau} \sum_{i=1}^{n} \ln \left(d_{i}^{(1)} d_{i}^{(2)} \ldots d_{i}^{(p)}\right)
\]

Вычитая $\sigma_{1}^{(p-1)}$ из $\sigma_{1}^{(p)}$ и используя (5.2.15), находим $p$-й показатель Ляпунова
\[
\sigma_{p}=\lim _{n \rightarrow \infty} \frac{1}{n \tau} \sum_{i=1}^{n} \ln d_{i}^{(p)} .
\]

Это соотношение фактически и используется при численном моделировании.

Бенеттин и др. [20] нашли таким способом все показатели Ляпунова для нескольких гамильтоновых систем, включая 4- и 6-мерные отображения. Мы приведем их результаты для системы с тремя степенями свободы, которая исследовалась Контопулосом и др. [93]:
\[
H=\frac{1}{2}\left(p_{1}^{2}+q_{1}^{2}\right)+\frac{\sqrt{2}}{2}\left(p_{2}^{2} \div q_{2}^{2}\right)+\frac{\sqrt{3}}{2}\left(p_{3}^{2}+q_{3}^{2}\right)+q_{1}^{2} q_{2}+q_{1}^{2} q_{3} .
\]

Движение этой системы ограничено при $H<0,097$. Для $H=0,09$ фазовое пространство разделяется, грубо говоря, на три области: большая область стохастичности с $\sigma_{1} \approx 0,03 ; \sigma_{2} \approx 0,008$ и $\sigma_{3} \approx 0$; область регулярного движения ( $\sigma_{1}=\sigma_{2}=\sigma_{3}=0$ ) и промежуточная область. На рис. 5.9 представлены результаты вычисления первых трех показателей Ляпунова для начальных условий в промежуточной области. Видно, что сходимость имеет место для $\sigma_{1} \approx 3 \times 10^{-3}$ и $\sigma_{3}=0$, соответствующего направлению вдоль траектории. Поведение показателя $\sigma_{2}$ не вполне ясно.

Как мы увидим в § 6.2, эти результаты на самом деле обманчивы. Действительно, в системе с тремя степенями свободы первая и третья области должны быть связаны слабой диффузией Арнольда, благодаря которой траектория переходит из одной области в другую. Поэтому, по-видимому, и для промежуточной области $\sigma_{1} \approx 0,03$, а $\sigma_{2} \approx 0,008$, что противоречит данным на рис. 5.9. Это еще раз указывает на основную трудность численного определения показателей Ляпунова: не существует априорного условия для определения достаточного числа итераций $n$. Поэтому при численном моделировании необходимо использовать и другие методы, такие, например, как метод сечения Пуанкаре ${ }^{1}$ ).