CC BY
44
Изв. вузов «ПНД», т. 16, № 3, 2008
УДК 517.9
ХАОТИЧЕСКАЯ ДИНАМИКА МОДЕЛИ ХАНТА -ИСКУССТВЕННО СКОНСТРУИРОВАННОЙ ПОТОКОВОЙ СИСТЕМЫ С ГИПЕРБОЛИЧЕСКИМ АТТРАКТОРОМ
Ю.С. Айдарова, С.П. Кузнецов
Исследуется численно хаотическое поведение, обусловленное присутствием гиперболического странного аттрактора типа Плыкина в модели Ханта - искусственно сконструированной динамической системе с непрерывным временем. Приводятся портреты аттрактора, графики реализаций порождаемого системой хаотического сигнала, иллюстрируется присущая хаосу чувствительная зависимость траекторий от начальных условий. Представлены также количественные характеристики аттрактора - показатели Ляпунова и оценка размерности. Обсуждается символическая динамика на аттракторе, найдены и проанализированы некоторые неустойчивые периодические орбиты, принадлежащие аттрактору.
Введение
В математической теории динамических систем вводится в рассмотрение класс однородно гиперболических хаотических аттракторов [1-9]. В таком аттракторе все траектории являются седловыми, обладают устойчивыми и неустойчивыми многообразиями, размерность которых одна и та же для всех точек на аттракторе, и эти устойчивые и неустойчивые многообразия не имеют касаний между собой. Динамика на странном гиперболическом аттракторе демонстрирует сильные хаотические свойства и допускает подробный математический анализ. Такие аттракторы обладают свойством грубости или структурной устойчивости, которое заключается в нечувствительности устройства аттрактора по отношению к вариациям функций и параметров в определяющих уравнениях. Основные положения соответствующей математической теории были сформулированы около 40 лет назад, однако до последнего времени считалось, что однородно гиперболические странные аттракторы не характерны для реальных систем со сложной динамикой. В учебниках и монографиях в качестве примеров традиционно фигурируют искусственные математические конструкции, такие как аттрактор Плыкина и соленоид Смейла-Вильямса. Названные примеры относятся к системам с дискретным временем - итерируемым отображениям специального вида. Аттрактор Смейла-Вильямса получается при отображении тороидальной области в себя в фазовом пространстве с размерностью 3 или более. Аттрактор Плыкина имеет место в некотором сложно организованном отображении, которое определено в ограниченной области на плоскости с тремя дырками.
176
В прикладных дисциплинах, физике и технике чаще всего приходится иметь дело с системами, функционирующими в непрерывном времени.
Переход от системы с дискретным временем, описываемой отображением xn+i = f (xn), к системе с непрерывным временем - потоковой системе, или просто потоку, математики называют надстройкой (suspension) [2-7]. К набору переменных х, определяющих состояние исходной системы, добавляется еще одна переменная t, отвечающая непрерывному времени, и вводится расширенное фазовое пространство {x, t}. Далее строится векторное поле v, зависящее от исходных переменных x и дополнительной переменной t. При этом за некоторый определенный интервал времени At = T состояние, эволюционирующее в соответствии с дифференциальным уравнением X = v(x, t), должно претерпевать трансформацию, в точности соответствующую применению отображения xn+i = f(xn). За пределами интервала T векторное поле v(x, t) считается периодической функцией времени. Процедура применима, если отображение xn+1 = f (xn) обратимое. С физической точки зрения, это следует трактовать как построение неавтономной системы, в которой за один период внешнего воздействия T реализуется преобразование состояния, заданное исходным отображением [8-12].
Недавно предложена и исследована, в том числе экспериментально, система с аттрактором, который отвечает надстройке над соленоидом Смейла-Вильямса [13, 14]. Система конструируется на основе двух неавтономных осцилляторов ван дер Поля, которые становятся активными попеременно и передают возбуждение друг другу так, что преобразование фазы колебаний за полный цикл описывается растягивающим отображением окружности (отображение Бернулли). Результаты проверки на компьютере условий теоремы, гарантирующей существование гиперболического аттрактора, представлены в работах [15, 16]. Другие варианты подобных схем на основе неавтономных и автономных осцилляторов обсуждаются в работах [17-19].
Для аттрактора Плыкина пока не известны конкретные физически реализуемые примеры, хотя в работе [20] аргументируется возможность существования такого аттрактора в отображении Пуанкаре автономной трехмерной системы обыкновенных дифференциальных уравнений, описывающих модель нейрона. С другой стороны, сконструированный в явном виде пример динамической системы, реализующей надстройку над аттрактором типа Плыкина1), был представлен в диссертации Т. Ханта [21], выполненной в Кембриджском университете под руководством профессора Роберта Мак Кэя. Модель Ханта описывается достаточно громоздкими выражениями, которые различаются в зависимости от рассматриваемой области фазового пространства, она содержит также множество искусственно введенных сглаживающих функций, и не является аналитической. (Построенный поток принадлежит классу С1, то есть решения непрерывны вместе с первыми производными.) Довольно трудно се-
х)Как сейчас известно, можно сконструировать множество различных примеров гиперболических странных аттракторов в ограниченной двумерной области с дырками, отличающихся друг от друга по топологической структуре. Поскольку лежащая в основе идея принадлежит Плыкину, представляется правомерным допустить двоякое употребление термина «аттрактор Плыкина» - в широком смысле, в отношении гиперболических аттракторов указанного типа в двумерном фазовом пространстве, и в узком смысле, когда имеется в виду конкретный пример из его исходной работы [22]. Объект, реализующийся в сечении Пуанкаре модели Ханта, может быть назван аттрактором Плыкина в широком смысле. Далее, для краткости мы будем пользоваться указанным термином без дополнительных оговорок.
177
бе представить, чтобы эту модель можно было реализовать физически, например, в виде радиотехнической схемы.
Тем не менее разработку модели Ханта как системы с непрерывным временем следует приветствовать как определенный шаг в направлении построения примеров реальных систем с гиперболическим странным аттрактором. В данной работе мы намерены воспроизвести конструкцию Ханта и выполнить численное исследование динамики этой системы с привлечением арсенала хорошо отработанных методов нелинейной динамики, включая построение фазовых портретов, анализ реализаций, вычисление показателей Ляпунова и размерностей. Одна из главных целей состоит в том, чтобы отладить методику и накопить опыт обращения с гиперболическими аттракторами, в том числе в отношении применения численной процедуры проверки гиперболичности, основанной на выполнении определенных условий вложения растягивающихся и сжимающихся конусов в пространстве векторов возмущения (касательном пространстве).
1. Качественное описание модели Ханта
Модель Ханта представляет собой неавтономную систему, описываемую дифференциальными уравнениями относительно двух переменных x и у, с правыми частями, зависящими от х, у и времени t:
dx/dt = f*(x,y,t), dy/dt = g* (x,y,t). (1)
Фигурирующие здесь функции f * и g* непрерывные и дифференцируемые и имеют период 2п по аргументу t. Формальное описание модели Ханта и математические соотношения для вычисления функций f* и д* приведены в Приложении A.
На рис. 1, а на плоскости (x, у) показана начальная (при t = 0) конфигурация области, динамика точек которой отслеживается в рамках предлагаемой конструкции. Эта область U содержит три подобласти Ui, в которых определены специальные криволинейные системы координат (r, 0)i, i = 1, 2, 3 (соответствующие семейства координатных линий приведены на рисунке). По ходу эволюции во времени положение точек, представляющих начала отсчета этих специальных координат, определенным образом изменяется, а сами области и координатные линии непрерывно деформируются.
Динамика переменных x, у на интервале времени 2п рассматривается на трех последовательных стадиях, протяженность каждой из которых 2п/3. На первой стадии область U претерпевает вертикальное сжатие и горизонтальное растяжение, как показано на рис. 1, б. (Говоря более аккуратно, сжатие и растяжение осуществляется по координатным линиям введенных криволинейных координат, соответственно, по линиям 0 = const и r = const.) На второй стадии расположенная в правой части фигуры трансформированная подобласть U2 плавно изгибается вверх, влево и вниз (рис. 1, в), так что в итоге оказывается расположенной вдоль границы области U3 (рис. 1, г). На третьей стадии подобласть Ui в левой части фигуры претерпевает аналогичную деформацию, изгибаясь вверх, вправо и вниз (рис. 1, д), так что в итоге ее нижняя граница располагается по краю преобразованной на предыдущей стадии области U2 (рис. 1, е). Фигура воспроизводит по форме исходную область,
178
Рис. 1. Эволюция области U = U1 U U2 U U3 на плоскости x, у за период At = 2п в модели Хан-
та. Жирные точки на диаграмме (а) отвечают точкам начала отсчета специальных криволинейных координат
причем координатные линии 0 = const и r = const в их новом положении направлены вдоль исходных координатных линий криволинейных координат. Описанная деформация составляет, так сказать, основное содержание динамики. Кроме того, в модель заложены дополнительные модификации, благодаря которым точки начала отсчета криволинейных координат становятся отталкивающими, а поперечное сжатие элементов фигуры несколько более сильным, так что в итоге трансформации она оказывается внутри исходной области (при совпадении нижних границ исходной области и ее образа). Степень сглаживания определяется параметром, величина которого в описанных ниже расчетах принята равной е = 0.17. Также в модели предусмотрено определение векторного поля, отвечающего правым частям уравнений (1), за пределами деформирующейся области U(t).
2. Аттрактор модели Ханта и его свойства
Если задать начальные условия xo = (xo,yo) и проинтегрировать уравнения (1) на временном интервале At = 2п, то получим новый вектор состояния xi = (xi, yi), который можно трактовать как результат применения некоторого отображения T к исходному вектору xo. Это отображение Пуанкаре или стробоскопическое отображение для построенного потока. При задании начального вектора в области U образ также попадает в U. Иными словами, образ T(U) принадлежит U, а аттрактор отображения T можно определить как пересечение множеств, получаемых многократным применением отображения: A = n^=iTk(U). Для системы с дискретным временем, заданной отображением T, это гиперболический странный аттрактор Плыкина [1-7, 22] (см. сноску 1). В контексте конкретной анализируемой модели, математическое обоснование свойства гиперболичности аттрактора представлено в
179
работе Ханта. Альтернативный подход к проверке этого свойства на основе компьютерных вычислений рассматривается в Приложении B.
Для численного решения уравнений мы использовали метод Рунге-Кутты 4-го порядка. Поскольку функции, фигурирующие в определении модели (1), не относятся к классу гладких с большим числом производных, метод высокого порядка неизбежно теряет в точности. Тем не менее эмпирически получаемые результаты с применением данного метода лучше, чем при использовании методов более низкого порядка. Если требовалась повышенная точность, это достигалось уменьшением шага интегрирования.
На рис. 2 показан портрет аттрактора в сечении Пуанкаре, полученный для модели Ханта в результате численного решения уравнений (1). Это набор точек представительной траектории на аттракторе в моменты времени, кратные 2п, в достаточно большом количестве (несколько десятков тысяч). Аттрактор характеризуется присутствием поперечной фрактальной структуры. Хорошо видны несколько первых уровней этой структуры - объект как бы составлен из полосок, каждая из которых содержит полоски следующего уровня, и так до бесконечности.
Поскольку основной мотив разработки модели Ханта заключался в построении системы с непрерывным временем, имеющей гиперболический аттрактор, естественно обратиться к иллюстрациям, которые передавали бы именно отличительные особенности потоковой системы.
На рис. 3 показан аттрактор потоковой системы в проекции на плоскость (x, у). Интересно, что, подобно аттрактору Лоренца, он составлен из двух «крыльев бабочки», которые, однако, асимметричны. За формирование этих «крыльев» отвечают вторая и третья стадия построения динамики модели Ханта, когда осуществляется деформация поворота боковых подобластей (см. рис. 1). На рис. 4 показаны зависимости переменных x и у от времени. Отметим нерегулярное поведение этих зависимостей при наблюдении на достаточно больших временах, что свидетельствует о хаотической динамике. Зависимости непрерывные, но на каждом периоде протяженностью 2п характеризуются присутствием стыкованных друг с другом участков разного поведения - горизонтальных плато и резких всплесков. В этом выражается специфика модели Ханта, динамика которой строится в виде последовательных стадий.
Рис. 2. Аттрактор модели Ханта в сечении Пуанкаре
Рис. 3. Портрет аттрактора модели Ханта в проекции на плоскость x, у
180
Рис. 4. Зависимости динамических переменных x и у от времени в модели Ханта
На рис. 5 приведен спектр мощности сигнала для переменной x при динамике на аттракторе в модели Ханта. По оси ординат использован логарифмический масштаб. В соответствии с рекомендациями теории численного спектрального анализа случайных процессов [23], имеющаяся реализация разбивалась на отрезки протяженностью То, на каждом отрезке производилось умножение на функцию sin2(nx/To) (так называемое «окно»), после чего выполнялось преобразование Фурье, и результат усреднялся по всем отрезкам реализации. Величина Дю = п/Т0 определяет разрешающую способность анализа. Чем больше количество имеющихся в наличии отрезков разбиения, тем меньше среднеквадратичная погрешность оценки спектральной плотности мощности.
Реализация была представлена выборкой отсчетов сигнала с шагом ДЬ = 2п/120, каждый отрезок Т0 содержал 213 отсчетов, усреднение проведено по 64 отрезкам. Из графика видно, что спектр сплошной, то есть имеет такой же характер, как спектр стационарного случайного процесса. Присутствующие пики на частотах 1, 2, 3... связаны с присущей модели Ханта периодичностью во времени правых частей определяющей системы дифференциальных уравнений.
Как известно, один из главных атрибутов динамического хаоса - неустойчивость фазовых траекторий или наличие высокой чувствительности движения по отношению к малым вариациям начальных условий. Для демонстрации этого свойства применительно к модели Ханта проводилось многократное численное решение уравнений (1) при задании в качестве начальных условий каждый раз некоторой определенной точки на аттракторе с малой случайной добавкой к координатам x и у. Полученные результаты представлены в виде наложения 20 реализаций на одном графике на рис. 6. Как видно из рисунка, на начальном участке реализации повторяют друг друга с высокой степенью точности, но с течением времени, ближе к право-
Рис. 5. Спектр Фурье сигнала, порождаемого динамикой модели Ханта на аттракторе
181
Рис. 6. Иллюстрация неустойчивости траекторий на аттракторе модели Ханта по отношению к малому возмущению начальных условий. На графике показаны наложенные друг на друга зависимости динамических переменных от времени для 20 реализаций со слегка различающимися начальными условиями
му краю диаграммы, они расходятся друг относительно друга все сильнее, и картина «замазывается». Сохраняющаяся на этом фоне периодическая структура всплесков обусловлена уже упоминавшейся спецификой модели Ханта как неавтономной динамической системы с периодическими во времени правыми частями уравнений.
Для количественной характеристики присущей хаосу неустойчивости траекторий используются показатели Ляпунова. В нашем случае имеется два нетривиальных показателя, один из которых, положительный, отвечает за неустойчивость движения на аттракторе, а другой, отрицательный - за приближение траекторий к аттрактору.
Для расчета показателей Ляпунова был применен алгоритм Бенеттина [24, 12]. Из-за сложности формального определения модели Ханта, выбрана версия алгоритма, не требующая вывода уравнений в вариациях. Производилось численное решение трех комплектов уравнений (1). Один из них отвечал движению изображающей точки по основной траектории, а два других - по близким возмущенным траекториям с фиксированной начальной нормой вектора возмущения л/ (хг — x0)2 + (У г — y0)2 = = £о << 1, где i = 1, 2. В конце каждого временного интервала протяженностью 2п проводилась ортогонализация и перенормировка векторов возмущения, с тем чтобы нормы снова равнялись е0. Показатели Ляпунова определялись как средние (за достаточно большое количество этих интервалов) скорости нарастания для накапливающихся сумм логарифмов отношения норм возмущений в конце и начале интервала. В результате усреднения по 22 реализациям длительностью 1000 периодов 2п получены следующие значения показателей: Xi = 0.1532 ± 0.0002, Х2 = —0.1930 ±0.0020, где в качестве погрешности указана величина среднеквадратичного отклонения. Показатели Ляпунова для отображения T определяются как Лц2 = 2nXi 2, и составляют, соответственно, Л1 = 0.9625, Л2 = —1.213. Заметим, что больший показатель хорошо согласуется с оценкой Л = ln[(3 + \/б)/2] = 0.9624, которая получается при описании динамики на аттракторе с использованием одномерного отображения (см. раздел 3). Сумма показателей Ляпунова отрицательная. Это выражает то обстоятельство, что в процессе динамики с выходом на аттрактор облака начальных точек его объем убывает во времени экспоненциально и стремится в итоге к нулю.
182
Оценка размерности аттрактора в сечении Пуанкаре по формуле Каплана-Йорке дает dL = 1 + Xi/|X2| ~ 1.793. Для аттрактора как объекта в расширенном фазовом пространстве размерность на единицу больше.
3. Символическая динамика и периодические орбиты
Обратимся к подходу, известному как символическая динамика. Для этого область в сечении Пуанкаре, посещаемую траекториями, разобьем определенным образом на подобласти, и каждую траекторию будем кодировать последовательностью символов, приписанных этим подобластям. Доказано, что применительно к гиперболическому аттрактору разбиение можно построить так, чтобы имело место взаимно однозначное соответствие между траекториями на аттракторе и множеством бесконечных символических последовательностей, составленных по некоторым правилам из набора имеющихся символов (так называемое марковское разбиение).
Для аттрактора модели Ханта разбиение области U, введенной в рассмотрение в разделе 1, показано на рис. 7. Напомним, что область U есть объединение трех областей U1, U2, U3. Подобласть U1 разбивается натри элемента p, a, q. Они определены так, что образы точек из p и q при однократном применении отображения Пуанкаре T(x) попадают снова в U1, тогда как образы точек, принадлежащих a -в область U3. Границы определяются координатными линиями криволинейной системы координат в области Ui (см. определение (A.1) в Приложении А), а именно уравнениями 0(1) = ±(X1 + п)(1 — 2/Х). Подобласть U2 разбита на три элемента x, с, у. Образы точек из x и у при действии отображения Пуанкаре T(x) попадают в U2, а точек из элемента с - в область U1. Границы заданы уравнениями
Рис. 7. Марковское разбиение области, в которой расположен аттрактор в сечении Пуанкаре (а). Представление динамики с использованием разветвленного одномерного многообразия (б), где на двух диаграммах, отвечающих исходному объекту и его трансформации, внутренность соответствующих друг другу петель закрашена одинаковым образом. Представление динамики на разветвленном одномерном многообразии с помощью итерационной диаграммы (в). Показан участок траектории с кодом xcqab...
183
0(2) = ±(X2 + п)(1 — 2/А) и отвечают координатным линиям криволинейной системы координат в области U2. Наконец, область U3 состоит из одного элемента b. Последовательности, кодирующие принадлежащие аттрактору траектории, строятся из семи символов (p, a, q, x, c, y, b) с применением правила, сформулированного ниже.
Полезное представление динамики получается, если принять соглашение не различать точки, которые в одной и той же области Ui имеют одинаковую координату 0(i). Модель Ханта устроена так, что в сечении Пуанкаре их образы будут оставаться на общей координатной линии 0 ив процессе последующих итераций. В такой трактовке фазовое пространство выступает как одномерное: его можно мыслить как резиновую нить с тремя петлями, натянутую на три гвоздика, как показано на рис. 7, б (так называемое разветвленное одномерное многообразие). Одна итерация отображения Пуанкаре соответствует тому, что мы определенным образом деформируем эту нить и вновь натягиваем ее на те же гвоздики.
Графически динамику можно иллюстрировать с помощью диаграммы, показанной на рис. 7, в. На каждой из координатных осей выделены три отрезка, в пределах каждого из которых отсчитывается координата 0(i) для i = 1, 2, 3. Итерации представляются в виде традиционной «лесенки Ламерея», причем посещение того или иного куска графика отвечает определенному символу в последовательности, кодирующей данную траекторию. (Отображение выглядит разрывным, но это дефект представления: на разветвленном многообразии оно на самом деле непрерывное.) Заметим, что все ветви графика имеют одинаковый наклон, выражающийся по абсолютной величине константой А = (3 + \/б)/2. Отсюда можно заключить, что отображение Пуанкаре должно иметь показатель Ляпунова, равный ln[(3 + л/5)/2] = 0.9624.2)
Из рассмотрения диаграммы на рис. 7, в видно, что «грамматика» символического описания траекторий выражается следующими правилами:
Р ^ P,a,q; a ^ b; q ^ p,a,q; x ^ x,c,y;
b 2)
c ^ p,a,q; y ^ x,c,y; b ^ x,c,y,
где для каждого символа указаны символы, которые после него разрешается ставить.
Периодические символические последовательности ассоциируются с периодическими орбитами или циклами на аттракторе модели Ханта. Используя правила (2), можно перечислить циклы всевозможных периодов. В частности, имеется три цикла периода 1 (коды p, q, y), два цикла периода 2 (коды pq, xy), пять циклов периода 3 (коды ppq, pqq, abc, xxy, xyy). Имеется десять циклов периода 4 (коды pppq, ppqq, pqqq, pabc, qabc, xxxy, xxyy, xyyy, xcab, ycab), двадцать четыре цикла периода 5, пятьдесят циклов периода 6 и т.д.
Поскольку аттрактор странный гиперболический, все эти циклы неустойчивы. Чтобы их найти в численных расчетах, приходится использовать специальные процедуры поиска. Один возможный подход состоит в том, чтобы итерировать представленное на рис. 7, в одномерное отображение в обратном времени, выбирая на каждом шаге нужную ветвь в соответствии с заданной заранее периодической символической последовательностью. Полученное в результате значение 0 дополняется
2) Разумеется, это только один, наибольший по величине из двух показателей. Второй при переходе к описанию динамики на одномерном разветвленном многообразии из рассмотрения исключается.
184
произвольно выбранным начальным значением г, и производится решение уравнений (1) в прямом времени за период цикла. Далее найденная точка на цикле уточняется, например, по методу Ньютона.
На рис. 8 показаны фазовые портреты некоторых циклов модели Ханта с периодами от 2 до 8 в проекции на плоскость (x, y). Для каждого цикла указана соответствующая ему символическая последовательность. Ее можно восстановить, внимательно рассматривая рисунок: жирные точки на фоне изображенного серым фоном
Рис. 8. Фазовые портреты неустойчивых циклов для модели Ханта в проекции на плоскость (x, y). Стрелкой показано направление движения по траекториям. Жирные точки на траектории отвечают моментам времени t = 2пи, где и - целое, то есть сечению Пуанкаре. Серым изображен аттрактор в сечении Пуанкаре. Для каждого цикла указан соответствующий ему символический код, количество символов которого есть период цикла в единицах 2п
185
Таблица
Символические коды, мультипликаторы и показатели Ляпунова для некоторых периодических орбит модели Ханта с периодами от 2 до 9
Период Код Hi И-2 Л1 Л2
2 pq -6.857 -0.1456 0.9626 -0.9634
3 pqq 17.928 0.04953 0.9621 -1.0017
4 pabc 46.977 0.004880 0.9624 -1.3307
4 ppqq 47.061 0.004976 0.9628 -1.3258
5 pqqpp 123.17 0.0004596 0.9627 -1.5370
5 pabcp 123.11 0.002027 0.9626 -1.2403
5 pabxc 122.89 0.003478 0.9623 -1.1323
6 ppqqqq 322.27 0.0008787 0.9626 -1.1728
6 pabcqq 321.90 0.0005048 0.9624 -1.2652
6 pabxcp 322.09 0.001021 0.9625 -1.1478
7 ppqqabc 844.70 0.0001869 0.9627 -1.2264
7 ppabxcp 845.55 0.0001006 0.9629 -1.3149
8 ppqqabcp 2203.1 0.00001750 0.9622 -1.3695
8 pabxcpqq 2201.8 0.0001814 0.9621 -1.0768
9 ppqqabyyc 5777.9 0.00002745 0.9624 -1.1670
9 pabcpqppq 5763.7 0.00002374 0.9621 -1.1831
9 ppabxyycp 5802.6 0.00001318 0.9628 -1.2486
аттрактора в сечении Пуанкаре располагаются как раз в областях, определенных как на рис. 7, а.
Для этих же и некоторых других циклов в таблице приводятся рассчитанные численно мультипликаторы - собственные числа матриц, описывающих преобразование вектора малого возмущения за период цикла в линейном приближении. Видно, что для каждого цикла из двух мультипликаторов один по модулю всегда больше единицы, а другой меньше. Так что эти циклы, в самом деле, неустойчивые, притом седловые, как это и должно быть для всех траекторий, принадлежащих гиперболическому странному аттрактору. В таблице приводятся также показатели Ляпунова, вычисленные для этих циклов. Больший показатель практически одинаков для всех циклов, и с точностью до погрешности вычислений совпадает с величиной Л = ln[(3 + \/б)/2] = 0.9624. Второй показатель для разных циклов варьируется в некоторых пределах, но всегда отрицателен, притом по абсолютной величине превышает Л.
Заключение
В данной работе мы воспроизвели построение искусственной модельной системы Ханта с непрерывным временем, обладающей странным аттрактором гиперболического типа. Было предпринято достаточно подробное численное исследование динамики этой системы, включая построение фазовых портретов, анализ реализаций, вычисление показателей Ляпунова и размерностей, спектра Фурье. Представлены результаты проверки в численных расчетах условий теоремы, гарантирующей гиперболическую природу аттрактора, что позволило устранить допущенную в ра-
186
боте Ханта неточность. Полученные результаты могут использоваться в учебном процессе в качестве примера динамики на гиперболическом странном аттракторе, и как иллюстрация соответствующей общей математической теории. Выполнение данного исследования способствовало накоплению опыта работы с гиперболическими аттракторами, что становится актуальным по мере появления примеров физически реализуемых систем, имеющих аттракторы подобного типа. Одним из интересных вопросов для дальнейших исследований остается поиск или конструирование реальных систем с аттракторами типа Плыкина, поскольку модель Ханта сама по себе вряд ли допускает реализацию в силу сложности и громоздкости соотношений, положенных в ее основу.
Работа выполнена при поддержке РФФИ (грант № 06-02-16619).
Приложение A
Формальное описание динамики модели Ханта
Динамика переменных х, у на временном интервале At = 2п рассматривается на трех последовательных стадиях, протяженность каждой из которых 2п/3. При формулировке соотношений используется фиктивное время s, которое в пределах каждой стадии меняется от 0 до 1. С целью обеспечить гладкость потока, связь между s и t выбирается так, чтобы на стыках стадий скорость движения изображающей точки на плоскости х, у обращалась в ноль.
При описании потока используются три специальных криволинейных системы координат (r, 0)\ i = 1, 2, 3, связанных с прямоугольными координатами х, у соотношениями
х = X + Si(n/2 — |0|), у = Y + rsgn0,
|0| > п/2,
х = X + Sir cos fr (0),
у = Y + r sin fr(0), |0| < n/2,
(A.1)
где fr(0) = 2n2r0/[n?(r+1)+4(r —1)02], Si = —1 и S2>3 = 1. Величины X и Y определяют положение центра i-й системы в исходных прямоугольных координатах. Для краткости переход к новым координатам и обратно обозначим, соответственно, (r, 0)i = (Е(х,у, X,Y, Si), 0(х,у, X,Y, Si)) и
(x, у) = (Fi(r, 0, X, Y, Si), Gi(r, 0, X, Y, Si)). Явные выражения для функций R, 0, Fi, R1 легко выводятся из формул (A.1).
Положим X = (3 + \/б)/2 и ^ = (3 — л/5)/2. Введем, кроме того, константы, определяющие геометрическое расположение элементов конструкции на плоскости (х,у):
187
Xi = 1 (6^5 - 6 - (3 -V5)k) , X2 = 3, X3 = 1 (18 - 6^5 - (V5 - 1)п) ,
Yi = 1(1 + V5), Y2 = 1, Y3 = 1(3 + V5),
Ri = 2(1 + V5), R2 = 1, R3 = 1(V5 - 1).
(A.2)
1. На первой стадии t e [0, 2п/3] полагаем s = sin2((3/4)t). Точки начала отсчета специальных координат перемещаются в зависимости от s в соответствии с уравнениями
Xji(s) = Х(Хг + п/2) - п/2, Yii = vsY, i = 1, 2, 3. (A.3)
Кроме того, полагаем Rii(s) = psRi. Задавшись определенными значениями s, x, y, находим векторное поле f(x,y,s) = (f (x,y,s), g(x,y,s)) посредством следующей процедуры.
а. Если x < 0, выполняем переход в систему координат с индексом 1
(r, 0)i = (R(x,y, -Xii(s),Yii(s), -1), &(x, y, -Xii(s),Yn(s), -1)) (A.4)
и, обозначая точкой производную по s, полагаем
0 = (ln 1)0, Г = (у(0, Xii(s))^i2(r, Rii(s)) + (1 - y(0, Xii(s)))) (ln p)r. (A.5)
Здесь введены функции
1,
y(0,x ) =
cos“ 24я (2|0| - п) X i^
0,
cos(n(1 - r/e)), 1,
hi2(r, R) = <
5 1 1
—I— cos —n(r - R) 4 4 V 2e v !
3
1 2,
1 I П
|e| < 2.
п п
- < 0 < - + X,
2 - ~2
|0|> 2+ x,
r < e,
e < r < R - e, ^ - R| < e, r > R + e.
(A.6)
Далее находим компоненты векторного поля (fi, gi) в прямоугольных координатах согласно формуле
188
X — Ssgn(0)0, |0| > п/2,
X + Sr cos fr(0) — Sr sin fr(0) (fdfr(0)/dr + 0dfr(0)/d0^ , |0| < n/2,
Y + sgn(0)f, |0| > n/2,
I Y + f sin fr(0) + rcos fr(0) (fdfr(0)/dr + 0dfr(0)/30^j , |0| < n/2,
. (A-7)
где следует положить i = 1, S = Si = —1, X = —(lnХ)(Хц + n/2), Y = (ln ^)Yn.
б. Если x > 0, то аналогичные вычисления выполняем в системах координат 2 и 3. Находим сначала
(r, 0)2 = (R(x,y,X2i(s),Y2i(s), 1), 0(x,y,X2i(s), Y2i(s), 1)). (A.8)
Полагаем 0 = (ln X)0,
r = (y(0, X2i(s) — X3i(s))^i2(r, R2i(s)) + (1 — y(0, X2i(s) — X3i(s)))) (ln ц)г
(A.9)
и по формуле (A.7), в которой нужно принять i = 2, S = S2 = 1,
X = (lnX)(X2i + n/2), Y = (lnp,)Y2i, находим компоненты векторного поля (f2,g2). Далее вычисляем
(r, 0)3 = (R(x,y,X3i(s),Y3i(s), 1), 0(x,y,X3i(s),Y3i(s), 1)) (A.10)
и полагаем
0 = (lnX)0, r = (y(0, X3i(s))^3(r,R3i(s)) + (1 — y(0, X3i(s)))) (ln ц)г, (A.11)
где
{cos n(1 — r/e), r < e,
(A.12)
1, r > e.
По формуле (A.7), в которой нужно взять i = 3, S = S3 = 1,
X = (lnX)(X3i + n/2), Y = (ln ^)Y3i, находим компоненты векторного поля (f3,g3).
в. Определим теперь векторное поле
{ (fi,gi), x < 0,
f (x,y,s) = < (A.13)
I w(d3,d2) ■ (f2,g2)+ w(d2,d3) ■ (f3,g3), x> 0,
где введена функция
w = w(u, v)
sin2 f 2nu(u + v) ij
(A.14)
а ее аргументы определены как da = max {Ra(x,y, Xai(s),Yai(s)) — цsRa, 0}, a = 2, 3.
189
г. В качестве заключительного шага вычислений на стадии 1 положим
3
f(x,y,s) = f(x,y,s) + ^ P(||x - Xn\\) • (Xji(s) -f(x,y,s) + (lnX)(x - Xii^.
i=i
. . . (A.15)
Здесь x = (x,y), Xii(s) = (SiXii, Уц), Xii(s) = (SiXii, Yi), и введена
функция
1, P < е/4,
в(р) = s cos2(n(2p/e — 1/2)), е/4 < p < е/2,
, 0, p > e/2.
(A.16)
Примечание. В этом месте в работе Ханта имеется неточность. В то время как определения в тексте его диссертации соответствуют формулам (A.15) и (A.16), в коде программы в среде «Mathematical в качестве аргумента функции в фигурирует квадрат нормы ||x — Xii ||2. Наши расчеты с использованием процедуры, описанной в Приложении B, показали, что это момент принципиальный: при подстановке квадрата нормы гиперболичность аттрактора нарушается! С другой стороны, если следовать определениям (A.15) и (A.16), то аттрактор действительно гиперболический. Однако при выбранной у Ханта слишком малой величине е = 0.05 фрактальная поперечная структура аттрактора в исправленной версии модели оказывается плохо различимой на иллюстрациях. Поэтому при вычислениях, представленных в настоящей статье, мы берем увеличенное значение е = 0.17, оставаясь в то же время внутри указанного в работе Ханта допустимого интервала.
2. На второй стадии t е [2п/3, 4п/3] полагаем s = sin2((3/4)t — (1/2)п). Берем за начало отсчета специальных координат точку
X22(s) = X2 + (1 - s)D,
Y22(s) = Y2 + (1 - s)D, (A.17)
D = (X - 1)(X2 + n/2) + ^^2.
Задавшись мгновенными значениями (x, y), перейдем к новым координатам
(r, 0) = (R(x,y,X22(s),Y22(s), 1), ®(x, y, X22(s), У22(s), 1)) (A.18)
и примем, что
r = -D, 0 = D. (A.19)
Производим обратный переход к прямоугольным координатам согласно формуле (A.7), где полагаем S = 1, X = -D, У = -D, и находим компоненты векторного поля (fo,go). Окончательно принимаем
f(x,y,s) = w(a2,b2) • (fo,go), (A.20)
190
где функция w определена формулой (A.13), а ее аргументы выражаются формулами
(A.21)
й2 = max{R(x,y,X22(s), 2(Y2 + Y22(s)), 1) - ^(-Yi + Y22(s)) - mR, 0}
b2 = max { Y22(s) - 2R2M- - R(x,y, X22(s),Y22(s), 1), 0} .
,3
3. На третьей стадии t £ [4п/3, 2п] полагаем s = sin2 (4t - л). Берем за начало отсчета специальных координат точку
Xi3(s) = Xi + (1 - s)D,
Yi3(s) = Yi + (1 - s)D, (A.22)
D = (X - 1)(X1 + п/2) + y,R1.
Задавшись исходными значениями (x, y), выполним переход к координатам (r, 0)
(r, 0) = (R(x, y, -Xi3(s),Yi3(s), -1), &(x,y, -Xi3(s),Yi3(s), -1)). (A.23)
Далее, определяем поток формулой (A.19) и производим переход к прямоугольным координатам согласно (A.7), считая, что S = -1, X = -D, Y = -D. Полагаем
f(x,y,s) = w(a3,b3) • (fo,go), (A.24)
где функция w дается выражением (A.14), и
аз = max{R(x, y, -XVi(s), 2(Yi + YVi(s)), -1) - 2(-Yi + Yi3(s)) - y,R2,0},
Ьэ = max{Yi3(s) - 2Ri^ - R(x, y, -Xi3(s),Yi3(s), -1), 0}.
4. В итоге, имея в виду, что f(x, y, s) = (f (x, y, s), g(x, y, s)), получаем уравнения, описывающие динамику модели с непрерывным временем на всех трех стадиях, в виде
dx/dt
33
41 sin2t\f (x,У,s(t)),
dy/dt
33
4\ sin2t\g(x,y,s(t)),
(A.25)
где в правых частях учтен множитель s'(t) = (3/4) \ sin(3/2)t\, возникающий при переходе от фиктивного времени s к естественному времени t.
Приложение B
Достаточные условия гиперболичности и методика их проверки
В основу проверки гиперболичности аттрактора положим известную из математической литературы теорему [1-7, 21]. Нам достаточно сформулировать ее применительно к диффеоморфизму класса с1 в двумерном евклидовом пространстве R2
191
{xi,x2}, в качестве какового будет выступать отображение Пуанкаре модели Ханта T(x). При введении малого возмущения исходного состояния x ^ x + 8x преобразование вектора возмущения при действии отображения x' = T(x) дается в линейном приближении соотношением 8x' = DTx8x, где DTx - матрица Якоби, вычисленная в точке x. Обозначение DT- отвечает матрице производных обратного отображения T-l(x).
Теорема [1-7, 21]. Пусть A есть подмножество пространства М2, инвариантное по отношению к действию диффеоморфизма T, отображающего ограниченную область U С М2 в себя. Множество A с Int U будет однородным гиперболическим, если существует некоторая константа у > 1 и выполнены следующие условия:
1) Для каждой точки x е A в пространстве Vx 2-мерных векторов bx определены расширяющийся и сжимающийся конусы, соответственно, SX и CX такие, что ||DTXu|| > у ||u|| для всех u е SX и ||DT-lv|| > у ||v|| для всех v е CX, причем для всех x е A выполнено SX П CX = 0 и SX + CX = Vx.
2. Конусы SX инвариантны по отношению к действию DT, а C'X по отношению к действию обратного отображения DT_l, то есть для всех x е A DTx(SX) С ST(X) и DT-l(CX) С CT_1(X).
Определение расширяющихся и сжимающихся конусов, вообще говоря, зависит от системы координат. Их невыполнение при попытке проверить условия теоремы может быть связано не с нарушением гиперболичности, а с неудачным выбором координат. Для модели Ханта неподходящими в этом смысле оказываются прямоугольные координаты (x,y), тогда как к успеху приводит использование криволинейных координат (0, т)г, введенных в Приложении А, см. формулы (A.1). Заметим, что эти координаты заданы по-своему в каждой из подобластей иг, i = 1, 2, 3, составляющих вместе область U.
Процедура проверки условий вложения для расширяющихся конусов состоит в следующем. Сначала, задавшись начальной точкой x = xo = (x0,y0) E иг, определим ее координаты (0о,то)г. Кроме того, возьмем точки xl = (0o + Н,т0)г и x2 = (0о, т0 + Ь)г, где h << 1, и вычислим для них согласно (А.1) прямоугольные координаты (xl,yl) и (x2,y2).
Выполним на интервале времени t от 0 до 2п численное решение уравнений (4) с начальными условиями (xk ,yk), k=0,1,2. Результирующие векторы состояния обозначим x' = x'o, x'l, x'2, и пусть x'o E Uj. Перейдя к криволинейным ко-
ординатам, найдем (0k,Tk)j, k=0,1,2, вычислим нормированные на h компоненты возмущения и составим матрицу
^ Ull Ul2 \ ( (0l - 0O)h-1
U ' 1 (Tl — т'0 )h~l
V
u2l u22
)
V
(02 0O)h l (т2 — т0 )h_l
)
DTX
(B.1)
Аналогичным образом при старте из точки x' = x'0, получаем матрицу U' DTX/.
Условие, что вектор u = (%, п) принадлежит конусу T(SX), можно представить в виде неравенства | u || > у ||U-lu| или
?2 + п2 > y2[(Ull? + Ul2n)2 + (U2l4 + U22n)2], (B.2)
192
U11 U\2
где Uj - элементы матрицы U"
U21 U22
D = u11u22 — u12u21. Это можно переписать в виде
a-2 + 2b-n + СП2 < 0,
U22D 1 —U12D 1
—u21 D-1 u11D-1
(B.3)
где a = U21 + U21 — y-2, b = U11U12 + u21u22, c = U22 + u22 — у-2. Условие, что тот
же вектор принадлежит конусу У^(x), представляется неравенством ||U'u|| > у ||и||. Его можно переписать в виде
(и11- + и1_2Л)2 + (и21 - + и>22П) > Г(1 + Л2)
(B.4)
или
a'-2 + 2b%n + с'ц2 > 0,
(B.5)
где a' = иП + U221 — У2. b' = U11U12 + U21U22’ c = U122 + U22 — У2-
Как устанавливается в ходе вычислений, для всех точек на аттракторе при надлежащем выборе константы у (из некоторого ограниченного интервала (у1, У2), где у1 < 1 и у2 > 1), справедливы неравенства a < 0, с > 0 и а' > 0, с' < 0, а также b2 — ac > 0 и b'2 — a'с' > 0. При этом соотношение (B.3) выполняется, если
k(1) < пр- < k(2),k(1,2) = —b/c ± д/(b/c)2 — a/c,
что определяет на плоскости (-, п) конус DT(SX)- С другой стороны, неравенство (B.5) имеет место, если
l(1) < п/- < l(2\l(1’2) = —b'/c' ± д/(b'/c')2 — a'/c',
что определяет конус У^(x)- Выполнение вложения DTx(5X) С У^(x) гарантировано, если l(1) < k(1) и l(2) > k(2). Достаточным условием этого служит неравенство, проверка которого и проводилась в ходе вычислений:
H>Q, (B.6)
где
H = д/(b'/c')2 — a'/c' — ^(b/c)2 — a/c, Q = \b'/c' — b/c\. (B.7)
Можно показать, что при выборе у < 1 та же процедура, выполняемая для точек на аттракторе, обеспечивает проверку условия для сжимающихся конусов с параметром у' = 1/у > 1: DT^C^7) С C1-1(x). Это связано с тем, что конусы УY и C 1/y являются друг для друга дополнениями, то есть УY U C 1/y = V.3)
Описанные расчеты были выполнены для множества из 20000 точек в сечении Пуанкаре аттрактора модели Ханта при е = 0.17.
3)При одном и том же у конусы Sy и Су имеют общую границу только при у = 1, тогда как при Y > 1 они не пересекаются: SY ^ СУХ = 0, как и требуется условием теоремы.
193
Для у2 = 2 обнаружено, что max a & —0.304 < 0, min c & 3.96 > 0, min(62 — ac) & 1.898 > 0, min a1 & 3.337 > 0, max c' & —0.311 > 0, min(&'2 — a'c') & 1.827 > 0, 1.1 < H < 4.4, max Q/H & 0.359 < 1. Для у2 = 1/2
найдено, что max a & —0.294 < 0, min c & 3.85 > 0, min(&2 — ac) & 1.932 > 0, min a' & 3.346 > 0, max c' & —0.322 > 0, min(&'2 — a'c') & 1.731 > 0, 1.1 < H < 4.4, max Q/H & 0.300 < 1.
Таким образом, проведенные вычисления свидетельствуют, что для значений у, равных \[2 и 1/\/2, предполагаемые неравенства действительно выполнены на аттракторе отображения Пуанкаре. Поэтому справедливы и соответствующие условия вложения для расширяющихся и сжимающихся конусов, откуда вытекает заключение о гиперболической природе аттрактора.
Библиографический список
1. Синай Я.Г. Стохастичность динамических систем // В кн. Нелинейные волны / Ред. А.В. Гапонов-Грехов. М.: Наука, 1979. С. 192.
2. Современные проблемы математики. Фундаментальные направления // Итоги науки и техники. Т. 2 / Под ред. Р.В. Гамкрелидзе. М.: Изд. ВИНИТИ АН СССР, 1985.
3. Eckmann J.-P. and Ruelle D. Ergodic theory of chaos and strange attractors // Rev. Mod. Phys. 1985. Vol. 57. P. 617.
4. Devaney R.L. An Introduction to Chaotic Dynamical Systems. Addison-Wesley, New York, 1989.
5. Shilnikov L. Mathematical Problems of Nonlinear Dynamics: A Tutorial // Int. J. of Bifurcation and Chaos. 1997. Vol. 7, № 9. P. 1353.
6. Каток А.Б., Хасселблат Б. Введение в современную теорию динамических систем. Пер. с англ. М.: Изд. «Факториал», 1999, 768 c.
7. Afraimovich V and Hsu S.-B. Lectures on chaotic dynamical systems, AMS/IP Studies in Advanced Mathematics, 28, American Mathematical Society, Providence, RI; International Press, Somerville, MA, 2003.
8. Рабинович М.И., Трубецков Д.И. Введение в теорию колебаний и волн. М.: Наука, 1984, 432 с.
9. Анищенко В.С. и др. Нелинейные эффекты в хаотических и стохастических системах. Институт Компьютерных исследований. Москва - Ижевск, 2003.
10. Андронов А.А., Витт А.А., Хайкин С.Э. Теория колебаний. М.: Наука, 1981, 568 с.
11. Неймарк Ю.И. Метод точечных отображений в теории нелинейных колебаний. М.: Наука, 1972, 472 с.
12. Кузнецов С.П.Динамический хаос. 2-е изд. М.: Физматлит, 2006.
13. Kuznetsov S.P. Example of a Physical System with a Hyperbolic Attractor of the Smale-Williams Type // Phys. Rev. Lett. 2005. Vol. 95. P. 144101.
14. Кузнецов C.П., Селезнев Е.П. Хаотическая динамика в физической системе со странным аттрактором типа Смейла-Вильямса // ЖЭТФ. 2006. Vol. 129, № 2. P. 400.
194
15. Кузнецов С.П., Сатаев И.Р. Проверка условий гиперболичности хаотического аттрактора в системе связанных неавтономных осцилляторов ван дер Поля // Изв. вузов. Прикладная нелинейная динамика. 2006. Т. 14, № 5. C. 3.
16. Kuznetsov S.P. and Sataev I.R. Hyperbolic attractor in a system of coupled nonautonomous van der Pol oscillators: Numerical test for expanding and contracting cones // Physics Letters A. 2007. Vol. 365, № 1-2. P. 97.
17. Isaeva O.B., Jalnine A.Yu. and Kuznetsov S.P. Arnold’s cat map dynamics in a system of coupled nonautonomous van der Pol oscillators // Phys. Rev. E. 2006. Vol. 74. P. 046207.
18. Kuznetsov S.P. and Pikovsky A. Autonomous coupled oscillators with hyperbolic strange attractors // Physica D. 2007. Vol. 232. P. 87. Grants of DFG and RFBR 04-02-04011 and 06-02-16619.
19. Кузнецов А.П., Кузнецов С.П., Пиковский А.С., Тюрюкина Л.В. Хаотическая динамика в системах связанных неавтономных осцилляторов с резонансным и нерезонансным механизмом передачи возбуждения // Известия вузов. Прикладная нелинейная динамика. 2007. Т. 15, № 6. P. 75.
20. Belykh V., Belykh I. and Mosekilde E. The hyperbolic Plykin attractor can exist in neuron models // Int. J. of Bifurcation and Chaos. 2005. Vol. 15, № 11. P. 3567.
21. Hunt TJ. Low Dimensional Dynamics: Bifurcations of Cantori and Realisations of Uniform Hyperbolicity. PhD Thesis, Univ. of Cambridge (2000).
22. Плыкин Р.В. Источники и стоки А-диффеоморфизмов поверхностей // Матем. сб. 1974. Т. 94(136), № 2(6). C. 243.
23. Свешников А.А. Прикладные методы теории случайных функций. М.: Изд-во «Наука», Глав. ред. физ.-мат. лит., 1968.
24. Benettin G., GalganiL., GiorgilliA., Strelcyn J.-M. Lyapunov characteristic exponents for smooth dynamical systems and for Hamiltonian systems: A method for computing all of them. Part I: Theory. Part II: Numerical application. Meccanica, 15,1980,9-30.
Саратовский филиал Института Поступила в редакцию 19.05.2008
радиотехники и электроники РАН
Саратовский государственный
университет
CHAOTIC DYNAMICS OF HUNT MODEL - ARTIFICIALLY CONSTRUCTED FLOW SYSTEM WITH A HYPERBOLIC ATTRACTOR
Yu.S. Aidarova, S.P. Kuznetsov
We study numerically chaotic behavior associated with the presence of a hyperbolic strange attractor of Plykin type in the model of Hunt, that is an artificially constructed dynamical system with continuous time. There are presented portraits of the attractor, plots of realizations for chaotic signal generated by the system, illustrations of the sensitive dependence on initial conditions for the trajectories on the attractor. Quantitative characteristics of the attractor are estimated, including the Lyapunov exponents and the attractor
195
dimension. We discuss the symbolic dynamics on the attractor, find out and analyze some unstable periodic orbits belonging to the attractor.
Айдарова Юлия Сериковна - родилась в Саратове (1983), окончила Саратовский государственный университет имени Н.Г. Чернышевского (2006). В данное время является аспирантом кафедры динамических систем на факультете нелинейных процессов. Является соавтором нескольких научных статей.
Кузнецов Сергей Петрович - родился в Москве (1951). Окончил Саратовский государственный университет (1973). С 1988 - сотрудник Саратовского филиала Института радиотехники и электроники РАН, в настоящее время -заведующий лабораторией теоретической нелинейной динамики. Защитил диссертацию на соискание ученой степени кандидата физико-математических наук в СГУ (1977) и доктора наук (1987) по специальности радиофизика. Профессор кафедры динамических систем СГУ. Автор учебно-научной монографии «Динамический хаос» и «Нелинейные колебания» (в соавторстве с А.П. Кузнецовым и Н.М. Рыскиным). Опубликовал свыше 150 научных статей по нелинейной динамике, радиофизике и электронике. Под его руководством защищены семь кандидатских диссертаций. Лауреат государственной научной стипендии для ученых России (1994-1996), Соросовский доцент (1998), Соросовский профессор (2000, 2001). Член редакционной коллегии журнала «Известия вузов. Прикладная нелинейная динамика».
E-mail: [email protected]
196
