Движение падающей пластины в жидкости: конечномерные модели и феномены сложной нелинейной динамики Текст научной статьи по специальности «Физика»

Нелинейная динамика. 2015. Т. 11. № 1. С. 3-49. Полнотекстовая версия в свободном доступе http://nd.ics.org.ru

УДК: 517.9, 532.3, 534.14, 519.622 М8С 2010: 34С15, 76Б99, 37Е99

Движение падающей пластины в жидкости: конечномерные модели и феномены сложной нелинейной динамики

С. П. Кузнецов

Представлен обзор результатов исследования плоской задачи о падении пластинки в сопротивляющейся среде на основе моделей в виде обыкновенных дифференциальных уравнений относительно небольшого числа переменных. Введена в рассмотрение обобщенная модель, в рамках которой с использованием одной и той же системы безразмерных переменных и параметров удается провести сравнительный анализ динамического поведения для моделей Козлова, Танабе-Канеко, Бельмонте-Айзенберга-Мозеса и Андерсена-Пе-савенто-Ванга. Показано, что общая структура устройства пространства параметров для разных моделей имеет определенное сходство, обусловленное, очевидно, одинаковой присущей симметрией и общей природой вовлеченных феноменов нелинейной динамики (неподвижные точки, предельные циклы, аттракторы, бифуркации). Для задачи о движении тела эллиптического профиля в вязкой среде в присутствии циркуляции вектора скорости и приложенного постоянного вращающего момента обнаружено присутствие странного аттрактора Лоренца в трехмерном пространстве обобщенных скоростей.

Ключевые слова: движение тела в жидкости, автоколебания, авторотация, флаттер, аттрактор, бифуркация, хаос, показатель Ляпунова

Получено 22 декабря 2014 года После доработки 16 января 2015 года

Работа выполнена при частичной поддержке гранта Президента РФ для ведущих научных школ НШ-1726.2014.2 «Фундаментальные проблемы нелинейной динамики и их приложения».

Кузнецов Сергей Петрович [email protected]

Институт радиотехники и электроники им. В. А. Котельникова РАН, Саратовский филиал 410019, г. Саратов, ул. Зеленая, д. 38

Содержание

Введение ....................................

1. Плоская задача о падении тела в идеальной жидкости . . . .

2. Консервативная динамика без гравитации............

3. Модель Козлова: падение тела при наличии вязкого трения

4. Модель с циркуляцией и аттрактор Лоренца..........

5. Модель Танабе — Канеко........................

6. Модель Бельмонте — Айзенберга — Мозеса ............

7. Модель Андерсена — Песавенто — Ванга..............

8. Обобщенная модель..........................

31 31 33 36 39

25

24

18

15

11

8

6

4

8.1. Уравнения и нормировка.....................

8.2. Мультистабильность и хаос в модели Козлова.........

8.3. Модифицированная модель Танабе-Канеко..........

8.4. Модифицированная модель Бельмонте - Айзенберга - Мозеса

8.5. Модель Андерсена-Песавенто-Ванга в сравнении с другими моделями . . 39

Заключение

45

Введение

Падение плоского листа в сопротивляющейся среде составляет содержание одной из классических задач гидро- и аэродинамики, анализ которой восходит к работам Максвелла, Кельвина, Кирхгофа, Жуковского, относящимся еще к XIX и началу XX века [1-6].

Элементарные эксперименты показывают, что возможны разнообразные режимы движения в зависимости от параметров задачи и начальных условий:

• простое равномерное падение,

• падение с колебаниями из стороны в сторону, которые могут быть как регулярными, так и нерегулярными (флаттер),

• падение с кувырканиями, периодическими или нерегулярными (авторотация).

Естественно, возникает желание выяснить условия реализации этих режимов, разобраться в их природе в контексте теории динамических систем, исследовать бифуркации, приводящие к возникновению тех или иных типов движения, и т. д.

Полное и корректное описание задачи о движении тела в вязкой несжимаемой жидкости подразумевает нахождение переменного во времени поля скоростей в окружающей области на основании уравнений Навье-Стокса [7-12], что требует непростых ресурсоемких компьютерных расчетов. Получаемые на этом пути результаты заведомо труднообозримы, и, с учетом ожидаемого многообразия режимов, зависящих от большого числа параметров, все равно требуют для своего понимания качественной интерпретации на физическом уровне.

Первый разумный шаг к упрощению анализа состоит в том, чтобы ограничиться плоской задачей и считать, что существенны две пространственные координаты X и У, тогда

как по третьей оси Z система протяженная, никаких существенных зависимостей от этой координаты нет и движений вдоль Z не происходит.

Далее, можно обратиться к приближенному описанию с привлечением моделей в виде обыкновенных дифференциальных уравнений для небольшого числа переменных, то есть к динамическим системам с относительно небольшой размерностью фазового пространства. Хотя правомерность такого подхода не очевидна, подчеркнем тот чрезвычайно ценный момент, что в его рамках для анализа задачи можно привлечь мощный и обширный концептуальный инструментарий современной теории динамических систем. В частности, при численном исследовании моделей будут полезны такие хорошо отработанные и апробированные методики, как построение портретов аттракторов [13-15], карт динамических режимов в пространстве параметров [15-17], вычисление показателей Ляпунова [14, 15, 18], бифуркационный анализ [19].

Принципиальной основой для обращения к конечномерному описанию может служить то обстоятельство, что в случае идеальной невязкой несжимаемой жидкости система уравнений для обобщенных координат и скоростей собственно твердого тела отделяется от уравнений динамики жидкости. Соответствующие уравнения в свое время выведены и изучены Кирхгофом [2, 7-10, 20, 21]. Влияние жидкости на динамику тела выражается в том, что модифицируются инерционные свойства системы — помимо массы и момента инерции тела в уравнениях присутствуют присоединенные массы и моменты инерции, обусловленные вовлечением в движение прилегающих областей жидкости. Действие силы на тело со стороны среды определяется величиной циркуляции вектора скорости вокруг профиля тела, причем в случае идеальной жидкости циркуляция выступает просто как постоянный параметр (в частности, она может быть нулевой).

Собственно уравнения Кирхгофа относятся к ситуации, когда потери механической энергии исключены и соответствующая динамическая система является консервативной. В этих предположениях, за рамками рассмотрения остаются многие принципиальные с практической точки зрения свойства динамического поведения при падении тел, обусловленные присутствием диссипации, в том числе возникновение установившихся режимов, таких как стационарное падение, регулярные и хаотические самоподдерживающиеся колебания и вращения [20, 21]. Для их изучения кажется естественным обратиться к феноменологическому учету диссипации, добавив в уравнения Кирхгофа выбранные подходящим образом дополнительные члены [22-28].

Настоящий обзор посвящен моделям на базе обыкновенных дифференциальных уравнений, которые формулируются в рамках плоской гидродинамической задачи, и сравнению этих моделей. Одна из главных целей в том, чтобы наполнить выявляемую картину феноменов нелинейной динамики конкретным содержанием и иллюстративным материалом, полученным в численных расчетах. Вопросы сравнения с результатами, базирующимися на уравнениях Навье-Стокса [29-32], а также с экспериментальными данными [33-38], затронуты не будут. Также не будут обсуждаться ситуации движения тел в сопротивляющейся среде, выходящие за рамки плоской задачи [36-38], и обобщения, в том числе касающиеся вопросов управления движением тела в жидкости [39-43].

В разделе 1 рассмотрена плоская задача о падении тела с профилем эллипса в идеальной жидкости и приведены уравнения движения с учетом эффекта присоединенных масс. В разделе 2 обсуждается частный случай, когда гравитация компенсирована выталкивающей силой, и описание динамики редуцируется к уравнению типа математического маятника с нелинейностью синуса. В разделе 3 рассматривается модель Козлова, которая получается из консервативной системы уравнений учетом только вязких (линейных

по обобщенным скоростям) сил и момента силы сопротивления. Рассмотрен частный случай, когда модель интегрируема, сводящийся к уравнению маятника с затуханием, и анализ потери устойчивости режима стационарного равномерного падения, что может приводить к переходу в режим авто ротации. В разделе 4 рассмотрена модификация модели, состоящая в том, что исключается гравитация, принимается постоянной циркуляция вектора скорости вокруг профиля и к телу прикладывается постоянный момент внешней силы. Показано, что в этой ситуации возможна хаотическая динамика, ассоциирующаяся с присутствием странного аттрактора типа Лоренца. В разделе 5 представлена модель Танабе-Канеко, учитывающая действие на тело сил сопротивления и подъемной силы из-за наличия циркуляции, выражающейся через динамические переменные на базе постулата Кутты - Жуковского -Чаплыгина. Изложены критические замечания, высказанные в литературе в адрес этой модели, и обозначены пути модификации для их учета. Представлен приближенный анализ нестандартной бифуркации, сопровождающей переход от режима стационарного падения к колебательному режиму в модели Танабе-Канеко. В разделе 6 обсуждается модель Бель-монте-Айзенберга-Мозеса, в которой, в отличие от модели Танабе-Канеко, сила и момент силы сопротивления считаются зависящими от скорости по квадратичному закону. В разделе 7 рассмотрена модель Андерсена-Песавенто-Ванга, в которой заложены эмпирически подобранные авторами на основании численного решения уравнений Навье-Стокса зависимости циркуляции и сил сопротивления от динамических переменных, характеризующих движение тела. В разделе 8 формулируются обобщенные уравнения, позволяющие с использованием одной и той же системы безразмерных переменных и параметров провести сравнительный анализ динамического поведения перечисленных моделей.

1. Плоская задача о падении тела в идеальной жидкости

Будем использовать две системы координат — лабораторную (X, У), где положение центра масс тела дается декартовыми координатами X по горизонтали и У по вертикали, и подвижную систему отсчета (х,у), оси координат которой фиксированы относительно тела (рис. 1).

Рис. 1. Лабораторная и подвижная системы координат в плоской задаче о падении тела в сопротивляющейся среде.

¥к

У

Обратимся вначале к задаче о падении тела в идеальной невязкой жидкости. Кинетическая энергия системы с учетом инерционных свойств тела и среды записывается в виде

Т = \ (ш + ///.,-) г; + i (m + my)v2y + \{I + J)02, (1.1)

где vx и vy — компоненты скорости в проекции на оси подвижной системы отсчета, 0 — угловая скорость, представляющую собой производную по времени от угла поворота тела, m и I — масса и момент инерции тела. Добавки mx, my, J учитывают эффект присоединенных масс, обусловленный вовлечением в движение прилегающих объемов жидкости. Будем считать, что тело имеет плотность ps, а жидкость — плотность pf .В соответствии с двумерной постановкой задачи, плотность определена как масса на единицу площади в поперечном сечении.

Согласно классической гидродинамике [7-12], при движении тела в идеальной невязкой жидкости силы сопротивления отсутствуют (парадокс Даламбера-Эйлера), а к появлению действующей на тело силы приводит наличие ненулевой циркуляции поля скорости Г вдоль контура, охватывающего тело (теорема Жуковского). Величина Г не зависит от конкретного выбора контура и при этом не меняется во времени, будучи определена начальными условиями для поля скоростей в жидкости (теорема Кельвина - Гельмгольца). Составляющие действующей на тело силы даются соотношениями fx = -pfrvy и fy = pfrvx.

С учетом сказанного, в присутствии силы тяжести, направленной по вертикали вниз и характеризуемой ускорением свободного падения g, с поправкой на архимедову выталкивающую силу, уравнения для составляющих скорости и угловой переменной записываются в виде

(m + mx)vx = (m + my)vy0 — pf rvy — mg(l — pf p-1) sin 0,

(m + my )vy = —(m + mx)vx0 + pf rvx — mg(1 — pf p-1) cos 0, (1.2)

(I + J )0 = (mx — my )vx vy.

Изменение во времени координат центра масс в лабораторной системе отсчета определяется тогда дифференциальными уравнениями

X = vx cos 0 — vy sin 0, Y = vx sin 0 + vy cos 0. (1.3)

Для тела с профилем в форме эллипса с полуосями a и b масса и момент инерции выражаются как

m = TTpsab, I = - 7Tpsab(a2 + б2), (1.4)

а присоединенные массы и присоединенный момент инерции в случае идеальной невязкой жидкости в плоской задаче определяются соотношениями [9]

тх = 7г pfb2, ту = 7г pfü2, J = i 7Tpf(a2 — b2)2. (1.5)

8

Подстановка этих выражений в уравнения (1.2) дает

rv

Avx = BvyÓ — p¡3~1—I — gil — p) sinó1, na2

Bi>4 = -Avxé + р[Г1Щ -g(l -р)со&в, (L6)

na2

Q0 = a-2(A — B)vx vy,

где введена сокращенная запись коэффициентов

A = l + pß, B = l + pß-\ Q = -^l+ß2) + ^pß-l{l-ß2)2 (1.7)

и обозначения р = pf /ps для отношения плотностей жидкости и тела и ß = b/a для отношения полуосей эллипса.

В некоторых случаях удобно использовать дополнительно нормировку составляющих скорости u = vx/a, v = vy/а, что приводит уравнения (1.6) к виду

Ай = BvÓ — p¡3 1—о — 9a 1 (1 — р) sin па2

Bv = -Айв + - да~1{1 - p) cos в, (L8)

na2

Q0 = (A - B)uv. 2. Консервативная динамика без гравитации

Рассмотрим сначала случай равных плотностей жидкости и тела, когда действие силы тяжести компенсировано выталкивающей силой, то есть положим р = 1. Тогда уравнения (1.8) принимают вид

Aii = BvO - T'v,

BV = -Aud + T'u, (2.1)

QO = (A - B)uv, где Г' = Г/па2в = Г/па&. Подстановка

u = RA-1 cos р, v = -RB-1 sin p (2.2)

приводит уравнения к виду

A-B

R = Г R——— sin ш cos (p, AB ^ ^

p = O - Г' (A-1 cos2 p + B-1 sin2 p), (2.3)

я A - B 2 e = ~^BQR

Переменная p имеет смысл угла между вектором скорости поступательного движения тела и главной осью эллиптического профиля, а переменная O — это угол наклона главной оси эллипса к горизонтали в лабораторной системе отсчета.

Если умножить первое уравнение на R и сложить с третьим уравнением, то можно видеть, что обращается в нуль производная от величины

D2 = R2 + 2^Q O, (2.4)

то есть это интеграл движения. Выразив отсюда R, с помощью второго и третьего соотношений системы (2.3) получаем уравнение для угловой переменной р:

A - B ( , A + B - A - B

0=2QÄB\-d + г HB -r Q^ÄB cos2" Isi"<2'5>

При нулевой циркуляции Г' = 0 уравнения (2.1) для переменных и, V, и> = в принимают вид [7]

AU = Bvw, Bv = -Auw, QW = (A — B )uv.

(2.6)

Согласно (2.4), в качестве интеграла движения в этом случае выступает величина Я, что соответствует сохранению импульса поступательного движения тела с учетом присоединенных масс. Из второго уравнения (2.3) видно, что ф = в, так что, не нарушая общности, можно не различать ф и в. (Это связано с тем, что начало отсчета угловой переменной в можно брать произвольным в силу пространственной изотропии задачи без гравитации.) Тогда для переменной в имеет место уравнение, приведенное в книге Ламба,

9 + ^А. R2 cos 9 sin 9 = 0. ABQ

(2.7)

Замечательно, что для удвоенного угла § = 29 это уравнение, в силу тождества sin 9 cos 9 = = sin 29, совпадает с уравнением математического маятника. На рисунке 2 показан фазовый портрет на плоскости переменных (9,9). В силу того, что 9 имеет смысл угловой переменной, конфигурация фазового портрета периодически повторяется, поэтому картину можно мыслить как развертку поверхности цилиндра, получаемого склеиванием вертикальных краев изображенного прямоугольника.

-1.6

0 7Г в

Рис. 2. Фазовый портрет уравнения (2.7); в = 0.25, R = 1.

Неподвижные точки 9 = п/2 и 9 = 3п/2 представляют собой устойчивые состояния типа центр, которые отвечают равномерному движению тела широкой стороной вперед. Неподвижные точки 9 = 0 и 9 = п — это неустойчивые состояния типа седло, соответствующие равномерному движению профиля узкой стороной вперед. На фазовом портрете можно видеть пару особых кривых — сепаратрис, сходящихся вместе в седловых точках. Движение по

сепаратрисам можно реализовать, задав начальные условия 9о = tj, Wq = Эти

кривые разделяют области колебательных движений около точек-центров, изображаемых замкнутыми кривыми, и ротационных движений, изображаемых кривыми выше и ниже сепаратрис.

Диаграммы на рисунке 3 иллюстрируют движение тела в реальном пространстве и построены по результатам совместного численного решения уравнений (2.6) и уравнений для координат центра масс X = u cos 9 — v sin 9, Y = u sin 9 + v cos 9. Верхний рисунок отвечает

¿0 = 0

0о = О.4

Рис. 3. Движение тела эллиптического профиля в идеальной жидкости без циркуляции при в = 0.25, К = 1. Показано положение главной оси в последовательные моменты времени. Начальный угол во = п/2, а начальные значения угловой скорости указаны на рисунке.

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

При ненулевой циркуляции переменные ф (угол между вектором скорости и главной осью эллипса) и в (угол наклона профиля в лабораторной системе отсчета) ведут себя по-разному, и их следует различать. Для изображения двумерного фазового портрета при фиксированной величине интеграла движения Б по оси абсцисс в этом случае откладываем величину ф. На рисунке 4 приведен фазовый портрет при величине интеграла движения Б = 1 и циркуляции Г = 0.4. Хотя вид траекторий трансформирован в сравнении с рисунком 2, топологическое устройство осталось таким же. Имеются неподвижные точки — состояния типа центр, а также неустойчивые точки типа седло, располагающиеся на пе-

Рис. 4. Фазовый портрет в случае динамики с циркуляцией; Г = 0.4, в = 0.25, Б = 1.

ресечении особых кривых — сепаратрис, разделяющих области колебательных движений, изображаемых замкнутыми кривыми, и ротационных движений, представленных кривыми выше и ниже сепаратрис. Анализ движения в пространстве при наличии циркуляции показывает, что колебательное или ротационное движение тела имеет место на фоне совершаемого центром масс движения вокруг некоторого центра (рис. 5).

Рис. 5. Диаграммы, иллюстрирующие движение тела эллиптического профиля в идеальной жидкости с циркуляцией при Г = 1, в = 0.25 с разными начальными условиями.

Приведенные достаточно простые и наглядные результаты можно рассматривать как отправной пункт для последующего анализа режимов, имеющих место при наличии диссипации и гравитации. Это аналогично тому, как в теории колебаний линейный консервативный осциллятор служит изначальной моделью для последующих модификаций с вовлечением эффектов, ответственных, скажем, за затухание колебаний или возникновение автоколебаний, представленных замкнутыми притягивающими орбитами — предельными циклами [44, 45].

Заметим, что для систем, состояния которых представляются на фазовом цилиндре, в теории колебаний вводят понятия предельных циклов первого рода, не охватывающих цилиндр, и второго рода, обходящих вокруг цилиндра. В контексте задачи о падении пластинки первый случай будет отвечать автоколебаниям при движении без переворачивания (флаттер), а второй — падению с кувырканиями (авторотация).

3. Модель Козлова: падение тела при наличии вязкого трения

Простейшая модель падения тела в жидкости под действием силы тяжести с учетом трения, пропорционального скорости [22], основывается на уравнениях вида (2.1), где циркуляция принимается нулевой, а в уравнения вводятся дополнительные члены, пропорциональные компонентам скорости и угловой скорости.

При учете сопротивления среды, обусловленного вязким трением, кажется естественным отправляться от известной формулы Стокса [7, 9—11]. Для установившегося движения

шара радиуса Ro в вязкой среде со скоростью V эта формула имеет вид

F = -6nnRo V, (3.1)

где п — коэффициент вязкости, который, как полезно вспомнить, представляется как произведение плотности жидкости на коэффициент кинематической вязкости: п = Pf v. Для тела в форме эллипсоида с полуосями a, b, c имеется обобщение [7] — формула того же вида, но с модифицированным геометрическим параметром R*, который в случае движения вдоль полуоси эллипсоида a определяется выражением

оо

1 = [ 3(А + 2а,2) d\

R* J 8л/(А + а,2)3(А + 62)(Л + с2)'

В контексте трехмерной задачи о падении в жидкости это решение можно было бы в соответствующем пределе использовать для тел в форме «блина» (диска эллиптической формы). Однако в строго двумерной постановке оно неприменимо, что является содержанием так называемого парадокса Стокса [7, 10, 11]. Тем не менее, с физической точки зрения, поскольку геометрические размеры тела в реальности ограничены по всем трем измерениям, кажется уместным допущение, что составляющие силы вязкого сопротивления представляются выражениями, аналогичными по структуре формуле Стокса, а именно,

Fx = -cinavx, Fy = -c2 navy, (3.3)

а момент силы вязкого сопротивления вращательному движению — выражением

Me = -ев паЧ, (3.4)

где ci)2,3 — некоторые коэффициенты. Это соответствует добавлению в уравнения (1.2) членов, представляющих собой производные по обобщенным скоростям от функции Рэлея [22]

R = \ crr/avl + \ с2??а/ц2 + \ сз'Г]а392. (3.5)

2 2 у 2

При приведении к форме уравнений, аналогичной (1.8), имеем

AU = —/1и + BvQ — P sin в,

Bv = — Аив — P cos в, (3.6)

Q в = —/з в + (А — B )uv,

где /1,2,3 = П@-1п-1а-1С\,2,з = рв-1 п-1иа-1сл,2,з = рв-1к\,2,з, P = ga-1(1 — р).

В случае, когда отношение коэффициентов трения при движениях тела эллиптического профиля вдоль одной и другой главной оси равно отношению соответствующих эффективных масс, то есть коэффициентов А и B, уравнения (3.6) интегрируются аналитически. В самом деле, пусть / = А/, / = B/, и положим z = (Au + iBv)e-^t-i®. Тогда из первой пары уравнений (3.6) для переменной z находим Z = —iPe^, и z = —iPf e^1 dt. Таким образом, имеем

Au = —P/-1 (1 — esin в + ce-lJt sin(e + а),

i ( Á t (3.7)

Bv = —P/-1 (1 — e-^) cos в + ce-fa cos(e + а),

где с и а — постоянные величины, определяемые начальными условиями. Отсюда видно, что при больших временах Ь реализуется режим падения, в котором переменные связаны

соотношениями

и = -РА-1/л-1 ссе в, V = -РЕ-1/!-1 8Ш в.

-1 ..-1,

(3.8)

Угловая переменная в этом режиме будет подчиняться третьему уравнению системы (3.6), где нужно сделать подстановку (3.8). В результате оно принимает вид уравнения маятника с затуханием

(3.9)

Отметим очевидное соответствие с уравнением (2.3) для консервативного случая. Если начальные условия таковы, что движение тела исходно было ротационным, то со временем оно сначала превратится в колебания без кувырканий, а затем колебания затухнут и будет происходить стационарное падение широкой стороной вперед, что соответствует устойчивой неподвижной точке

и = 0, V = -Р//2, г = 0 = 0, в = 0. (3.10)

(Имеется второе, не отличающееся по своим свойствам, стационарное решение в = п, поскольку возможны две равноправные ориентации тела при падении — одной или другой широкой стороной вниз.)

В работе Козлова исследована устойчивость стационарного решения (3.10) для произвольного соотношения коэффициентов трения /1,2,3. Подставляя в (3.6) решение (3.10) с малыми добавками к переменным, обозначенными тильдой, получаем в первом порядке по возмущениям

1

Аи = — ЕР/2 г — Рв, Ев = —/2 в,

Qw = —ц3ги + (Е — А)Р/-1и, и= гв.

Подстановка в, в, гв, в ~ е^ приводит к характеристическому уравнению

(3.11)

,3 , I , мз'\2 , ми4мз + (В-А)ВР2 (В — А)Р2 А С}) + ^АЯ

0.

(3.12)

Среди четырех корней этого уравнения имеется один тривиальный, Зо = —/2/Е, а три остальные корня получаются из решения кубического уравнения с вещественными коэффициентами. Условие потери устойчивости — перехода действительной части пары комплексно-сопряженных корней через нуль — можно получить, если искать решение в виде = , что приводит к соотношениям

2 , /1/2/3 + (Е — А)ЕР2

—с2 +

ц22AQ

С = 0,

/1 из \ 2 _ (В - А)Р2 = А С?) /2 АЯ

1 +

/1/22/3

(3.13)

(3.14)

откуда в предположении ( = 0 имеем

/2 _ //1 , /з

В ~ \ А Я ) ' (В-А)ВР2

(Заметим, что в интегрируемом случае /12/Е = /1/А равенство невозможно: пока /1,2,3 > 0 и Е > А, правая часть с очевидностью больше левой, так что неподвижная точка (3.10) всегда устойчива.)

Чтобы перейти от тела с профилем эллипса к тонкой пластине, естественно рассмотреть предел в = Ь/а ^ 0.1 Из уравнений (3.6) видно, что в таком пределе содержательным оказывается случай, когда одновременно стремится к нулю также и отношение плотностей, в то время как величина г = р/в, оставаясь фиксированной, рассматривается как параметр задачи. Производя нормировку времени, составляющих скорости и коэффициентов трения согласно соотношениям

и = и' да~1, v = у' \/ да~1, I = д~1а, &1,2,з = да~1, (3.15)

приходим к уравнениям

и = —гк1и + (1 + г)ув — 8Ш в, (1 + г) = —гк2у — ив — сов в,

^ + ^ г ) в = —гкзв — гш

4 8 )

(3.16)

где штрихи для краткости опущены.

Заметим, что уравнения (3.16) соответствуют по форме (3.6), где А = 1, В = 1 + г,

= \ + г, Р = 1, ¿¿1 2 з = ткх 2 3) поэтому все приведенные выше результаты легко на них

4 8 ' ' ' '

переносятся. В частности, это касается редукции к уравнению маятника с затуханием (3.9) в интегрируемом случае, а именно, при к1 = к2/(1+г), а также условия потери устойчивости неподвижной точки, отвечающей режиму стационарного падения (3.14). Последнее, как нетрудно видеть, переписывается в виде

1 + г \ 2+ г) \ 1 + г )

На рисунке 6 иллюстрируется динамика модели (3.16) на фазовой плоскости (в, в), или, как можно интерпретировать, на фазовом цилиндре. В ситуации, представленной на диаграмме (а), движение, будучи вначале ротационным, становится затем колебательным, и затухание колебаний приводит в итоге к неподвижной точке, отвечающей устойчивому стационарному падению профиля широкой стороной вперед. На диаграмме (Ь) траектория уходит по спирали от неустойчивой неподвижной точки и приходит к охватывающему цилиндр предельному циклу второго рода, который соответствует режиму авторотации. На рисунке 7 можно видеть иллюстрации соответствующего пространственного движения пластинки.

Несмотря на простоту модели Козлова, в ней оказываются возможны нетривиальные феномены нелинейной динамики — переход к хаосу через удвоения периода, странные аттракторы, мультистабильность (см. параграф 8.2).

В заключение этого раздела упомянем кратко работу [26], где постановка задачи аналогична модели Козлова, но с добавлением членов, учитывающих нелинейный закон зависимости силы сопротивления от скорости. Автор задает нелинейное трение в таком виде, что

1 Строгости ради заметим, что такой предельный переход приводит на самом деле к пластинке с неоднородным распределением массы (с максимумом в центре и с уменьшением линейной плотности к краям по закону \/о2 — ж2). Количественное отличие от однородной пластинки в том,

что вместо момента инерции ^ та2 для эллиптического профиля в указанном пределе получается то2. Впрочем, пока речь идет о качественном анализе динамики, эту разницу можно считать не принципиальной.

Рис. 6. Фазовые траектории в модели Козлова, заданной уравнениями (3.16), для г = 2, = 1, кз = 0.04, к^ = ^ (а) и к^ = -рг (Ь). Красным цветом показаны аттракторы — неподвижная точка на

3 5

рисунке (а) и предельный цикл второго рода, отвечающий авторотации, на рисунке (Ь).

(а)

I м 1111111111 м 11 м) м м 1 ^ м м м V ^ г 11 +-\ \ м I м м 1111 f —4 ^—^—^——^—ч-—

(Ь)

Рис. 7. Пространственное движение тела в модели Козлова. Значения параметров соответствуют диаграммам (а) и (Ь) на предыдущем рисунке.

при поступательном движении у каждой компоненты силы коэффициент сопротивления зависит только от «своей» координатной составляющей скорости, что вряд ли правомерно. Впрочем, значительная часть приведенных материалов исследования относится к случаю отсутствия нелинейных поправок. Помимо того, в модели учитывается возможная асимметрия — смещение центра масс тела относительно геометрического центра сечения в форме эллипса, а также возможность наличия циркуляции вектора скорости по контуру сечения тела, которая считается заданной постоянным параметром Г. Приводимые аналитические результаты соответствуют результатам Козлова с тем дополнением, что рассмотрен асимметричный случай. Численных результатов представлено мало, и их не удается воспроизвести из-за очевидных ошибок в тексте статьи при указании параметров.

4. Модель с циркуляцией и аттрактор Лоренца

Аттрактор Лоренца [46-48], хорошо изученный математиками классический объект нелинейной динамики, относится к классу квазигиперболических (или сингулярно гиперболических) странных аттракторов. Хаотическая природа динамики на аттракторе Лоренца установлена и строго обоснована [49].

Покажем, что при движении тела эллиптического профиля в жидкости с вязким трением при определенных условиях может реализоваться хаос, ассоциирующийся с аттрактором лоренцевского типа.

В предположениях модели Козлова (3.6) рассмотрим случай равенства плотностей тела и жидкости р =1, когда действие силы тяжести оказывается исключенным (т.е. Р = 0), и будем полагать присутствующей отличную от нуля постоянную циркуляцию Г. В уравнение для угловой скорости добавим постоянный член, отвечающий действию на тело неизменного во времени вращающего момента. Используя обозначение и> = в, получаем систему,

похожую на модель Лоренца:

Ли = Bvw — rv — ¡iu,

Bv = —Auw + Ги — ¡2v, (4.1)

QW = —(B — A)uv — ¡3w + M.

Заменой переменных и параметров

и = хл/A~1BQ/(B - А), V = у \/AB~1Q/(B - A), w = M/j^1-z, (4.2) vi = ¡1Л-1, hi = M¡-1 — ГВ-1,

V2 = ¡2 В-i, h2 =ГЛ-1 — M¡-1, (4.3) V3 = ¡¡3Q-1.

Эта система приводится к виду

x = hiy — vix — yz,

y = h2X — V2y + xz, (4.4)

Z = —V3Z + xy,

что в случае hi = h2 в точности совпадает с уравнениями для действительных амплитуд в задаче о параметрическом возбуждении волн, рассмотренной в работе [50], где продемонстрировано наличие аттрактора типа Лоренца, в частности, при vi = 1, V2 = 4, V3 = 1, h = 5.875.

Приняв произвольно в = 0.25 и, соответственно, полагая Л = 1 + в = 1.25, В = 1 +

+ ¡3~l = 5, Q = i (1 + /32) + i ¡3~l{l — /З2)2 = 0.705, при пересчете параметров согласно (4.3) получаем

¡i = 1.25, ¡2 = 20, ¡з = 0.705, M = 6.903, Г = 19.583. (4.5)

На рисунке 8 показан портрет аттрактора в трехмерном фазовом пространстве системы (4.1), построенный по результатам численного решения уравнений при указанных параметрах. Также приводится график отображения, полученный в соответствии с процедурой, предложенной в исходной работе Лоренца, где по осям координат отложены значения минимумов переменной w, достигаемых по ходу временной эволюции системы. Вид графика с острым максимумом, напоминающий классическое отображение «зуб пилы» [14, 15, 46, 47, 50], свидетельствует, что наблюдаемый аттрактор квазигиперболический, как и классический аттрактор Лоренца.

На рисунке 9 приводятся портреты аттрактора на плоскости переменных и, v и ф = = arg(u + iv), w = в, которые можно сравнить с диаграммами для других моделей, рассмотренных в настоящем обзоре. Из диаграммы (b) видно, что динамику можно интерпретировать как хаотическую ротацию.

Рисунок 10 иллюстрирует реальное пространственно-временное движение, совершаемое телом, когда в подпространстве обобщенных скоростей имеет место динамика на данном аттракторе. Для построения этой диаграммы численное решение системы (4.1) проводится совместно с уравнениями для угловой скорости и координат центра масс

в = w, X = и cos в — v sin в, Y = и sin в + v cos в. (4.6)

М/ и)

-0.2

—ад

71+1

-4.2

-4.2

-0.2

Рис. 8. Портрет аттрактора Лоренца в трехмерном пространстве состояний (а) и отображение для последовательных минимумов переменной и> в процессе временной эволюции (Ь). Графики построены по результатам численного решения уравнений (4.1) при = 1.25, = 20, рз = 0.705, М = 6.903, Г = 19.583.

Рис. 9. Портрет аттрактора системы (4.1) в проекции на плоскости переменных (и, у) и Параметры те же, что на предыдущем рисунке. То обстоятельство, что траектории проходят с пересечением левой и правой границ прямоугольника на диаграмме (Ь), то есть охватывают фазовый цилиндр, свидетельствует о наличии ротации (кувыркания) тела.

Из рисунка видно, что движение сопровождается хаотическими колебаниями и кувырканиями тела, в соответствии с хаотической природой аттрактора типа Лоренца.

Для количественной характеристики хаоса уместно воспользоваться показателями Ляпунова. Спектр показателей Ляпунова полной системы (4.1)-(4.6) включает в себя три нулевых показателя, отвечающих возмущениям типа изменения угла поворота и сдвига центра масс по двум осям координат, и три нетривиальных показателя, относящихся к аттрактору Лоренца подсистемы (4.1). Их вычисление путем совместного численного решения уравнений (4.1) и соответствующих уравнений в вариациях с применением алгоритма Бе-неттина [15, 18] дает

А: = 0.390 ± 0.003, А2 = -0.0002 ± 0.0007, Аз = -6.390 ± 0.003. (4.7)

Рис. 10. Иллюстрация пространственного движения тела в ситуации, когда в пространстве обобщенных скоростей (и, V, ад) имеет место аттрактор Лоренца. Параметры те же, что на предыдущих рисунках.

Присутствие в списке (4.6) положительного показателя свидетельствует о наличии хаоса, характеризуемого экспоненциальным ростом возмущений, при динамике на аттракторе. Второй показатель, равный нулю с точностью до ошибки вычисления, связан с возмущением типа сдвига вдоль траектории. Третий показатель отрицательный и отвечает за приближение траекторий к аттрактору. Тот факт, что сумма показателей отрицательна, свидетельствует о сжатии фазового объема в фазовом подпространстве (и^у'ш). Величина суммы показателей согласуется с аналитическим вычислением дивергенции векторного поля, компоненты которого заданы правыми частями уравнений (4.1): ё1у Е = ди/и + / + + дм= —¡\/Л — ¡12¡В — (при принятых параметрах она равна -6).

5. Модель Танабе — Канеко

Для правильного описания движения тела в вязкой жидкости принципиален учет зависимости циркуляции вектора скорости вокруг тела от динамических переменных и параметров. Для тела в виде плоской пластины величину циркуляции можно определить с привлечением постулата Кутты - Жуковского - Чаплыгина [7-12], заключающегося в предположении об отсутствии сингулярности поля скорости на задней кромке обтекаемого профиля. Тогда становится возможным найти по теореме Жуковского действующую на тело подъемную силу, а также составляющую силы, ответственную за лобовое сопротивление. В работе Танабе и Канеко [23] делается заключение, что учет этих эффектов приводит к возможности возникновения сложной динамики и хаоса при падении тела в жидкости под действием силы тяжести.

При использовании в качестве динамических переменных компонент скорости в проекции на оси системы координат, связанной с телом, уравнения Танабе - Канеко упрощаются в сравнении с оригинальной записью [23] и представляются в виде

11х + к\\пх = виу - д0 вт в + 7грп'у sgnих, уу + к±иу = —9ух — до ссе в — ж~рухИу sgn ух,

(5.1)

в + к±9 = —3-7Г р1 11>хУу,

где р = pfl/m, I = 2а — ширина пластины, до — ускорение свободного падения, т — масса пластины. Коэффициенты kц и k± характеризуют вязкое трение пластины в жидкости при движении в продольном и поперечном направлениях.2 Для определения координат центра масс в лабораторной системе отсчета систему (5.1) нужно дополнить уравнениями

X = vx cos в — vy sin в, Y = vx sin в + vy cos в. (5.2)

В отсутствие трения система (5.1), очевидно, соответствует по структуре уравнениям (1.6), в которых следует положить

А — В Зтг~б

Г = -irlvysgavx, А = В = 1, =--—, д{ 1 - р) = до- (5.3)

Тот факт, что второе и третье равенства (5.3) противоречат друг другу, связан с определенной некорректностью постановки задачи Танабе и Канеко, отмеченной в критическом комментарии [24, 25], а именно: из-за того, что авторы не принимают во внимание эффект присоединенных масс, коэффициенты A и B получаются равными друг другу, что ведет к обращению в нуль коэффициента в третьем уравнении (5.1), являющегося существенным для наблюдаемых эффектов сложной динамики.3 Далее, формула для циркуляции Г согласно (5.3) принимает во внимание вклад только поступательного движения пластины, тогда как присутствует еще и вклад вращательного движения [9, 28]. Кроме того, Танабе и Канеко не учитывают архимедову выталкивающую силу.4

Несмотря на эти, казалось бы, существенные недочеты, модель в качественном отношении дает разумное представление о возможных режимах сложной динамики при падении пластины в жидкости. Это подтверждается также сравнительным анализом скорректированной версии модели (см. параграф 8.3).

Приведем некоторые численные результаты для модели (5.1), приняв выбранные авторами параметры р = 0.1, I = 1, д = 9.8. Величины Щ и к± будем варьировать, анализируя реализующиеся в зависимости от этого динамические режимы.

Для построения отображения Пуанкаре в качестве секущей оказывается удобным взять в четырехмерном фазовом пространстве трехмерную поверхность, заданную условием

5 = sin в = 0. (5.4)

Вычисление отображения Пуанкаре реализуется в виде специальной подпрограммы, выполняющей численное решение системы дифференциальных уравнений методом Рунге-Кут-ты четвертого порядка. Для построения сечения Пуанкаре в соответствии с условием (5.4) используется метод Эно [15, 51]. Согласно этому методу, численное интегрирование дифференциальных уравнений производится до обнаружения такой ситуации, что в процессе вычислений на очередном шаге величина S меняет знак. Тогда последний шаг аннулируется

2Постулируемое Танабе и Канеко равенство коэффициентов трения для поперечного поступательного и вращательного движения вряд ли можно считать всегда оправданным.

3 Происхождение этого члена в предложенных авторами уравнениях иное и определяется допущением о наличии определенного смещения центра приложения подъемной силы относительно геометрического центра пластины по аналогии с теорией аэродинамического обтекания тонкого профиля.

4Последнее, впрочем, без труда исправляется введением в качестве параметра до эффективного ускорения свободного падения д(1 — р), ав некоторых случаях, таких как падение листа в воздухе, может быть вообще несущественным.

и делается дополнительный шаг с помощью того же разностного метода, но за независимую переменную принимается 5, а величина шага дается предыдущим значением этой переменной с обратным знаком. Это обеспечивает возвращение изображающей точки на секущую поверхность 5 = 0, притом согласованное по точности с используемой разностной схемой. Аналогичная процедура выполняется при пересечении фазовой траекторией трехмерной гиперповерхности и = 0, чтобы во избежание потери точности используемой разностной схемы аккуратно локализовать момент прохода и сделать его совпадающим с узлом разностной сетки.

Для построения карты динамических режимов производится сканирование плоскости параметров к± и / = к±/Щ путем перебора узлов сетки с некоторым шагом по двум параметрам. В каждой точке выполняется порядка 103 итераций отображения Пуанкаре, и по результатам последних шагов итераций проводится анализ на предмет наличия периода повторения величины, задающей квадрат угловой скорости, от 1 до 14 с некоторым заданным изначально уровнем допустимой погрешности. При обнаружении периодичности соответствующий пиксель на диаграмме обозначается цветом, определяемым числом изменения знака величины 5 на периоде, и производится переход к анализу следующей точки на плоскости параметров. При этом в качестве начальных условий в новой точке разумно задавать состояние, полученное в итоге итераций в предыдущей точке («сканирование с наследованием»), что в большинстве случаев способствует ускорению сходимости к установившемуся режиму динамики.

На рисунке 11 в центре показана карта режимов, а по периферии — портреты аттракторов, отвечающие представительным точкам на плоскости параметров. Изображение фазовых портретов на плоскости (9, 9) позволяет по их виду легко делать заключение о природе режима падения пластины. Если траектория или множество траекторий, принадлежащих аттрактору, охватывает фазовый цилиндр, то режим классифицируется как авторотация — периодическая (РИ, диаграммы (а), (Ь), (с)) или хаотическая (СИ, диаграммы (ё), (Ь)). Если же охвата цилиндра нет, то это отвечает падению с колебаниями без кувыркания — флаттеру, который может быть периодическим (PF, диаграмма (е)) или хаотическим (CF, диаграмма (£)). Черная область в нижней части карты соответствует режиму стационарного падения без колебаний (SPF).

Спектр показателей Ляпунова полной системы (5.1), (5.2) включает в себя два нулевых показателя, отвечающих возмущениям сдвига координат центра масс X и У, и четыре показателя, относящихся к аттрактору подсистемы (5.1). Для аттракторов, отличных от неподвижной точки, один из этих четырех показателей всегда нулевой; он ассоциируется с возмущением типа сдвига вдоль фазовой траектории.

Наличие разрывной функции sgn в уравнениях (5.1) делает предпочтительным вариант методики расчета показателей Ляпунова, не использующий линеаризованных уравнений в вариациях [15]. При заданных параметрах и начальных условиях с применением подпрограммы вычисления отображения Пуанкаре выполняем совместно итерации для четырех состояний, одно из которых отвечает опорной траектории, а три других — слабо возмущенным относительно него состояниям. После каждой итерации для трех векторов возмущения производится ортогонализация по Граму-Шмидту и приведение к заданной малой фиксированной норме, после чего для продолжения расчетов используются переопределенные векторы возмущения. Три нетривиальных показателя Ляпунова получаются как коэффициенты, характеризующие нарастание или убывание накапливающихся сумм логарифмов отношения норм (после ортогонализации, но до перенормировки) к исходным значениям.

Рис. 11. Карта режимов на плоскости параметров к±_, ] = /кц для модели Танабе-Канеко. Значения остальных параметров: ~р = 0.1, / = 1, д = 9.8. Цвета определяются из анализа периода повторения состояний в сечении Пуанкаре (синий цвет — режимы периода 1 с симметрией, а зеленый — без симметрии, остальные цвета соответствуют большим периодам). Белый цвет отвечает хаосу или нераспознанным регулярным режимам. Черный цвет — режим равномерного стационарного падения (SPF). По периферии показаны портреты аттракторов (РИ — периодическая ротация, СИ — хаотическая ротация, ИР — периодические колебания без кувыркания — флаттер, CF — хаотический флаттер). Сине-белая горизонтальная штриховка отвечает области сосуществования аттракторов, ассоциирующихся с хаотической ротацией и периодическими колебаниями, что иллюстрируется диаграммой (Ь).

При движении по параметру к±, то есть по горизонтальному пути на карте слева направо, в модели Танабе-Канеко можно наблюдать переход к хаосу через каскад бифуркаций удвоения периода.

Для систем с симметрией, к которым относится рассматриваемая модель, симметричный предельный цикл не может претерпеть бифуркацию удвоения периода, но может иметь место бифуркация потери симметрии [52]. На карте она соответствует переходу из синей в зеленую область, а за ней следует каскад бифуркаций удвоения периода с переходом к хаосу.

На рисунке 12 показаны однопараметрическая бифуркационная диаграмма («бифуркационное дерево») и зависимость старшего нетривиального показателя Ляпунова от параметра, иллюстрирующие переход к хаосу через бифуркации удвоения периода. Диаграммы отвечают движению на карте по горизонтали при f = 50 в указанных на графиках пределах по к±. Графики имеют хорошо узнаваемый характерный вид, свойственный сценарию

4.70

4.75

4.80

к± (а)

4.85 4.70

4.75

4.80

4.85

к± (Ь)

Рис. 12. Бифуркационное дерево и зависимость показателя Ляпунова от параметра, иллюстрирующие переход к хаосу через последовательность бифуркаций удвоения периода в модели Тана-бе-Канеко на горизонтальном участке ] = 50 на карте режимов (см. рис. 11) в указанных пределах изменения к±_.

перехода к хаосу через бесконечный каскад бифуркаций удвоения периода, подчиняющийся универсальным закономерностям Фейгенбаума [14, 15, 54, 55]. Это качественное заключение подтверждается полученными в расчетах численными оценками для константы сходимости точек бифуркации к точке накопления (5 ~ 4.67) и константы, характеризующей расщепление ветвей «дерева» (а ~ —2.50).

Рисунок 13 воспроизводит изображения фазовых портретов аттракторов модели Та-набе-Канеко, представленные в их оригинальной работе. По осям координат отложены компоненты скорости поступательного движения в лабораторной системе отсчета согласно (5.2). Показатели Ляпунова этих аттракторов (за исключением нулевых, относящихся к подсистеме (5.2)), полученные путем расчетов по описанной выше методике, приведены в подписи к рисунку. Аттрактор на диаграмме (а) регулярный и состоит из единственной замкнутой траектории, представляя собой предельный цикл, так что старший показатель нулевой, а остальные отрицательные. Аттракторы (Ь) и (с) характеризуются присутствием положительного показателя, что говорит о хаотической природе режима. Полная сумма показателей в каждом случае отрицательна, что свидетельствует о сжатии фазового объема в фазовом подпространстве (и, V, 9, в). Ее величина согласуется с результатом аналитического вычисления дивергенции векторного поля, компоненты которого задаются правыми частями уравнений (5.1): ё1уЕ = —Щ — 2к±.

Обсудим вопрос о природе перехода от устойчивого стационарного падения к падению с колебаниями (флаттеру) в модели Танабе-Канеко. Нестандартная природа бифуркации обусловлена тем, что уравнения содержат разрывную функцию sgn.

Рассматривая решения, близкие к режиму стационарного падения 9 = 0, vx = 0, Vy = —до/к±, положим 9 = §, vx = и, Vy = —д0к-1 + V, где ^ 1, |и| ^ 1, |V| ^ 1. Если пренебречь возмущениями поперечной скорости vy (что при относительно большом коэффициенте к±, по-видимому, разумно), уравнения для остальных переменных можно

-18

Рис. 13. Фазовые портреты аттракторов модели Танабе-Канеко при к = 4.84 (a), 4.9 (b), 5 (с). Остальные параметры: ~р = 0.1, / = 1, g = 9.8, k\\ = lOOfc. По осям координат отложены составляющие скорости поступательного движения в лабораторной системе отсчета X = = vx cos в — vy sin 0, Y = vx sin в + vy cos 0. Показатели Ляпунова: {0, -0.085, -2.026, -7.62} (a), {0.162 ± 0.005,0, —2.23 ± 0.01, —7.78 ± 0.1} (b), {0.345 ± 0.007, 0, —2.36 ± 0.01, —7.86 ± 0.05} (с).

записать в виде

й + k\\u + go£ = тг рд^к±2 sgnи,

(5.5)

где £ = к—1{} + Отсюда для величины £ получается уравнение, совпадающее по форме с уравнением осциллятора с сухим трением [53]

£ + кц£ + ш2£ = Р sgn (5.6)

где и2 = Зтг'рд^-1^2, Р = Зтг2~р2 д^-1^4. При этом, однако, коэффициент, отвечающий за «сухое трение», имеет знак, противоположный стандартной задаче.

Примем сначала для простоты кц = 0. На фазовой плоскости (£, £) в верхней полуплоскости фазовые траектории представлены семейством эллипсов с центром на оси абсцисс в точке £ = Ри)~2 = тт'рдк^2 (точка А на рис. 14), а в нижней полуплоскости — эллипсами с центрами в точке £ = -Рш~2 (точка В на рис. 14). Поскольку Р > 0, движение имеет вид нарастающих колебаний.

Пусть теперь коэффициент кц — положительный. Пока он невелик, вместо эллипсов в верхней полуплоскости будем иметь семейство спиральных траекторий, скручивающихся к точке А, а в нижней полуплоскости — к точке В. Соответствующие частота колебаний и коэффициент затухания определяются мнимой и действительной частями решения

£

£ е

Рис. 14. Фазовые траектории уравнения (5.6) на плоскости (£, £) в случае кц =0 (а) и 0 < к^ < к* (Ь). Горизонтальный отрезок АВ, помеченный синим цветом, состоит из неустойчивых неподвижных точек. Красным цветом на диаграмме (Ь) показан устойчивый предельный цикл. Траектории пересекают ось абсцисс без разрыва.

характеристического уравнения Л2 + к^Л + ш2 = 0: Л\,2 = кц/2 ± ^у ш2 _ к|/4. В отличие

от случая к\\ = 0, теперь осциллирующее движение изображающей точки, посещающей по очереди верхнюю и нижнюю полуплоскость, приводит к предельному циклу, охватывающему отрезок АВ. Ситуация изменится, когда характер приближения к неподвижной точке станет не осциллирующим, а монотонным, что имеет место в случае

к\\

> 4ш2. Возвращаясь

к параметрам модели Танабе-Канеко, видим, что этому соответствует условие

На карте режимов (рис. 11) соответствующая область показана черным цветом. Численные расчеты подтверждают, что внутри этой области реализуется стационарное падение, а при выходе из нее развиваются автоколебания (флаттер).5

6. Модель Бельмонте —Айзенберга —Мозеса

Отправляясь от экспериментов с падающими в жидкости пластинами, авторы работы [27] на основе обобщения и обработки полученного материала обращаются к постановке задачи, в рамках которой силы сопротивления поступательному и вращательному движению зависят от скорости по квадратичному закону, а именно:

^ = ~ а\\р11Уух, Еу = а±р/1Уиу, ^ = ашр^4\9\9, (6.1)

5Упрощенные уравнения (5.1) можно считать приближенно справедливыми, пока < 1 . Поскольку |г?| ~ это подразумевает | < 1- При параметрах ~р = 0.1, / = 1, д = 9.8 имеем А* = 4.36 и = ъ'рдк'^2 « ЗА^2, так что должно быть выполнено условие ЗА+2 < 1. Черная область на приведенной карте режимов располагается в той части плоскости параметров, где это условие с запасом выполнено.

где V = ^и'х + 1'у, скц, а_\_ и аш — постоянные коэффициенты. Уравнения, аналогичные уравнениям Танабе-Канеко, принимают в этом случае вид

ух = --а\\рУих - д* эш в + ттру2 sgn их + вьу, V = +

:ух — — — их — у ОИ1 и ~г /I у и у 0511 их ~г 1>Чу,

1

1

в = — Заш~р1\0\0 + Ы~1ж~рп1хуу,

уу = — - а±рУиу — д* соэ в — жрухуу sgnг^I. — 9ух

(6.2)

где использованы такие же обозначения, как в уравнениях (5.1), причем константа д* характеризует ускорение свободного падения, скорректированное с учетом архимедовой выталкивающей силы.

Особо авторы подчеркивают то обстоятельство, что, благодаря квадратичной зависимости сил сопротивления, при приведении уравнений к безразмерному виду в них появляется известный в гидродинамике характеристический безразмерный параметр — число Фруда Fr: в результате замены переменных и параметров

1- Ег- —

уравнения (6.2) приводятся к виду

ух = (-а\\Уих - эш в + 47гу2 sgn г^/Ег + виу, V = +г Ъу = (—а±Уъу — сов 9 — 4пъхъу sgn vx)/Fr — 9ъх, (6.4)

9 = —12 • аш|9|9/РГ2 + 12пЪхЪу,

где штрихи для краткости опущены. Авторы принимают значения коэффициентов сопротивления

а\\ = 0.88, а± = 4.1, аш = 0.0674. (6.5)

На рисунке 15 показаны диаграммы, полученные при численном решении уравнений (6.4) и воспроизводящие результаты работы [27], — иллюстрации движения падающей пластинки, отвечающие режимам падения с колебаниями (флаттер) и с кувырканием (авторотация). На рисунке 16 приведены соответствующие этим режимам портреты аттракторов в проекции на плоскости (Х,У) и (9,9). Из диаграмм в правой колонке видно, что в первом случае наблюдаемому режиму отвечает предельный цикл первого рода, а во втором — второго рода (охватывающий фазовый цилиндр).

В системе (6.2), как и в исходной версии модели Танабе-Канеко, не учитываются присоединенные массы и момент инерции, а также поправка к циркуляции, обусловленная вращательной составляющей движения. Поэтому имеет смысл скорректировать уравнения в том же стиле, как это было сделано для модели Танабе-Канеко (см. раздел 8).

7. Модель Андерсена — Песавенто — Ванга

Среди обсуждаемых в данном обзоре систем модель Андерсена-Песавенто-Ванга [28] представляется наиболее проработанной попыткой конечномерного описания задачи о падении в жидкости плоской пластины или тела с профилем в виде эллипса. Авторы исходят

Рис. 15. Иллюстрации пространственного движения пластины в жидкости в режиме автоколебаний при Fr = 0.45 (а) и авторотации при Fr = 0.89 (Ь) по результатам численного решения уравнений (6. 4).

-1.4

-1.4

Рис. 16. Фазовые портреты аттракторов, отвечающих автоколебаниям при Fr = 0.45 (а) и авторотации при Fr = 0.89 (Ь) по результатам численного решения уравнений (6.1). Черным и красным цветом показаны сосуществующие аттракторы.

из уравнений (1.2), куда добавляются члены, учитывающие сопротивление среды поступательному движению, Fx, Fy, и момент силы трения Fg, учитывающий сопротивление вращению:

(m + mx)vx = (m + my)vy9 — pf Tvy — mg(1 — pf p-1) sin 9 — Fx,

(m + my)vy = —(m + mx)u9 + pf Tvx — mg(1 — pf p-1) cos 9 — Fy, (7.1)

(I + J )9 = (mx — my )vxVy — Fg.

На основании обработки данных, полученных при численном решении двумерных гидродинамических задач обтекания эллиптических профилей с привлечением уравнений Навье-

Стокса [28, 30], для циркуляции Г и сил сопротивления авторы предлагают следующие аппроксимирующие выражения:

Г = -2 СТа—^= + 2 CRa29

\JV'Í+Vy

Fg = npf a4 + й) •

При конкретных вычислениях для коэффициентов принимаются численные значения

6 7

Ст = -=, Сц = ж, Са = ~, Св = 1, (7.3)

5 5

а параметры /1,2 варьируются.

Подставим в уравнения (7.1) выражения для массы, момента инерции, присоединенных масс и присоединенного момента инерции (1.4), (1.5), введем обозначения

Г = (3 = - (7.4)

ар$ а

и нормируем составляющие скорости и время так, чтобы обратился в единицу коэффициент при членах, отвечающих за действие силы тяжести:

u

, t'= ty/gl*(l-p). (7.5)

(7-6)

(7-7)

Vff1*^ - P)' а-vV*(l - P)

В результате получаем

(Г + 02)й = (F + l)v9 - Tv - sin 9 - Ku, (F + 1 )v = -(/* + l32)uÓ + Tu -cosd- Kv,

\ (1 + ¡32)F + i (1 - /З2)2) 9 = (f32 - 1 )uv - KoO.

В пределе в 0, соответствующем переходу к задаче о падении тонкой пластины,

Гй = (Г + l)vÓ - Tv - sin 9 - Ku, (I* + l)v = -I*uÓ + Tu - cos 9 - Kv,

\ I* + i) <9 = -uv - Keé. В соответствии с (7-2) в уравнениях (7-6) и (7-7) следует положить

т = 1(-сТП^== + Снё), к = ±(са-св14^УЛ2^2, , ч

vr V Vи + V2 ) IT V u2 + v2J (7.8)

Kg = ^ + ¡J,2191, Aii = Ail - P).

Для определения безразмерных координат центра масс падающего тела совместно с (7-6) или (7-7) решаются уравнения

Х = u cos 9 — v sin 9, Y = u sin 9 + v cos 9. (7-9) _НЕЛИНЕЙНАЯ ДИНАМИКА. 2015. T. 11. №1. С. 3-49_^

vy

v

x

Андерсен, Песавенто и Ванг, особо не акцентируя этот момент, ограничивают рассмотрение специальным случаем, когда между введенными коэффициентами, случайно или нет, имеет место соотношение Ca + Cb = 2Ct . Благодаря этому некоторые члены в уравнениях сокращаются, и они приводятся к виду

273

Гй = (I* - 1)7)0--/ - sin 0,

5^Vu2 + v2

1472 + 1 2v;2

(Г + 1)7' = (2 - />0 - — -—^=====v - cos 0, (7.10)

5^Vu2 + v2

'1.1

а вывод условия потери устойчивости равномерного стационарного падения существенно упрощается. При линеаризации системы (7.10) в окрестности неподвижной точки u = 0, v = —V = — у57г/12, 0 = 0, 0 = 0 первое и третье уравнение дают

1*и = -(/* - 1)V9 - 0, Qr + 0 0 = Vu-^в, (7.11)

и подстановка u, в ~ esí приводит к характеристическому уравнению

s3 Q I* + ^ /* + sViГ + ,(Г - 1)F2 + V = 0. (7.12)

Условие потери устойчивости, состоящее в переходе действительной части пары комплексно-сопряженных корней через нуль, можно найти, если искать решение в виде s = iZ, где ( — действительное число. Это приводит к соотношениям

(3 Г + 1)Г = (Г " 1)У2(' = V' (7ЛЗ) откуда для критического значения коэффициента трения получаем [28]

21* + 1 1 / 3 21* + 1

Потеря устойчивости неподвижной точки при переходе в область /1 < /* сопровождается суперкритической (нормальной) бифуркацией Андронова-Хопфа [19, 44, 45], в результате которой рождается предельный цикл, отвечающий периодическим автоколебаниям пластины из стороны в сторону (флаттер). Заметим, что такая бифуркация имеет место именно и только при учете вязкой составляющей силы сопротивления вращению тела. Если зависимость момента силы сопротивления от угловой скорости квадратичная, то есть /1 = 0, то вблизи состояния стационарного равномерного падения в линейном приближении сопротивление отсутствует; в этом приближении динамика вблизи неподвижной точки консервативная и, формально говоря, анализ устойчивости должен проводиться с учетом нелинейных членов. При этом авторотация не имеет порога по коэффициенту силы сопротивления.

На рисунке 17 показана карта режимов модели (7.10) на плоскости параметров I* и /, построенная в предположении /1 = /2 = /, что соответствует двумерной фазовой

Рис. 17. Карта режимов на обсуждавшейся в работе Андерсена-Песавенто-Ванга плоскости параметров I* и р для модели (7.10), где ¡1 = ¡2 = ¡. По периферии рисунка показаны портреты аттракторов в представительных точках (РИ — периодическая ротация, СИ, — хаотическая ротация, ИР — периодические колебания без кувыркания (флаттер), CF — хаотический флаттер). Карта с раскраской, определяемой периодом при итерациях отображения Пуанкаре, получена сканированием в направлении снизу вверх с наследованием. Синий цвет обозначает режимы периода 1 с симметрией, а зеленый — без симметрии, остальные цвета отвечают большим периодам. Белый цвет — хаос или нераспознанные регулярные режимы. Черный цвет — режим равномерного стационарного падения (БРР). Пунктирная линия отвечает порогу потери устойчивости неподвижной точки (7.14).

диаграмме в работе Андерсена-Песавенто-Ванга [28]. Карта построена по методике, описанной выше для модели Танабе-Канеко, при сканировании в направлении снизу вверх с наследованием. Для раскраски областей на карте анализируется периодичность или отсутствие таковой при численных итерациях отображения Пуанкаре, заданного условием прохождения переменной 9 нулей функции sin 9. Способ реализации отображения Пуанкаре в численных расчетах аналогичен описанному выше для модели Танабе-Канеко, но проще, поскольку уравнения не содержат разрывных коэффициентов.

Пунктирная линия на карте обозначает порог потери устойчивости неподвижной точки (7.14). Переход через эту границу из черной области SPF, отвечающей падению пластины без колебаний, в темно-синюю область сопровождается возникновением автоколебаний (флаттер, область PF). Области РЯ соответствуют периодическим режимам авторотации. Синий цвет обозначает режимы, для которых система демонстрирует однооборотный предельный цикл с симметрией, а зеленый — без симметрии. Остальные цвета отвечают большим периодам. Белый цвет обозначает режимы, где периодичность не обнаруживается, соответствующая область по большей части занята режимами с хаотической динамикой.

На рисунке 18 показаны фазовые портреты аттракторов модели (7.10) в нескольких представительных точках плоскости параметров (см. пояснения к рисунку). По осям координат отложены составляющие безразмерной скорости в лабораторной системе отсчета, чтобы можно было качественно сравнить структуру наблюдаемых аттракторов с портретами, приведенными выше для модели Танабе - Канеко. На диаграммах (с)-(е) можно видеть сосуществующие аттракторы, представленные черным и красным цветом.

Рис. 18. Фазовые портреты аттракторов модели Андерсена-Песавенто-Ванга (7.10) при следующих значениях /1,2: а) 0.35, Ь) 0.25, с) 0.23, ё) 0.13, е) 0.12, Г) 0.1. Остальные параметры: I * = 1.3, Ст = 1.2, С я = п, С а = 1.4, С в = 1.0. На диаграммах (а)-(е) можно видеть сосуществующие аттракторы, представленные черным и красным цветом. По осям координат отложены составляющие безразмерной скорости в лабораторной системе отсчета (7.9). Показатели Ляпунова: {0, -0.135, -0.910,-1.075} (а), {0, -0.0159, -0.871, -0.934} (Ь), {0, -0.158, -0.792, -0.827} (с), {0, -0.0353, -0.362, -0.766} (ё), {0.051,0, -0.435, -0.735} (е), {0.073,0, -0.474, -0.690} (Г).

Показатели Ляпунова этих аттракторов, полученные в расчетах, приведены в пояснении к рисунку (помимо двух нулевых, относящихся к подсистеме (7.9)). В случаях (е) и (£) аттрактор хаотический, что отвечает положительному старшему показателю Ляпунова. Имеется также нулевой показатель, а остальные — отрицательные. В случаях (а)-(ё) аттракторы представляют собой предельные циклы, у которых старший показатель нулевой, а остальные — отрицательные. Сумма показателей во всех случаях отрицательна, что свидетельствует о сжатии фазового объема, сопровождающем приближение к аттрактору.

8. Обобщенная модель

Сопоставление результатов, относящихся к рассмотренным выше конечномерным моделям падения тела в жидкости, затруднено тем, что авторы используют разную нормировку переменных и параметров. Ниже представлены уравнения, которые охватывают все эти модели и записаны в единообразно определенных безразмерных переменных и параметрах для случая предельного перехода от эллиптического профиля к тонкой пластинке. Принятая нормировка соответствует использованной при записи модели Козлова в версии (3.16). При рассмотрении моделей Танабе - Канеко и Бельмонте-Айзенберга-Мозеса учтены и скорректированы моменты, подвергавшиеся критике [24, 25], так что речь будет идти о модифицированной версии этих моделей.

8.1. Уравнения и нормировка

За исходную форму уравнений примем записанную Андерсеном, Песавенто и Вангом систему (7.1), которая относится к ситуации падения в жидкости тела эллиптического профиля с полуосями a и b при наличии циркуляции Г и сил сопротивления, характеризуемых коэффициентами Kx y g. Полагая

Г = тга2Г, (8.1)

Fx = npf a3Kx u, Fy = npf a3Kyv, Fg = npf a4Kg Q, (8.2)

подставляем в уравнения (7.1) выражения для масс и моментов инерции эллиптического профиля (1.4), (1.5) и, вводя обозначения р = pfр-1 в = b/a, u = vx/a, v = vy/a,

(1 + рр)й = (1 + pfrl)v9 - pf3~lTv - ga~l{ 1 - p) sin9 - р[Г1Кхи,, (1 + pj3~l)v = -(1 + pf3)u9 + p!3~lTu - ga~l{ 1 - p) cos 9 - pfS^KyV,

имеем

\ (1 + ñ + I Pß~1(1 - /З2)2) e = P(ß - ß~l)uv - pß-lKeö.

При переходе к пределу плоской пластины ß ^ 0, р ^ 0, r = pß-1 это дает

й = (1 + r)v9 — гГг> — да-1 sin 9 — гКхи, (1 + r)v = —ив + гТи — да~1 cos в — rKyv,

7 + ^ г ) @ = ~ruv ~ fKgé. 4 8 )

(8.3)

(8-4)

В качестве заключительного шага, перенормировкой времени и составляющих скорости устраняем из уравнений параметр д,

и = и' да~1, v = у' \/ да~1, I = д~1а, (8.5)

и, опуская для краткости штрихи, приходим к системе

й = (1 + r)v9 — гГг> — sin в — гКхи,

(1 + r)v = —ив + гТи — cos в — rKyV, g-j

\ + - г ) в = —гш — гКдв.

4 8 /

Безразмерные координаты центра масс падающей пластины определяются путем решения совместно с (8.6) уравнений

X = u cos в — v sin в, Y = u sin в + v cos в. (8.7)

Если коэффициенты сил сопротивления в (8.6) постоянные,

Kx,y,e = k1,2,3, (8.8)

то при нулевой циркуляции Г = 0 получаем модель Козлова (3.16).

Введем теперь циркуляцию, определяемую движением пластины, на основе постулата Кутты - Жуковского - Чаплыгина. Для этого используем формулу из книги Седова [9], где учитывается вклад поступательной и вращательной компонент движения: Г = —2navy sgn vx + + пва2. В принятой здесь нормировке это дает

Г = —2v sgn u + в, (8.9)

и тогда при задании коэффициентов сопротивления согласно (8.8) уравнения (8.6) будут соответствовать модифицированной модели Танабе — Канеко.

Далее, если при задании циркуляции выражением (8.9) положить

Кх = аЫи2 + V2, Ку = а\\/и2 + V2, Кв = а' \в\,

11 У ^ (8.10) = сц/2п = 0.14, = а±/2п = 0.65, а'ш = 4аш /п = 0.086,

то приходим к модифицированной модели Бельмонте —Айзенберга —Мозеса.

Наконец, уравнения модели Андерсена — Песавенто — Ванга в принятой нами здесь нормировке получаются, если положить

— 2 / т> Л

Г = — ( —Ст—/=== + Сяв ), Ст = 1.2, Ск = тт,

п \ У и2 + V2 )

Кх = Ку = -(сА-Св14^)л/^2Т^, Сл = 1.4, Св = 1, (8-П) п \ и2 + V2 )

Кв = Ж + ц21 в I, = /л 1 \/д~1а.

Все модели характеризуются наличием симметрий, которые полезно принимать во внимание при анализе динамического поведения:

• (Я1). Зеркальная симметрия: и — -и, V — V, в — -в, X — —X, Y — Y (при этом Г — —Г).

• (Я2). Переворот пластины: и — —и, V — —V, в — в+п, X — X, Y — Y (при этом Г — Г).

Как обычно, для любого объекта в фазовом пространстве системы с симметрией, например, регулярного или хаотического аттрактора, имеется две возможности [17]. Первая состоит в том, что объект симметричен, то есть при соответствующем преобразовании переменных переходит сам в себя. Вторая — объект несимметричен, тогда при преобразовании получается другой объект той же природы — симметричный партнер.

8.2. Мультистабильность и хаос в модели Козлова

Хотя на возможность сложной динамики при конечномерном описании падающей в жидкости пластинки первыми указали на основе своей модели Танабе и Канеко, такого рода динамика реализуется, как оказывается, в рамках существенно более простой постановки задачи, без учета циркуляции и связанных с ней действующих на тело сил.

Как отмечалось в разделе 3, в интегрируемом случае к\ = ^/(1 + г) модель Козлова (3.16) редуцируется к уравнению маятника с затуханием и демонстрирует только регулярную динамику. При уходе по параметрам от интегрируемости, можно наблюдать потерю устойчивости неподвижных точек, отвечающих стационарному падению пластины, с возникновением предельного цикла первого рода через бифуркацию Андронова-Хопфа. Он отвечает периодическим автоколебаниям (флаттеру) до тех пор, пока располагается, грубо говоря, внутри сепаратрисы (в терминах редуцированной модели), но с изменением параметров может приблизиться к ней, проходя вблизи седел, что создает предпосылки для возникновения сложной динамики и хаоса.

На рисунке 19 показана карта режимов модели (3.16) на плоскости параметров (г, кз) при фиксированных значениях к1 = 0.06 и к2 = 0.3. Карта получена сканированием снизу вверх с наследованием. Обозначения на карте следуют номенклатуре Танабе и Канеко, использованной выше. По периферии представлены фазовые портреты аттракторов в проекции на плоскость переменных (в, в). Аттракторы, сосуществующие при одних и тех же параметрах, обозначены разными цветами. Пунктирная линия отвечает выполнению условий порога потери устойчивости неподвижной точки (3.17).

Двигаясь на плоскости параметров снизу вверх, наблюдаем сначала бифуркацию потери симметрии режима авторотации (переход из синей в зеленую область), а затем последовательность бифуркаций удвоения периода с переходом к хаосу. Вначале появляется хаотический аттрактор, не обладающий симметрией. Поэтому здесь сосуществует четыре аттрактора, являющихся симметричными партнерами исходного объекта относительно операций Я1 и Я2. Далее происходит слияние в единый симметричный хаотический аттрактор. При движении по параметру наряду с хаосом наблюдаются узкие окна периодичности. После выхода из области сложной динамики попадаем в режим периодических автоколебаний, обозначенный на карте темно-синим цветом, а затем в режим стационарного падения, отвечающий черной области. Следует заметить, что граница раздела этих областей в правой части карты не совпадает с линией потери устойчивости, что, очевидно, связано с изменением характера бифуркации Андронова-Хопфа (суперкритическая, сопровождающаяся мягким рождением или исчезновением предельного цикла, в левой части карты, и субкритическая, отвечающая жесткому переходу и гистерезису, в правой части карты).

На рисунке 20 показаны фазовые портреты аттракторов на плоскости, где по осям координат отложены составляющие поступательной скорости пластины, допускающие сравнение с портретами модели Танабе-Канеко, имеющими аналогичный вид (см. рис. 13). Аттрактор на диаграмме (а) отвечает притягивающему циклу, возникшему после бифуркации потери симметрии и двух бифуркаций удвоения периода. На диаграмме (Ь) — хаотический аттрактор, сформировавшийся в результате каскада бифуркаций удвоения периода. Он несимметричен, то есть наряду с ним имеются еще партнеры, получаемые операциями симметрии Я1 и Я2. На диаграмме (с) можно видеть хаотический аттрактор, который уже симметричен и может интерпретироваться как результат объединения партнеров в единый объект. Показатели Ляпунова этих аттракторов, рассчитанные путем совместного решения

Рис. 19. Карта режимов модели (3.16) при фиксированных значениях к\ = 0.06 и &2 = 0.3, полученная сканированием плоскости параметров (г,кз) в направлении снизу вверх с наследованием, и фазовые портреты аттракторов в проекции на плоскость (в, 0) в представительных точках. Аттракторы, являющиеся симметричными партнерами относительно 81 и 82, обозначены разными цветами. Область РИ отвечает периодическим режимам авторотации, PF — периодическим колебаниям (флаттер), 8PF — падению пластины без колебаний. Обозначения на фазовых портретах СИ и CF относятся к хаотическим режимам ротации и флаттера. Пунктирная линия соответствует порогу потери устойчивости неподвижной точки (3.17).

системы (3.16) и соответствующих уравнений в вариациях с помощью алгоритма Бенетти-на [14, 15, 18], приведены в пояснении к рисунку.

Для выявления мультистабильности и гистерезиса обратимся к однопараметрическим бифуркационным диаграммам. Эффект гистерезиса состоит в том, что, медленно изменяя один или более параметров, с приходом в одну и ту же финальную точку в пространстве параметров получаем разные режимы динамики, то есть разные аттракторы. На рисунке 21 в качестве переменной, откладываемой по вертикальной оси, используем абсолютную величину угловой скорости. При проведении вычислений параметр кз изменяется шаг за

Рис. 20. Фазовые портреты аттракторов модели Козлова при k3 = 0.423 (a), 0.424 (b), 0.437 (c). Остальные параметры: r = 0.9, ki =0.06, k2 = 3. По осям координат отложены составляющие безразмерной скорости поступательного движения в лабораторной системе отсчета X = u cos в — v sin 0, Y = u sin в + v cos 0. Показатели Ляпунова: {0, -0.089, -0.326, -2.106} (a), {0.091,0, —0.326, —2.272} (b), {0.186,0, —0.341, —2.398} (c).

го

(а)

....................

0.45

го

0

||| (Ь)

Ч м

fc3

0.65 0.35

0.55

Рис. 21. Однопараметрические бифуркационные диаграммы, построенные численно для модели Козлова (3.16) для к\ = 0.06, к2 = 3, г = 0.4 (а) и г = 0.9 (Ь). Направление движения по параметру вдоль ветвей диаграммы показано стрелками. Стрелкой ниже оси абсцисс отмечена точка потери устойчивости режима стационарного падения согласно (3.17).

шагом, и при каждом новом его значении начальные условия для итераций задаются состоянием, находимом из финального состояния в предыдущей точке. Данные, полученные при сканировании слева направо и справа налево, помечены стрелками. Несовпадение режимов, наблюдаемых в одном и другом случае, свидетельствует о присутствии при одних и тех же параметрах системы сосуществующих аттракторов, каждый из которых имеет свой бассейн в пространстве состояний — множество начальных точек, стартуя из которых

фазовые траектории приходят к данному аттрактору. Подчеркнем, что эти аттракторы не являются симметричными партнерами по отношению к операциям S1 и S2.

На диаграмме (a) можно видеть единственную область бистабильности, где один из сосуществующих аттракторов — предельный цикл, отвечающий периодическому флаттеру, а второй соответствует периодической (в левой части интервала) или хаотической (у правого края интервала) авторотации. На правом краю диаграммы имеет место бифуркация рождения предельного цикла — суперкритическая (нормальная) бифуркация Андронова-Хопфа.

На диаграмме (b) имеется две области бистабильности. В одной из них, в левой части рисунка, сосуществуют аттракторы, отвечающие периодическому флаттеру (нижняя ветвь) и хаотической авторотации (верхняя ветвь). Во второй области, расположенной справа, сосуществуют предельный цикл, отвечающий периодическим автоколебаниям, и устойчивая неподвижная точка, соответствующая стационарному падению. Это указывает на субкритический характер бифуркации Андронова-Хопфа [19, 44, 45] на пороге потери устойчивости неподвижной точки.

8.3. Модифицированная модель Танабе — Канеко

Модификация модели Танабе - Канеко, обсуждавшаяся в параграфе 8.1, приводит курав-нениям

U = v6 + 2rv2 sgn и — sin в — rkiu,

(1 + r)V = (r — 1)ив — 2ruv sgn и — cos в — rk2v, (g

i + i r J в = —ruv — rk¿9. 4 8 J

Подбором параметра r оказывается возможным добиться неплохого соответствия динамики этой системы и оригинальной модели (5.1).6 На рисунке 22 показана карта режимов модели (8.12) на плоскости параметров (k2, f = k2fk\) при кз = k2/4 и r = 0.1. Ее можно сравнить с картой исходной системы Танабе - Канеко (см. рис. 11: хотя нормировка в уравнениях разная, масштаб подобран так, чтобы сделать сравнение визуально наглядным). Как можно видеть, карты выглядят весьма похожими. В обеих версиях модели наблюдаются аналогичные типы динамического поведения, и имеет место очевидное соответствие взаимного расположения областей различных режимов на плоскости параметров. В частности, на рисунке 23 показаны фазовые портреты аттракторов модели (8.12) при f = k2/k\ = 100, k3 = k2/4 для нескольких значений k2, которые аналогичны аттракторам оригинальной модели на рисунке 13.

Имея в виду последующее сравнение с моделью Андерсена-Песавенто-Ванга, рассмотрим карту режимов модели (8.12) на плоскости параметров (r,k3) при фиксированных значениях ki и k2, которая показана на рисунке 24. Карта получена сканированием снизу вверх с наследованием. По периферии представлены фазовые портреты аттракторов

бУравнения (5.1) заменой переменных и параметров v:x = uy/2gl/3, vy = vy/2gl/3, t = T^/2l/'ig, rk\ = А;ц y/2l/3g, rk>2 = k±y/2l/3g, = k±y/2l/3g приводятся к виду ú + k\u = 6v — sin в +

o _ . • o _ •• • 4 _ 1_

+ ^ Tiplv sgnu, v + koi' = —ви — cosO — ^ npluv sgnu, в + к-^в = — — nphiv. Если положить г = ^

что при принятых в работе [23] параметрах дает r « 0.1, то приведенные выражения весьма близки к уравнениям (8.12), поскольку отличаются только отсутствием недалеких от единицы коэффициентов 1 + r, 1 — r, 1+r/2 . Поэтому динамика моделей (5.1) и (8.12) будет сходной, что и наблюдается в вычислениях.

7.8 к2 26.0

Рис. 22. Карта режимов модифицированной системы Танабе -Канеко (8.12) для г = 0.1 при к\ = = к2//, кз = к2/4 на плоскости параметров (к2,/). Цвета определяются периодом повторения значений квадрата угловой скорости в сечении Пуанкаре (синий цвет — режимы периода 1 с симметрией, а зеленый — без симметрии; остальные цвета соответствуют большим периодам). Белый цвет — хаос или нераспознанные регулярные режимы. Помечены области периодической ротации (РИ), хаотической ротации (СИ) и периодических колебаний (PF). Черный цвет — режим равномерного стационарного падения (SPF).

Рис. 23. Фазовые портреты аттракторов модифицированной модели (8.12) при к2 = 10.86 (а), к2 = 11 (Ь), к2 = 11.6 (с), к\ = к2/100 для г = 0.1. По осям координат отложены составляющие безразмерной скорости в лабораторной системе отсчета для удобства сравнения с диаграммами из работы Танабе-Канеко. Показатели Ляпунова: {0, —0.029, —0.359, —1.906} (а), {0.043, 0, —0.435, —1.862} (Ь), {0.098,0, —0.488, —1.757} (с).

Рис. 24. Карта режимов модифицированной модели Танабе -Канеко (8.12) на плоскости параметров (г, кз) при к\ = 0.06, = 3. По периферии рисунка показаны портреты аттракторов в представительных точках (РИ — периодическая ротация, СИ, — хаотическая ротация, ИР — периодические колебания без кувыркания (флаттер), CF — хаотический флаттер). Карта с раскраской, определяемой периодом изменения квадрата угловой скорости при итерациях отображения Пуанкаре, получена сканированием в направлении снизу вверх с наследованием. Синий цвет обозначает режимы периода 1 с симметрией, а зеленый — без симметрии; остальные цвета отвечают большим периодам. Белый цвет — хаос или нераспознанные регулярные режимы.

в проекции на плоскость переменных (9, 9). Аттракторы, сосуществующие при одних и тех же параметрах, обозначены разными цветами.

Двигаясь на плоскости параметров снизу вверх, то есть при увеличении параметра диссипации, в левой части карты наблюдаем сначала бифуркацию потери симметрии режима авторотации (переход из синей области в зеленую), а затем последовательность бифуркаций удвоения периода с переходом к хаосу. Вначале появляется хаотический аттрактор, не обладающий симметрией, вместе с симметричными партнерами относительно операций Я1 и Я2. Далее происходит слияние в единый симметричный хаотический аттрактор. В области

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

В правой части карты выход из области периодической динамики при движении по параметру снизу вверх сопровождается жестким переходом и гистерезисом; в частности, в точке (0.95, 0.1) сосуществуют периодические режимы флаттера и авторотации. При движении по карте в противоположном направлении, сверху вниз, сначала наблюдаем периодические автоколебания (флаттер), которым отвечает предельный цикл первого рода (темно-синяя область в правом верхнем углу карты). При движении по параметру вниз этот цикл увеличивается в размере и в какой-то момент, связанный, очевидно, с приближением к седловой точке, его топологическая природа меняется, и он превращается в предельный цикл второго рода, соответствующий авторотации (область, показанная более светлым оттенком синего цвета). Далее этот цикл теряет симметрию (переход из синей области в зеленую), после чего имеет место каскад бифуркаций удвоения периода режимов авторотации с переходом к хаотической авторотации.

8.4. Модифицированная модель Бельмонте —Айзенберга —Мозеса

Система уравнений для модели Бельмонте-Айзенберга-Мозеса, модифицированной, как описано в параграфе 8.1, с использованием введенной нормировки переменных и параметров, имеет вид

v, = v0 + 2tv2 Sgn u — Sin 9 — 0!||Г '\/и'2 + v2u,

(1 + r)v = (г — u)é — 2ruv sgn и — cos 9 — а±г\/и2 + v2v, ^g

- + - г | в = —ruv — aur\é\é. 4 8 J

На рисунке 25 показана карта режимов модели (8.13) на плоскости параметров (г,аш), полученная сканированием слева направо с наследованием. По периферии приводятся фазовые портреты аттракторов в представительных точках в проекции на плоскость переменных (9,9). Аттракторы, сосуществующие при одних и тех же параметрах, обозначены разными цветами. Обозначения на карте и портретах соответствуют номенклатуре Танабе и Канеко.

Аттракторы, отвечающие периодическим режимам флаттера (PF) и авторотации (PR), демонстрируют определенное визуальное сходство с портретами на рисунке 16. Однако количественное сопоставление модифицированной модели (8.13) и исходной модели Бельмонте - Айзенберга - Мозеса (6.4) не представляется правомерным. В самом деле, как можно показать, определенное согласно [26] число Фруда выражается через параметры модели (8.13) соотношением Fr = \/ж/г, так что принятым на рисунках 15 и 16 числам 0.89 и 0.45 отвечают большие значения r (приблизительно 4 и 15). Это соответствует ситуации весьма существенного эффекта присоединенных масс, который не учитывается в уравнениях (6.4).

8.5. Модель Андерсена — Песавенто — Ванга в сравнении с другими моделями

Как уже отмечалось, модель Андерсена-Песавенто-Ванга [28] выступает как наиболее проработанный вариант конечномерной модели падения пластины в сопротивляющейся среде. При использовании нормировки предыдущего параграфа уравнения (7.10) записы-

Рис. 25. Карта режимов модели (8.13) при а\\ = 0.14, а± = 0.65 и фазовые портреты аттракторов в проекции на плоскость (в, 0) в представительных точках. Аттракторы, являющиеся симметричными партнерами относительно 81 и 82, обозначены разными цветами. Области РИ отвечают периодическим режимам авторотации, а области PF — периодическим колебаниям (флаттер). Обозначения на фазовых портретах СИ и CF относятся к хаотическим режимам ротации и флаттера.

ваются в виде

и = (1 — т)ув —

2ти3

Ъжл/и'2 + V2

— 8Ш в,

л л 14и2 + 12^2 (1 + г)у = (2г — 1)ив — г- , V — сов в,

5^уи2 + V2

I + I Г) ^ = ~гиг> ~ г +

(8.14)

где т = 1/1 *. Имея в виду сопоставление с другими моделями, естественно исследовать динамику в зависимости от тех параметров, которые присутствуют как в системе (8.14), так и в подвергаемых сравнению моделях: в качестве одного варьируемого параметра будем

Рис. 26. Фазовые портреты аттракторов модели (8.13) при r = 0.6 и различных значениях aw: 0.157 (a), 0.16 (b), 0.2. По осям координат отложены составляющие безразмерной скорости поступательного движения в лабораторной системе отсчета X = u cos в — v sin в, Y = = u sin в+v cos в. Показатели Ляпунова: {0, —0.042, —0.337, —1.493} (a), {0.037,0, —0.414, —1.516} (b), {0.133,0, —0.526, —1.726} (c).

использовать величину г, а в качестве другого — коэффициент линейной или квадратичной силы сопротивления вращению в третьем уравнении. Что касается первых двух уравнений, то их вид и фигурирующие в них коэффициенты будем полагать неизменными.

Потере устойчивости режима стационарного падения отвечает критическое значение коэффициента трения в принятой нормировке

1 / 3 2+Г

На рисунке 27 показана карта режимов на плоскости, где по осям координат отложены коэффициент вязкого трения и параметр г в отсутствие квадратичной силы сопротивления при вращении (/2 = 0). На рисунке 28 приводится карта на плоскости коэффициента квадратичного трения /2 и параметра г в отсутствие вязкого трения (/ = 0). Обе карты получены сканированием снизу вверх с наследованием. Обозначения на картах и портретах аттракторов следуют номенклатуре Танабе и Канеко (РИ и СИ — периодический и хаотический режимы авторотации, PF — периодические автоколебания, SPF — падение пластины без колебаний). По периферии рисунка представлены фазовые портреты аттракторов в проекции на плоскость переменных (9,9). Разными цветами обозначены аттракторы, сосуществующие при одних и тех же параметрах и являющиеся симметричными партнерами по отношению к операциям S1 и S2.

Рис. 27. Карта режимов модели (8.14) при = 0, полученная сканированием плоскости параметров (г, ц{) в направлении снизу вверх с наследованием, и фазовые портреты аттракторов в проекции на плоскость (в, 0) в представительных точках. Аттракторы, являющиеся симметричными партнерами относительно 81 и 82, обозначены разными цветами. Область РИ отвечает периодическим режимам авторотации, PF — периодическим колебаниям (флаттер), 8PF — падению пластины без колебаний. Обозначения на фазовых портретах СИ и CF относятся к хаотическим режимам ротации и флаттера. Наклонными полосами в правой части области SPF обозначена зона сосуществования авторотации или флаттера и режимов стационарного падения.

Карту на рисунке 27 уместно использовать при сопоставлении результатов с моделями Козлова и Танабе-Канеко, а карту на рисунке 28 — при сравнении с моделью Бельмонте-Айзенберга - Мозеса.

В дополнение к модели Андерсена-Песавенто-Ванга интересно рассмотреть приближение, используя выражения для циркуляции и сил сопротивления в предположении, что продольное поступательное движение пластины характеризуется скоростью много большей, чем поперечное. В соответствующем предельном случае |и| ^ |г| для циркуляции получаем

Рис. 28. Карта режимов модели (8.14) при /1 = 0, полученная сканированием плоскости параметров (г, /2) в направлении снизу вверх с наследованием, и фазовые портреты аттракторов в проекции на плоскость (9, 9) в представительных точках. Обозначения те же, что и на предыдущем рисунке.

формулу

- 12

Г =--sgn и + 29 = -0.76 sgn и + 29. (8.16)

5п

По структуре она аналогична формуле Седова (8.9), выведенной в рамках постулата Кут-ты - Жуковского - Чаплыгина и использованной нами в рамках модифицированных моделей Танабе-Канеко и Бельмонте-Айзенберга-Мозеса (8.12) и (8.13). При этом имеет место отличие в числовых коэффициентах: вместо чисел 2 и 1 в выражении (8.9) формула (8.16) содержит константы 0.76 и 2.

Для силы сопротивления поступательному движению в этом приближении имеем

2

кх,у = т~ vii2 + v2 = 0.13 \/и2 + v2, (8.17)

что отличается от модели Бельмонте-Айзенберга-Мозеса числовыми константами (вместо констант 0.14 и 0.65 для продольного и поперечного движения имеет место один и тот же коэффициент 0.13). Таким образом, уравнения, отвечающие данному специальному предельному случаю модели Андерсена-Песавенто-Ванга, записываются в виде

• 2 2

й = (1 — r)v6--ru sgn и — sin 0,

5п

• 14

(1 + r)vy = (2г — 1 )ив — — ruv sgn и — eos 0, (8.18)

5п

+ ^ г) в = -ruv - Г (Vi + /¿2и)

1 1

4

с использованием приближенных выражений для циркуляции и силы сопротивления (8.16) и (8.17).

На рисунке 29 показаны карты режимов на плоскости, где по горизонтали отложен параметр г, а по вертикали — коэффициент вязкого сопротивления вращательному движению тела. Панели (а) и (Ь) отвечают модели Андерсена-Песавенто-Ванга (8.14) и предельному случаю этой модели (8.18) в предположении |и| ^ |г|. Карта, представленная на панели (с), соответствует модифицированной модели Танабе - Канеко (8.12).

На рисунке 30 показаны карты режимов, где по вертикали отложен коэффициент квадратичной силы сопротивления вращательному движению, а вязкое трение отсутствует (/1 = 0). Как и на предыдущем рисунке, панели (а) и (Ь) относятся к модели Андерсена-Песавенто-Ванга (8.14) и предельному случаю этой модели (8.18). Карта на панели (с) соответствует модифицированной модели Бельмонте-Айзенберга-Мозеса (8.13).

Как показывает сравнение панелей (а) и (Ь) на рисунках 29 и 30, общая структура областей в пространстве параметров остается такой же, хотя относительные размеры и взаимное расположение областей на карте несколько изменяется. Это наблюдение свидетельствует в пользу разумности использования описания циркуляции в рамках постулата Кутты - Жуковского - Чаплыгина (по крайней мере, при качественном анализе). Меньшее, хотя заметное, сходство наблюдается на обоих рисунках между диаграммами (Ь) и (с).

Можно обратить внимание на наличие областей PF в правом верхнем и РИ в левом нижнем углу карт, хаотической динамики в центральной части карт, переходы к хаосу через

Рис. 29. Карты режимов, представленные для сравнения модели Андерсена - Песавенто - Ванга (а), предельного случая этой модели (Ь) и модифицированной модели Танабе-Канеко (с) на плоскости параметра г и параметра вязкого сопротивления вращательному движению тела.

(а) (Ь) (с)

Рис. 30. Карты режимов, представленные для сравнения модели Андерсена - Песавенто - Ванга (а), предельного случая этой модели (Ь) и модифицированной модели Бельмонте - Айзенберга - Мозеса (с) на плоскости параметра г и параметра квадратичной силы сопротивления вращательному движению тела.

потерю симметрии и каскад удвоений периода в левой нижней части карты, и т.д. Таким образом, несмотря на определенное количественное различие рассмотренных моделей, в целом картина динамики падающей пластины в сопротивляющейся среде, вырисовывающаяся на их основе, выглядит достаточно определенной и согласованной. Это обстоятельство, очевидно, определяется общностью консервативной динамики, соответствующей уравнениям Кирхгофа, а также общими свойствами симметрии и универсальной природой вовлеченных феноменов нелинейной динамики (неподвижные точки, предельные циклы, аттракторы, бифуркации).

Заключение

В статье воспроизведены и подвергнуты сравнительному анализу результаты исследования плоской задачи о падении пластинки в сопротивляющейся среде на основе моделей в виде обыкновенных дифференциальных уравнений относительно небольшого числа переменных. Методической основой для такого подхода служит то обстоятельство, что в случае идеальной невязкой несжимаемой жидкости для обобщенных координат и скоростей твердого тела применимы уравнения Кирхгофа, отделенные от уравнений динамики жидкости.

В рамках обзора мы старались провести линию рассуждений, отправляясь от ситуации, когда уравнения Кирхгофа редуцируются к уравнению консервативного маятника с нелинейностью синуса к моделям, демонстрирующим автоколебательные и авторотационные периодические и хаотические режимы благодаря учету в рамках тех или иных предположений о силах, действующих на тело со стороны вязкой среды. Подобный подход плодотворен в теории колебаний, когда консервативный осциллятор используется как исходная модель для модификаций с вовлечением эффектов, ответственных за затухание колебаний, возбуждение автоколебаний, возникновение хаотической динамики.

В задаче о падении тела в среде периодические колебания (флаттер) и периодическое движение с кувырканием (авторотация) интерпретируются как ассоциирующиеся с предельными циклами (соответственно, первого и второго рода) в фазовом пространстве системы (точнее, в подпространстве обобщенных скоростей).

Для исследования задачи о падении тела в жидкости привлечен инструментарий численного анализа феноменов нелинейной динамики, включающий в себя построение портретов аттракторов, анализ отображений возврата Пуанкаре, построение карт динамических режимов в пространстве параметров, анализ бифуркационных диаграмм («деревьев»), вычисление спектра показателей Ляпунова. Это позволило наполнить выявляемую на основе конечномерных моделей картину динамики конкретным содержанием и иллюстративным материалом.

Выявлено и заслуживает внимания неожиданное богатство динамического поведения, демонстрируемого простейшей моделью [22], где присутствуют только линейные по обобщенной скорости силы вязкого сопротивления поступательному и вращательному движению пластины, а циркуляция и связанные с ней эффекты типа подъемной силы не учитываются. Ранее внимание обращалось только на простые режимы динамики этой модели, тогда как в ней в определенной области параметров оказываются возможными достаточно сложные феномены, такие как хаос, каскад бифуркаций удвоения периода, мультистабильность.

Обнаружено присутствие странного аттрактора лоренцевского типа в трехмерном пространстве обобщенных скоростей для задачи о движении тела эллиптического профиля в условиях скомпенсированной гравитации, при наличии вязкого трения, приложении к телу постоянного вращающего момента и постоянстве циркуляции вектора скорости вокруг профиля.

В отношении модели Танабе-Канеко, которая учитывает действие подъемной силы на падающую в среде пластинку в рамках постулата Кутты - Жуковского - Чаплыгина, проведен анализ, включающий в себя воспроизведение результатов оригинальной работы, модификацию модели с учетом критики, которой она была подвергнута (учет присоединенных масс), а также сопоставление исходной и модифицированной моделей. Обнаруживается, что, несмотря на казалось бы существенные недочеты, модель в качественном отношении дает разумное представление о возможных режимах сложной динамики при падении пластины в жидкости (по крайней мере, в определенной области параметров). Аналогичное исследование проведено в отношении модели Бельмонте-Айзенберга-Мозеса, в которой вместо вязкого трения введено нелинейное трение, квадратичное по обобщенным скоростям.

Для модели Андерсена-Песавенто-Ванга, опирающейся на предложенные авторами эмпирические формулы для сил сопротивления и подъемной силы, вытекающие из результатов численного решения плоской задачи с привлечением уравнений Навье - Стокса, представлен обширный материал численных расчетов, включая карты режимов в пространстве параметров.

Введена в рассмотрение обобщенная модель, в рамках которой с использованием одной и той же системы безразмерных переменных и параметров удается провести сравнительный анализ динамического поведения для моделей Козлова, Танабе-Канеко, Бельмонте-Айзенберга - Мозеса и Андерсена-Песавенто-Ванга.

Отмечается существенная роль симметрии задачи о падении пластинки в жидкости для понимания наблюдаемых феноменов сложной динамики. В частности, обращает на себя внимание возможность сосуществования аттракторов, представляющих собой взаимно симметричные объекты, и возможность их объединения в единый аттрактор при движении по параметрам.

Обнаруживается, что в целом структура устройства пространства параметров для разных моделей демонстрирует определенное сходство. Таким образом, несмотря на количественное различие в динамическом поведении рассмотренных моделей, вырисовывающаяся на их основе картина движения падающей пластины в сопротивляющейся среде выглядит

достаточно согласованной. Это обстоятельство, очевидно, определяется общностью динамики всех моделей в консервативном пределе (уравнения Кирхгофа), а ее модификация при феноменологическом учете эффектов вязкости имеет общие черты благодаря одинаковой присущей симметрии и универсальной природе вовлеченных феноменов нелинейной динамики (неподвижные точки, предельные циклы, аттракторы, бифуркации).

Выражаю благодарность А.В.Борисову, привлекшему мое внимание к данной проблеме.

Список литературы

[1] Maxwell J.K. On a particular case of descent of a heavy body in a resisting medium // Cambridge and Dublin Math. Journ., 1854, vol. 9, pp. 145-148.

[2] Кирхгоф Г. Р. О движении тела вращения в жидкости // Густав Роберт Кирхгоф. Избранные труды / Е. И. Погребысская, Л. С.Полак. Москва: Наука, 1988. С. 280-300.

[3] Жуковский Н. Е. О парении птиц // Жуковский Н. Е. Собрание сочинений: В 7 тт.: Т. 4: Аэродинамика / под ред. Л. С. Лейбензона, С. А. Христиановича и др. Москва - Ленинград: Гостех-издат, 1949. С. 3-34.

[4] Belmonte A., Moses E. Flutter and tumble in fluids // Physics World, 1999, vol. 12, no. 4, pp. 21-25.

[5] Finn D.L. Falling paper and flying business cards // SIAM News, 2007, vol.40, no.4, 3pp.

[6] Ern P., Risso F., Fabre D., Magnaudet J. Wake-induced oscillatory paths of bodies freely rising or falling in fluids // Annu. Rev. Fluid Mech., 2012, vol.44, pp. 97-121.

[7] Ламб Г. Гидродинамика. Москва-Ленинград: ОГИЗ, 1947. 929с.

[8] Борисов А. В., Мамаев И. С. Динамика твердого тела: Гамильтоновы методы, интегрируемость, хаос. Москва-Ижевск: Инст. компьютерн. исслед., 2005. 576с.

[9] Седов Л. И. Плоские задачи гидродинамики и аэродинамики. Москва: Наука, 1966. 448 с.

[10] Кочин Н. Е., Кибель И. А., Розе Н. В. Теоретическая гидромеханика: Т. 2. Москва: Физматлит, 1963. 728 с.

[11] Биркгоф Г. Гидродинамика: Методы, факты, подобие. Москва: ИЛ, 1963. 244 с.

[12] Аржаников Н. С., Садекова Г. С. Аэродинамика летательных аппаратов. Москва: Высшая школа, 1983. 359 с.

[13] Странные аттракторы / под ред. Я. Г. Синая, Л. П. Шильникова. Москва: Мир, 1981. 256 с.

[14] Шустер Г. Детерминированный хаос: Введение. Москва: Мир, 1988. 253 с.

[15] Кузнецов С. П. Динамический хаос. 2-е изд. Москва: Физматлит, 2006. 356 с.

[16] Kuznetsov A. P., Kuznetsov S.P., Sataev I.R., Chua L.O. Multi-parameter criticality in Chua's circuit at period-doubling transition to chaos // Internat. J. Bifur. Chaos Appl. Sci. Engrg., 1996, vol.6, no. 1, pp. 119-148.

[17] Borisov A. V., Jalnin A. Yu., Kuznetsov S.P., Sataev I.R., Sedova J.V. Dynamical phenomena occurring due to phase volume compression in nonholonomic model of the rattleback // Regul. Chaotic Dyn., 2012, vol. 17, no. 6, pp. 512-532.

[18] Benettin G., Galgani L., Giorgilli A., Strelcyn J.-M. Lyapunov characteristic exponents for smooth dynamical systems and for Hamiltonian systems: A method for computing all of them // Meccanica, 1980, vol. 15, pp. 9-30.

[19] Kuznetsov Yu. A. Elements of applied bifurcation theory. 3rd ed. (Appl. Math. Sci., vol. 112.) New York: Springer, 2004. 631 pp.

[20] Борисов А. В., Козлов В. В., Мамаев И. С. О падении тяжелого твердого тела в идеальной жидкости // Тр. ИММ УрО РАН, 2006, т. 12, №1, с. 25-47.

[21] Borisov A. V., Mamaev I. S. On the motion of a heavy rigid body in an ideal fluid with circulation // Chaos, 2006, vol. 16, no. 1, 013118, 7pp.

[22] Козлов В. В. К задаче о падении тяжелого твердого тела в сопротивляющейся среде // Вестн. Моск. ун-та. Сер. Матем. Механ., 1990, №1, с. 79-86.

[23] Tanabe Y., Kaneko K. Behavior of a falling paper // Phys. Rev. Lett., 1994, vol. 73, no. 10, pp. 13721375.

[24] Mahadevan L., Aref H., Jones S. W. Comment on «Behavior of a falling paper» // Phys. Rev. Lett., 1995, vol. 75, no. 7, p. 1420.

[25] Tanabe Y., Kaneko K. Tanabe and Kaneko Reply // Phys. Rev. Lett., 1995, vol. 75, no. 7, p. 1421.

[26] Mahadevan L. Tumbling of a falling card // C. R. Acad. Sci. Paris, Ser. 2b, 1996, vol. 323, pp. 729736.

[27] Belmonte A., Eisenberg H., Moses E. From flutter to tumble: Inertial drag and froude similarity in falling paper // Phys. Rev. Lett., 1998, vol.81, no. 2, pp. 345-348.

[28] Andersen A., Pesavento U., Wang Z.J. Analysis of transitions between fluttering, tumbling and steady descent of falling cards //J. Fluid Mech., 2005, vol. 541, pp. 91-104.

[29] Pesavento U., Wang Z.J. Falling paper: Navier - Stokes solutions, model of fluid forces, and center of mass elevation // Phys. Rev. Lett., 2004, vol. 93, no. 14, 144501, 4 pp.

[30] Andersen A., Pesavento U., Wang Z. J. Unsteady aerodynamics of fluttering and tumbling plates // J. Fluid Mech., 2005, vol.541, pp. 65-90.

[31] Noor D. Z., Chern M. J., Horng T. L. Study of a freely falling ellipse with a variety of aspect ratios and initial angles // Proc. of the 22nd Internat. Conf. on Parallel Computational Fluid Dynamics http://www1.math.fcu.edu.tw/~tlhorng/paper/Extended_abstract.pdf (2010).

[32] Caetano V. F.R. Calculation of the dynamic behavior of a falling plate or disk in a fluid. https://fenix.tecnico.ulisboa.pt/downloadFile/395142133553/resumo.pdf (2010).

[33] Dupleich P. Rotation in free fall of rectangular wings of elongated shape. http://digital.library.unt.edu/ark:/67531/metadc64688/ (1949, UNT Digital Library).

[34] Huang W., Liu H., Wang F., Wu J., Zhang H. P. Experimental study of a freely falling plate with an inhomogeneous mass distribution // Phys. Rev. E., 2013, vol.88, no. 5, 053008, 7pp.

[35] Mahadevan L., Ryu W. S., Samuel A.D.T. Tumbling cards // Phys. Fluids, 1999, vol.11, no. 1, pp. 1-3.

[36] Field S.B., Klaus M., Moore M.G., Nori F. Chaotic dynamics of falling disks // Nature, 1997, vol. 388, no. 6639, pp. 252-254.

[37] Heisinger L., Newton P., Kanso E. Coins falling in water // J. Fluid Mech., 2014, vol. 742, pp. 243253.

[38] Leweke T., Thompson M. C., Hourigan K. Motion of a Mobius band in free fall //J. Fluids Struct., 2009, vol. 25, no. 4, pp. 687-696.

[39] Lugt H. J. Autorotation // Annu. Rev. Fluid Mech., 1983, vol. 15, no. 1, pp. 123-147.

[40] Paoletti P., Mahadevan L. Planar controlled gliding, tumbling and descent //J. Fluid Mech., 2011, vol. 689, pp. 489-516.

[41] Рамоданов С.М., Тененев В. А. Движение тела с переменной геометрией масс в безграничной вязкой жидкости // Нелинейная динамика, 2011, т. 7, №3, с. 635-647.

[42] Fernandes A. C., Sefat S. M. Bifurcation from fluttering to autorotation of a hinged vertical flat plate submitted to a uniform current // Proc. of the 11th Internat. Conf. on the Stability of Ships and Ocean Vehicles (23-28 September 2012, Athens, Greece), pp. 1-12.

[43] Michelin S., Smith S.G. L. Falling cards and flapping flags: Understanding fluid-solid interactions using an unsteady point vortex model // Theor. Comp. Fluid Dyn., 2010, vol. 24, nos. 1-4, pp. 195200.

[44] Андронов А. А., Витт А. А., Хайкин С.Э. Теория колебаний. Москва: Физматгиз, 1959. 915 c.

[45] Бутенин Н. В., Неймарк Ю.И., Фуфаев Н.А. Введение в теорию нелинейных колебаний. 2-е изд., перераб. Москва: Наука, 1987. 382 с.

[46] Лоренц Э. Детерминированное непериодическое течение // Странные аттракторы / под ред. Я. Г. Синая, Л. И. Шильникова. Москва: Мир, 1981. С. 88-116.

[47] Sparrow C. The Lorenz equations: Bifurcations, chaos, and strange attractors. New York: Springer, 1982. 269 pp.

[48] Shilnikov L. Mathematical problems of nonlinear dynamics: A tutorial // J. Franklin Inst., 1997, vol. 334, no. 5, pp. 793-864.

[49] Tucker W. A rigorous ODE solver and Smale's 14th problem // Found. Comput. Math., 2002, vol. 2, no. 1, pp. 53-117.

[50] Пиковский А. С., Рабинович М. И., Трахтенгерц В. Ю. Возникновение стохастичности при рас-падном ограничении параметрической неустойчивости // ЖЭТФ, 1978, т. 74, с. 1366-1374.

[51] Henon M. On the numerical computation of Poincare maps // Phys. D, 1982, vol. 5, pp. 412-414.

[52] Swift J. W., Wiesenfeld K. Suppression of period doubling in symmetric systems // Phys. Rev. Lett., 1984, vol. 52, no. 9, pp. 705-708.

[53] Магнус К. Колебания: Введение в исследование колебательных систем. Москва: Мир, 1982.

[54] Feigenbaum M.J. Quantitative universality for a class of nonlinear transformations //J. Stat. Phys., 1978, vol. 19, no. 1, pp. 25-52.

[55] Фейгенбаум М. Универсальность в поведении нелинейных систем // УФН, 1983, т. 141, №10, с. 343-374.

Motion of a falling card in a fluid: Finite-dimensional models, complex phenomena, and nonlinear dynamics

Sergey P. Kuznetsov

Kotel'nikov's Institute of Radio Engineering and Electronics of RAS, Saratov Branch

410019 Saratov, Zelenaya 38, Russian Federation

[email protected]

Results are reviewed relating to the planar problem for the falling card in a resisting medium based on models represented by ordinary differential equations for a small number of variables. We introduce a unified model, which gives an opportunity to conduct a comparative analysis of dynamic behaviors of models of Kozlov, Tanabe-Kaneko, Belmonte-Eisenberg-Moses and Andersen - Pesavento - Wang using common dimensionless variables and parameters. It is shown that the overall structure of the parameter spaces for the different models shows certain similarities caused obviously by the same inherent symmetry and by universal nature of the involved phenomena of nonlinear dynamics (fixed points, limit cycles, attractors, bifurcations). In concern of motion of a body of elliptical profile in a viscous medium with imposed circulation of the velocity vector and with the applied constant torque, a presence of the Lorenz-type strange attractor is discovered in the three-dimensional space of generalized velocities.

MSC 2010: 34C15, 76D99, 37E99

Keywords: body motion in fluid, oscillations, autorotation, flutter, attractor, bifurcation, chaos, Lyapunov exponent

Received December 22, 2014, accepted January 16, 2015

Citation: Rus. J. Nonlin. Dyn., 2015, vol. 11, no. 1, pp. 3-49 (Russian)

304 с.