САНКТ-ПЕТЕРБУРГСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ
Направление — Астрономия (011501)
Профиль — звездная астрономия
Асипович Ася Дмитриевна
БЛИЗКИЕ К ПЕРИОДИЧЕСКИМ ОРБИТЫ В
РАВНОБЕДРЕННОЙ ЗАДАЧЕ ТРЕХ ТЕЛ
Дипломная работа
Научный руководитель:
доктор физико-математических наук, профессор
Орлов Виктор Владимирович
Рецензент:
кандидат физико-математических наук, доцент
Мартынова Алия Ибрагимовна
Санкт-Петербург
2016
SAINT PETERSBURG STATE UNIVERSITY
Main Field of Study — Astronomy (011501)
Area of Specialisation — stellar astronomy
Asipovich Asia
ORBITS CLOSE TO PERIODIC ONES IN THE
ISOSCELES THREE-BODY PROBLEM
Graduation Thesis
Scientific supervisor:
doctor of physics and mathematics, professor
Orlov Victor
Reviewer:
candidate of physics and mathematics, assistant professor
Martynova Aliyah
Saint-Petersburg
2016
3
Оглавление
Стр.
Введение . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4
Глава 1. Постановка задачи. . . . . . . . . . . . . . . . . . . . . .
7
1.1 Общие сведения о задаче . . . . . . . . . . . . . . . . . . . .
7
1.2 Основные формулы и обозначения . . . . . . . . . . . . . . .
7
1.3 Работы, посвященные данной задаче . . . . . . . . . . . . . .
10
1.4 Поиск орбит, близких к периодическим. . . . . . . . . . . . .
15
Глава 2. Методы решения задачи и детали их реализации .
16
2.1 Борьба с сингулярностью . . . . . . . . . . . . . . . . . . . .
16
2.1.1
Приближение задачей двух тел . . . . . . . . . . . . .
16
2.1.2
Приближение упругим отскоком . . . . . . . . . . . .
18
2.1.3
Регуляризация . . . . . . . . . . . . . . . . . . . . . .
20
2.2 Программные средства и детали реализации. . . . . . . . . .
20
2.2.1
Пакет MATLAB . . . . . . . . . . . . . . . . . . . . .
21
2.2.2
Программа TRIPLE . . . . . . . . . . . . . . . . . . .
23
2.2.3
Дополнительные программные средства . . . . . . . .
24
Глава 3. Результаты . . . . . . . . . . . . . . . . . . . . . . . . . .
25
3.1 Приближение задачей двух тел . . . . . . . . . . . . . . . . .
25
3.2 Приближение упругим отскоком . . . . . . . . . . . . . . . .
28
3.3 Регуляризация. . . . . . . . . . . . . . . . . . . . . . . . . . .
33
Заключение . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
38
Список литературы
40
. . . . . . . . . . . . . . . . . . . . . . . . . .
4
Введение
Задача трех тел — одна из классических проблем небесной механики.
Впервые поставлена она была Ньютоном, впоследствии ею занимались такие великие ученые, как Эйлер, Лагранж и Пуанкаре. Различным ее приложениям (например, поиску периодических орбит, хореографий, изучению
хаотизации движений и т.п.) посвящено множество статей. Помимо прочих
причин, это связано с тем, что именно подсистемы из трех тел встречаются
чаще всего в звездных системах, количество тел в которых больше двух.
Таким образом, данная задача крайне важна для исследования эволюции
звездных скоплений, галактик, групп и скоплений галактик.
Эйлер, в частности, нашел класс точных периодических решений задачи трех тел — когда в течение всей эволюции системы тела находятся
на одной прямой (сигизии), вращающейся вокруг общего центра масс. Эти
решения называют коллинеарными или эйлеровыми. Также именно Эйлер
первым ввел вращающуюся систему координат, а еще предложил идею регуляризации движений при двойном соударении в задаче двух тел [1].
Лагранж нашел еще один класс точных решений — когда три тела
постоянно находятся в вершинах равностороннего треугольника, обращающегося вокруг центра масс системы. Такие решения называют треугольными или лагранжевыми.
Одним из первых исследователей, достигших значительных результатов в исследовании задачи трех тел был А. Пуанкаре. В своей работе он
искал минимумы функционала действия для плоской задачи трех тел произвольных масс при некоторых наложенных на орбиты ограничениях. При
этом он, помимо прочих результатов, задал направление исследований, которое используется и в некоторых современных работах — рассмотрение
вместо движений тел деформаций треугольников, ими образованных. В
процессе работы им была обнаружена одна из ключевых проблем задачи
трех тел — то, что функционал действия сходится в окрестности столкновений для ньютонова потенциала, пропорционального 1/𝑟. Будучи не в
5
силах решить задачу для этого случая, Пуанкаре рассмотрел так называемое ”сильное взаимодействие“, пропорциональное 1/𝑟2 .
Несмотря на то, что задача трех тел является всего лишь частным
случаем более общей задачи 𝑁 тел, богатство возможных решений и вариантов подзадач поражает воображение. Одной из этих подзадач является так называемая обобщенная задача Ситникова. В своей работе Ситников рассматривал пространственную ограниченную равнобедренную задачу трех тел, когда тело нулевой массы движется вдоль прямой, проходящей
через центр масс двух других тел равных масс перпендикулярно их орбитальной плоскости. Ситников показал, что при подобной постановке задачи
существует несколько вариантов решения:
1. тело нулевой массы покоится в барицентре системы (частный случай решения Эйлера);
2. тело нулевой массы cовершает затухающие колебания около барицентра системы;
3. периодическое решение;
4. тело нулевой массы cовершает колебания с нарастающей амплитудой, завершающиеся его уходом;
5. тело нулевой массы совершает осциллирующие колебания — с нарастающей амплитудой, но без ухода на бесконечность.
Также Ситников установил, что если масса тела, движущегося по
прямой, не равна нулю, но мала по сравнению с массами крайних тел, то
осциллирующие движения тоже возможны.
Важный вклад в изучение обобщенной задачи Ситникова внес Алексеев. Применив методы символической динамики, он показал, что в зависимости от выбора начальных условий в этой задаче реализуются все
возможные типы финальных движений по классификации Шази.
В конце XIX в. Пуанкаре, обобщив теорему Брунса, показал, что общее решение задачи трех тел невозможно выразить через алгебраические
или трансцендентные функции координат и скоростей тел. В 1912 г. Сундман представил общее решение в виде рядов по вспомогательной переменной, однако оказалось, что для практического применения эти ряды непригодны, так как сходятся черезвычайно медленно. Таким образом, так как
6
приемлемое аналитическое решение задачи трех тел невозможно, важную
роль в ее исследовании играют вычислительные методы с использованием
ЭВМ. К сожалению, эти методы сами содержат некоторую погрешность,
а конечная машинная точность добавляет к ней еще и ошибки округления, поэтому следует помнить, что все результаты, полученные с помощью
вычислений на ЭВМ, являются приближенными.
В данной работе будут рассмотрены различные вычислительные методы решения плоской равнобедренной задачи трех тел, проведено их сравнение и представлены результаты поиска периодических орбит в этой задаче.
7
Глава 1. Постановка задачи.
1.1
Общие сведения о задаче
Одним из вариантов задачи трех тел, близким к обобщенной задаче Ситникова, является плоская равнобедренная задача трех тел. В этом
случае масса центрального тела не ограничивается (она может быть как
сравнима, так и больше массы крайних тел), движение всех трех тел происходит в одной плоскости. Таким образом, крайние тела равных масс в
системе координат, связанной с их центром масс, движутся строго по направлению друг к другу и обратно, сталкиваясь и расходясь. Центральное
тело движется по прямой, перпендикулярной соединяющей крайние тела
и проходящей через центр масс системы крайних тел. Все время эволюции
конфигурация тройной системы представляет собой равнобедренный треугольник. Этот вариант задачи трех тел проиллюстрирован на рисунке 1.1
Поиск периодических орбит в этой задаче важен, поскольку с помощью небольшой вариации начальных условий можно перейти к порожденным ими семействам периодических решений для пространственной равнобедренной задачи трех тел или же плоской неравнобедренной задаче трех
тел.
1.2
Основные формулы и обозначения
В данной работе используются следующие начальные условия:
1. Все три тела лежат на одной прямой.
2. Центральное тело движется вертикально вверх.
3. Крайние тела, находящиеся на некотором отдалении от него, движутся под некоторым углом к горизонтальной прямой, либо отдаляясь, либо сближаясь.
8
Рисунок 1.1 — Плоская равнобедренная задача трех тел. CM - центр масс
системы 𝑚1 и 𝑚2
Были введены следующие переменные: 𝑟 — расстояние между крайними
телами, 𝑅 — расстояние от центра масс крайних тел до центрального тела.
При этом массы крайних тел 𝑚1 и гравитационная постоянная 𝐺 принимались равными единице, отношение масс центрального и крайнего тела —
𝜀, а полная энергия системы 𝐸 = −1. Начальные условия при этом можно
задать двумя параметрами [2]:
1. вириальный коэффициент:
𝑟˙02
𝜀𝑅˙ 02 1 + 4𝜀 −1
𝑘=( +
)(
)
4
𝜀+2
𝑟0
2. отношение скоростей:
𝑟˙0
𝜇 = √︁
𝑟˙02 + 𝑅˙ 02
9
При этом 𝑘 ∈ (0, 1), 𝜇 ∈ (−1, 1). Границы интервалов соответствуют тройным соударениям и не рассматривались. В силу симметрии задачи можно
рассматривать только случай 𝑅˙ 0 > 0.
В таком случае начальные условия определяются через параметры
задачи следующим образом:
𝑟0 = (1 − 𝑘)(1 + 4𝜀),
𝑅0 = 0,
√︃
𝑟˙0 = 𝜇
4𝑘
,
4𝜀
(1 − 𝑘)(𝜇2 + 2+𝜀
(1 − 𝜇2 ))
√︃
𝑅˙ 0 =
4𝑘
2
𝜇
(1 − 𝑘)( 1−𝜇
2 +
4𝜀
2+𝜀 )
.
Уравнения движения выглядят так:
⎧
8𝜀𝑟
⎨𝑟¨ = − 22 −
,
𝑟
(𝑟2 +4𝑅2 )3/2
⎩𝑅
¨ = − 8(2+𝜀)𝑅 .
(𝑟2 +4𝑅2 )3/2
Интеграл энергии записывается следующим образом:
1
𝜀 ˙2 1
4𝜀
𝐸 = 𝑟˙ 2 +
𝑅 − − 2
.
4
2+𝜀
𝑟 (𝑟 + 4𝑅2 )1/2
Используются динамические единицы длины и времени:
(1.1)
10
𝐺
,
2|𝐸|
√︃
𝜏=
𝑚𝑖 𝑚𝑗
𝑖̸=𝑗
𝑑=
𝐺
∑︀
3
∑︀
𝑖=1
𝑚𝑖
∑︀
𝑚𝑖 𝑚𝑗
𝑖̸=𝑗
(2|𝐸|)3/2
,
где 𝑑 — средний размер системы, а 𝜏 — среднее время ее пересечения [3].
В условиях задачи они выражаются следующим образом [3]:
𝑑 = 1/2 + 𝜀,
√︀
𝜏 = (1/2 + 𝜀) 1 + 𝜀/2.
1.3
Работы, посвященные данной задаче
Исследованию данной задачи посвящено довольно большое количество статей. Отметим следующие из них:
– Работа [4] была одной из первых в исследовании плоской равнобедренной задачи трех тел. В частности, в ней рассматривался случай
системы с нулевой энергией, а также тройные столкновения. Была также построена регуляризация для данной задачи. Одним из
основных результатов этой работы является обнаружение так называемой орбиты Брука, являющейся периодической. Вид орбиты
в различных системах координат можно увидеть на рисунке 1.2.
– Статья [5] посвящена той же задаче, что и данная работа, но рассматривается только случай равных масс. В ней исследуется зависимость длины выброса центрального тела от начальных условий.
Также рассматриваются различные отклонения от условий задачи — во-первых, смещение центрального тела к одному из край-
11
а)
б)
Рисунок 1.2 — Вид орбиты Брука в различных системах координат.
них, во-вторых, отклонение вектора скорости центрального тела
на небольшой угол. В результате обнаружена заметная зависимость
между начальными условиями и длиной первого выброса центрального тела, проиллюстрированная рисунком 1.3. Также было замечено, что небольшие отклонения от условий задачи могут как увеличивать длину выброса для широких начальных конфигураций, так
и уменьшать ее для тесных. В частности, область ухода центрального тела деформируется и уменьшается при отклонении условий
задачи от равнобедренной.
– В статье [2] была предпринята попытка исследовать зависимость
типа орбиты от начальных данных в плоской равнобедренной задаче трех тел. В отличие от предыдущей работы, в ней рассматривались различные массы центрального тела, а интегрирование орбит
проводилось не до первого разворота центрального тела, а до его
ухода из системы (𝑅 > 100𝑑), либо длительное время (𝑡 ≤ 100𝜏 ).
В случае ухода центрального тела рассматривалось также количество предшествующих этому сближений тел. В результате исследования орбиты были разбиты на несколько классов:
1. орбиты с уходом центрального тела;
2. устойчивые орбиты — периодические и порожденные ими;
3. стохастические орбиты.
12
Рисунок 1.3 — Зависимость длины первого выброса центрального тела от
начальных условий
Примеры данных орбит можно видеть на рисунке 1.4. В результате исследования было обнаружено, что для любого значения массы
центрального тела в центре области начальных условий наблюдается область устойчивых орбит, максимальная при массе центрального тела, равной массе одного из крайних тел. Вокруг этой области расположена область стохастического движения, возрастающая при увеличении массы центрального тела. Области орбит, на
которых центральное тело уходит из системы после некоторого количества "пролетов-(ситуаций, когда центральное тело пересекает
прямую, соединяющую крайние тела), отделены от области устойчивых орбит областью стохастических орбит и уменьшаются с ростом массы центрального тела. Иллюстрацию данных результатов
можно видеть на рисунке 1.5.
13
a)
b)
c)
d)
Рисунок 1.4 — Примеры орбит: a) орбита с уходом, b) периодическая
орбита, c) порожденная периодической устойчивая орбита, d)
стохастическая орбита
14
Рисунок 1.5 — Области начальных условий 𝑘, 𝜇, относящиеся к
траекториям с уходом центрального тела после 𝑛 пролетов. Нижний ряд
демонстрирует области устойчивых орбит. 𝜀 — отношение массы
центрального тела к массе одного из крайних.
15
1.4
Поиск орбит, близких к периодическим.
Как и в других разновидностях задачи трех тел, особый интерес в
плоской равнобедренной задаче представляют периодические орбиты. Их
поиском, в частности, занимался Брук [4]. Как уже упоминалось, им была
обнаружена так называемая орбита Брука (1.2), являющаяся периодической и устойчивой, а также обнаружено несколько других периодических
орбит для систем с 𝐸 ̸= −1, где 𝐸 - полная энергия системы. Однако
вопрос о существовании других периодических орбит в плоской равнобедренной задаче трех тел остается открытым, и именно он рассмотрен в этой
работе.
К сожалению, в силу разностных (обусловленных приближенностью
интегрирования) и машинных (вызванных конечной разрядной сеткой компьютеров) ошибок, любая построенная нами орбита является лишь приближенной. Таким образом, даже получив точное совпадение начальных данных с данными в какой-то момент времени, мы можем утверждать лишь
то, что данная орбита близка к периодической.
Для определения, является ли данная орбита близкой к периодической, был введен параметр орбиты 𝐹𝑚𝑖𝑛 = min 𝐹 (𝑡), где
𝑡
√︃
𝐹 (𝑡) =
˙ − 𝑅˙ 0 )2
(𝑟(𝑡) − 𝑟0 )2 (𝑅(𝑡) − 𝑅0 )2 (𝑟(𝑡)
˙ − 𝑟˙0 )2 (𝑅(𝑡)
+
+
+
𝑑2
𝑑2
(𝑑/𝜏 )2
(𝑑/𝜏 )2
Несложно видеть, что чем меньше 𝐹𝑚𝑖𝑛 , тем точнее орбита в процессе своей
эволюции повторяет начальные данные. Соответственно, если мы найдем
орбиту, для которой в какой-то момент времени значение 𝐹 (𝑡) мало, она
будет являться близкой к периодической.
Следует учитывать, что в начальный момент времени 𝐹 (𝑡0 ) = 0, затем
значение 𝐹 (𝑡) возрастает вплоть до первого максимума в момент времени
𝑡1 . Вычислять 𝐹𝑚𝑖𝑛 следует начиная с этого момента времени и далее.
16
Глава 2. Методы решения задачи и детали их
реализации
2.1
Борьба с сингулярностью
Как видно из уравнений движения (1.1), при 𝑟 = 0 возникает деление
на 0. В отличие от общей задачи трех тел, где столкновения тел являются скорее исключением, чем правилом, в плоской равнобедренной задаче
трех тел двойные столкновения крайних тел случаются постоянно. Таким
образом, обычное интегрирование уравнений движения в данной задаче
неприемлимо, необходим метод устранения особенности.
В данной работе применялись три метода обработки двойных столкновений тел:
– приближение задачей двух тел;
– приближение упругим отскоком;
– регуляризация (с помощью программы TRIPLE).
Рассмотрим их подробнее.
2.1.1
Приближение задачей двух тел
При достижении 𝑟 < 𝑟𝑐𝑟𝑖𝑡 будем считать, что центральное тело слабо
влияет на крайние, и движение крайних тел можно приближенно рассматривать как прямолинейную задачу двух тел. При этом уравнение движения
центрального тела оставим прежним. Как только 𝑟 станет больше 𝑟𝑐𝑟𝑖𝑡 , мы
снова вернемся к рассмотрению задачи трех тел.
Уравнение относительного движения в прямолинейной задаче двух
тел выглядит так:
𝑟¨ = −
𝐺(𝑚1 + 𝑚2 )
𝑟2
17
В нашем случае
2
𝑟2
Интеграл энергии для задачи двух тел имеет следующий вид:
𝑟¨ = −
2
1
𝐸2 = 𝑟˙ 2 −
2
𝑟
Здесь по-прежнему присутствует особенность при 𝑟 = 0. Однако его можно
регуляризовать следующей заменой переменной времени [6]:
𝑑𝑡 = 𝑟𝑑𝜏
В таком случае уравнение движения примет вид
𝑑2 𝑟
= 𝑟′′ = 2 + 2𝐸2 𝑟
2
𝑑𝜏
Это уравнение не имеет особенностей при 𝑟 = 0 и может быть решено в
явном виде. Возможны три варианта:
1. 𝐸2 < 0
𝑟′′ + 𝑎2 𝑟 = 2,
𝑟(𝜏 ) =
2
𝑎2
√
−2𝐸2
+ 𝑐1 𝑐𝑜𝑠(𝑎𝜏 ) + 𝑐2 𝑠𝑖𝑛(𝑎𝜏 )
𝑐1 = 𝑟 |𝜏 =0 − 𝑎22 ,
𝑡 − 𝑡0 =
𝑎=
∫︀ 𝜏
0
𝑐2 =
𝑟(𝜉)𝑑𝜉 =
2𝜏
𝑎2
𝑟′ |𝜏 =0
𝑎
+
𝑐1
𝑎 𝑠𝑖𝑛(𝑎𝜏 )
2. 𝐸2 > 0
𝑟′′ − 𝑎2 𝑟 = 2,
𝑎=
√
2𝐸2
𝑟(𝜏 ) = − 𝑎22 + 𝑐1 𝑒𝑎𝜏 + 𝑐2 𝑒−𝑎𝜏
−
𝑐2
𝑎 (𝑐𝑜𝑠(𝑎𝜏 )
− 1)
18
𝑐1 =
𝑟|𝜏 =0 + 𝑎22 +
2
𝑡 − 𝑡0 =
∫︀ 𝜏
0
𝑟 ′ |𝜏 =0
𝑎
,
𝑐2 =
𝑟|𝜏 =0 + 𝑎22 −
2
𝑟(𝜉)𝑑𝜉 = − 2𝜏
𝑎2 +
𝑐1 𝑎𝜏
𝑎 (𝑒
𝑟 ′ |𝜏 =0
𝑎
− 1) −
𝑐2 −𝑎𝜏
𝑎 (𝑒
− 1)
3. 𝐸2 = 0
𝑟′′ = 2
𝑟(𝜏 ) = 𝜏 2 + 𝑐1 𝜏 + 𝑐2
𝑐1 = 𝑟′ |𝜏 =0 ,
𝑡 − 𝑡0 =
∫︀ 𝜏
0
𝑐2 = 𝑟 |𝜏 =0
𝑟(𝜉)𝑑𝜉 = 31 𝜏 3 +
𝑐1 2
2𝜏
+ 𝑐2 𝜏
Следует заметить, что данный метод применим не всегда. Если в момент сближения крайних тел центральное тело находится достаточно близко к их центру масс, то пренебрежение его влиянием может существенно
изменить орбиты тел. Поэтому для применения метода необходимо отслеживать параметр 𝑅/𝑟, и если его значение в момент 𝑟 = 𝑟𝑐𝑟𝑖𝑡 меньше некоторого значения, то данный метод неприменим. Для интегрирования подобных орбит нужен другой метод борьбы с сингулярностью.
2.1.2
Приближение упругим отскоком
При достижении 𝑟 < 𝑟𝑐𝑟𝑖𝑡 будем считать, что центральное тело слабо
влияет на крайние. Зафиксируем текущие ускорения крайних тел и будем
считать что вплоть до столкновения их ускорения остаются постоянными.
В момент столкновения происходит упругий отскок — скорости крайних
тел меняют направления на противоположные, ускорения остаются прежними. Уравнение движения центрального тела не меняется. Когда 𝑟 вновь
достигает значения 𝑟𝑐𝑟𝑖𝑡 , возвращаемся к задаче трех тел.
19
Таким образом, в нашей системе координат движение крайних тел
будет выглядеть следующим образом:
𝑟¨(𝑡0 )
(𝑡 − 𝑡0 )2 , (𝑡 − 𝑡0 ) ≤ 𝑡𝑐𝑜𝑙𝑙𝑖𝑠𝑖𝑜𝑛
2
𝑟¨(𝑡0 )
(𝑡 − 𝑡𝑐𝑜𝑙𝑙𝑖𝑠𝑖𝑜𝑛 )2 , (𝑡 − 𝑡0 ) > 𝑡𝑐𝑜𝑙𝑙𝑖𝑠𝑖𝑜𝑛
𝑟(𝑡) = −𝑟(𝑡
˙ 𝑐𝑜𝑙𝑙𝑖𝑠𝑖𝑜𝑛 )(𝑡 − 𝑡𝑐𝑜𝑙𝑙𝑖𝑠𝑖𝑜𝑛 ) +
2
(2.1)
𝑟(𝑡) = 𝑟(𝑡0 ) + 𝑟(𝑡
˙ 0 )(𝑡 − 𝑡0 ) +
Где 𝑡𝑐𝑜𝑙𝑙𝑖𝑠𝑖𝑜𝑛 - момент столкновения крайних тел. Его можно вычислить из условия 𝑟(𝑡𝑐𝑜𝑙𝑙𝑖𝑠𝑖𝑜𝑛 ) = 0. Подставляя в 2.1, получим квадратное
уравнение:
(𝑡𝑐𝑜𝑙𝑙𝑖𝑠𝑖𝑜𝑛 − 𝑡0 )2 +
2𝑟(𝑡
˙ 0)
2𝑟(𝑡0 )
(𝑡𝑐𝑜𝑙𝑙𝑖𝑠𝑖𝑜𝑛 − 𝑡0 ) +
=0
𝑟¨(𝑡0 )
𝑟¨(𝑡0 )
Его положительное решение выражается следующей формулой:
𝑡𝑐𝑜𝑙𝑙𝑖𝑠𝑖𝑜𝑛 − 𝑡0 =
√︀
1
(−𝑟(𝑡
˙ 0 ) + 𝑟˙ 2 (𝑡0 ) − 2¨
𝑟(𝑡0 )𝑟(𝑡0 ))
𝑟¨(𝑡0 )
Теперь несложно вычислить скорость в момент столкновения:
𝑟(𝑡
˙ 𝑐𝑜𝑙𝑙𝑖𝑠𝑖𝑜𝑛 ) = 𝑟(𝑡
˙ 0 ) + (𝑡𝑐𝑜𝑙𝑙𝑖𝑠𝑖𝑜𝑛 − 𝑡0 )¨
𝑟(𝑡0 )
Для данного метода характерны те же недостатки, что и в предыдущем случае — если в момент сближения крайних тел центральное тело
находится достаточно близко к их центру масс, то пренебрежение его влиянием может существенно изменить орбиты тел. Поэтому для применения
этого метода также необходимо отслеживать параметр 𝑅/𝑟, и если его значение в момент 𝑟 = 𝑟𝑐𝑟𝑖𝑡 меньше некоторого значения, применение данного
метода является некорректным. Для интегрирования подобных орбит необходимо использовать другой метод борьбы с сингулярностью.
20
2.1.3
Регуляризация
В случае регуляризации приближения не используются, вместо этого
в уравнениях движения тел делается переход к новой переменной, после
чего особенность в данном месте пропадает. После прохождения опасной
точки происходит возврат к нормальным переменным. В этой работе использовался метод регуляризации Арсета-Заре [7]. Это адаптация метода
KS-преобразований для задачи трех тел. Подробности метода можно посмотреть, например, в [6]. К сожалению, данная регуляризация неприменима для тройных столкновений, которые представляют собой существенную
особенность.
2.2
Программные средства и детали реализации.
В процессе поиска периодических орбит важную роль играет выбор
программных средств, а также алгоритмы, с помощью которых интегрируются уравнения движения и происходит оценка периодичности орбиты.
Необходимо достигнуть баланса между быстродействием, точностью интегрирования и шагом по параметрам, так, чтобы полученные результаты не
были чересчур приближенными, и при этом работа программы занимала
разумное время. Поскольку чем точнее интегрирование и чем мельче шаг
по параметрам, тем более точными и подробными являются результаты,
в первую очередь необходимо сконцентрироваться на увеличении быстродействия программы.
Используемые в данной работе программные средства можно разделить на три вида:
– приближение упругим отскоком и задачей двух тел моделировались с помощью пакета MATLAB;
– регуляризация задачи выполнялась с помощью программы
TRIPLE;
21
– также в процессе обработки данных использовались иные программные средства.
Рассмотрим их подробнее.
2.2.1
Пакет MATLAB
Выбор пакета MATLAB для работы с задачей был обусловлен
несколькими причинами:
1. Наличие различных встроенных алгоритмов интегрирования, реализованных профессиональными математиками и программистами.
2. Наличие встроенных формул интерполяции решения системы дифференциальных уравнений с сохранением точности.
3. Возможность менять метод интегрирования заменой всего
нескольких строк программы.
4. Возможность использовать параллельные вычисления на СPU без
необходимости изучения нюансов языка.
5. Наличие различных встроенных функций (например, поиск нуля
функции).
Также следует отметить следующие преимущества пакета MATLAB,
которые не были использованы в этой работе, однако могут быть применены для продолжения исследования этой задачи:
1. Возможность использовать параллельные вычисления на GPU без
переписывания программы на другом языке.
2. Возможность при необходимости подключать функции, написанные на C/FORTRAN.
В реализации использовался встроенный алгоритм решения систем
ОДУ — ode113. Это многошаговый метод Адамса–Башфорта–Мултона переменного порядка. Он был выбран, так как, во-первых, рекомендуется в
руководстве по MATLAB для решения задач, жестко критичных к ошибкам, и, во-вторых, показал лучшие результаты, чем встроенный метод
22
Рунге-Кутта 4-го порядка, на модельной задаче. С деталями работы алгоритма можно ознакомиться по ссылке [8].
Первоначальная версия программы, использовавшая интегрирование
уравнений движения с мелким шагом и проверявшая на каждом шаге значение параметра 𝐹 (𝑡) работала неудовлетворительно долго. С помощью
следующих модификаций удалось уменьшить время работы программы более чем в 20 раз:
1. Изменение шага интегрирования в зависимости от значения 𝐹 (𝑡).
При достижении порогового значения шаг интегрирования уменьшался в 100 раз до тех пор, пока значение 𝐹 (𝑡) вновь не становилось больше порогового.
2. Использование параллельных вычислений на четырехядерном
процессоре Intel(R) Core(TM) i5-4460 CPU @ 3.20GHz.
3. Рандомизация списка орбит, вычисляемых на каждом ядре процессора. При обычной параллелизации цикла по одному из параметров 𝑘, 𝜇 ядра процессора заканчивали свою работу не одновременно, продлевая время работы программы на существенный срок.
Для компенсации этого эффекта было найдено следующее решение: список орбит генерировался заранее, затем перемешивался в
случайном порядке и уже после этого вычисление орбит распараллеливалось. Это позволило получить существенный прирост скорости работы программы.
4. Прекращение интегрирования орбиты, если ее просчет занимал более 30 секунд (интегрирование обычной орбиты занимает 1-3 секунды).
Для проверки корректности интегрирования орбиты отслеживалось
значение интеграла энергии. Для большинства орбит отклонение не превышало 10−8 , однако присутствовали и орбиты с существенными отклонениями от начального значения интеграла энергии. Они отбрасывались как
некорректные.
23
2.2.2
Программа TRIPLE
Данная программа предназначена для интегрирования задачи трех
тел. Она написана и поддерживается астрономом Sverre Aarseth, научным
сотрудником Института Астрономии в Кембриджском Университете. Использование этой программы достаточно просто, однако в процессе ее эксплуатации требовалось решить несколько проблем:
– Проверить, корректно ли программа TRIPLE работает для плоской
равнобедренной задачи трех тел.
– Добавить код, позволяющий искать периодические орбиты без вмешательства в основную часть программы.
Данные проблемы были решены следующими способами:
1. Был введен параметр
∆𝑠𝑦𝑚 =
√︀
(𝑥1 + 𝑥2 )2 + (𝑦1 − 𝑦2 )2
для отслеживания нарушения симметрии орбит. Было обнаружено, что для некоторых орбит нарушение симметрии является существенным, что невозможно в условиях данной задачи и свидетельствует о некорректности построенной орбиты. Орбиты с
∆𝑠𝑦𝑚 > 10−4 и малым 𝐹𝑚𝑖𝑛 исследовались дополнительно, чтобы выяснить значение ∆𝑠𝑦𝑚 в момент достижения 𝐹𝑚𝑖𝑛 . Если оно
оказывалось больше 10−4 , то орбита отбрасывалась как некорректная.
2. Помимо этого, отслеживалось значение интеграла энергии. Для
большинства орбит отклонение не превышало 10−8 , однако присутствовали также орбиты с существенными отклонениями от начального значения интеграла энергии. Они также отбрасывались
как некорректные.
3. Шаг интегрирования программы при точности интегрирования
10−14 достаточно велик, в случае же перехода к 10−15 время работы программы скачкообразно возрастает. Однако для поиска 𝐹𝑚𝑖𝑛
24
необходим мелкий шаг прохода по орбите. Для того, чтобы преодолеть это затруднение, использовалась квадратичная интерполяция по трем точкам. Для каждых трех опорных точек орбиты,
вычисленных программой, строилась квадратичная интерполяция
˙ Далее на данной параболе бра(параболой) параметров 𝑟, 𝑅, 𝑟,
˙ 𝑅.
лось 1000 промежуточных точек и для каждой из них вычислялось
значение 𝐹 (𝑡).
Помимо этого, некоторые орбиты интегрировались на протяжении длительного времени. Интегрирование таких орбит прерывалось вручную, так
как их количество было мало.
Для ускорения работы программы так же, как и в предыдущем случае, использовались параллельные вычисления. Программа выполнялась
на 24-ядерном процессоре Intel Xeon CPU 2.00GHz.
2.2.3
Дополнительные программные средства
Необработанные файлы выходных данных занимают много строк,
вплоть до десятков гигабайт. Обработать подобные данные вручную, разумеется, невозможно. Таким образом, были необходимы программные средства, позволяющие отсеивать орбиты по некоторым критериям. Для этого
был написан ряд скриптов на языке Perl. Для визуализации результатов
использовался пакет GNUPLOT.
25
Глава 3. Результаты
3.1
Приближение задачей двух тел
Основным ограничением на подробность сканирования является процессорное время. Для оптимизации поиска работа проходила в несколько
этапов.
Первоначально была рассмотрена область 𝑘 ∈ [0.1, 0.9], 𝜇 ∈ [−0.9, 0.9]
с шагом 0.01 по каждой из осей. Несмотря на довольно крупный шаг,
это дало возможность выделить перспективные области для поиска периодических орбит (Рис. 3.1). Было выбрано весьма слабое ограничение на
𝐹 : как близкие к периодическим рассматривались орбиты с 𝐹 < 0.002
(Рис. 3.2). Для сравнения укажем, что 𝐹 для известной орбиты Брука
(𝑘 = 0.418, 𝜇 = 0) составило 0.0006.
Рисунок 3.1 — Полная карта области.
26
Рисунок 3.2 — Точки первого сканирования с 𝐹 ≤ 0.002.
Поскольку данное сканирование заняло более 11 часов, логично предположить, что сканирование с шагом 0.001 займет в 100 раз больше времени, то есть примерно полтора месяца. Так как данный срок не представлялся разумным, для подробного сканирования были выбраны наиболее
перспективные области:
1. 𝑘 ∈ [0.520, 0.700], 𝜇 ∈ [0.410, 0.800],
2. 𝑘 ∈ [0.400, 0.480], 𝜇 ∈ [0.220, 0.460],
3. 𝑘 ∈ [0.520, 0.700], 𝜇 ∈ [−0.800, −0.410],
4. 𝑘 ∈ [0.400, 0.480], 𝜇 ∈ [−0.460, −0.220],
5. 𝑘 ∈ [0.400, 0.480], 𝜇 ∈ [−0.220, 0.220].
Получившийся результат можно видеть на Рис. 3.3.
Начальные условия для орбит, удовлетворяющих условию 𝐹 ≤
0.0006, показаны на Рис. 3.4.
При этом было найдено довольно много орбит с соответствующими
значениями 𝐹 (Таб. 3.1).
27
Рисунок 3.3 — Области, просканированные с шагом 0.001.
Рисунок 3.4 — Начальные условия для орбит, близких к
периодическим(𝐹 ≤ 0.0006).
28
Таблица 3.1
Количество найденных орбит в зависимости от 𝐹 .
𝐹 (·10−4 ) 5 - 6 4 - 5 3 - 4 2 - 3 1 - 2 < 1 Всего
N
40
43
36
21
16
4
160
3.2
Приближение упругим отскоком
Последовательность действий, предпринятых в данном пункте, аналогична предыдущему. Первоначально была рассмотрена область 𝑘 ∈
[0.1, 0.9], 𝜇 ∈ [−0.9, 0.9] с шагом 0.01 по каждой из осей (Рис. 3.5)
Рисунок 3.5 — Полная карта области.
29
Было выбрано то же, что и в предыдущем случае, слабое ограничение
на 𝐹 : как близкие к периодическим рассматривались орбиты с 𝐹 < 0.002.
Как можно видеть на Рис. 3.6, перспективные области остались теми же.
Рисунок 3.6 — Точки первого сканирования с 𝐹 ≤ 0.002.
Для подробного сканирования были выбраны следующие области:
1. 𝑘 ∈ [0.520, 0.700], 𝜇 ∈ [0.410, 0.800],
2. 𝑘 ∈ [0.400, 0.480], 𝜇 ∈ [0.220, 0.460],
3. 𝑘 ∈ [0.520, 0.700], 𝜇 ∈ [−0.800, −0.410],
4. 𝑘 ∈ [0.400, 0.480], 𝜇 ∈ [−0.460, −0.220],
5. 𝑘 ∈ [0.400, 0.480], 𝜇 ∈ [−0.220, 0.220].
Получившийся результат можно видеть на Рис. 3.7.
Орбиты, удовлетворяющие условию 𝐹 ≤ 0.0006, показаны на Рис. 3.8.
Орбит с соответствующими значениями 𝐹 было найдено довольно
много (Таб. 3.2).
30
Рисунок 3.7 — Области, просканированные с шагом 0.001.
Рисунок 3.8 — Орбиты, близкие к периодическим(𝐹 ≤ 0.0006).
31
Таблица 3.2
Количество найденных орбит в зависимости от 𝐹 .
𝐹 (·10−4 ) 5 - 6 4 - 5 3 - 4 2 - 3 1 - 2 < 1 Всего
N
68
44
52
33
17
9
223
Следует отметить, что результаты, полученные этим методом, отлично согласуются с результатом, полученым методом приближения задачей
двух тел, что можно видеть при сравнении рисунков 3.4 и 3.8.
Для примера приведем несколько построеных орбит (Рис. 3.9.). Левая колонка соответствует методу приближения задачей двух тел, правая
- приближению упругим отскоком.
32
3
3
2
2
1
1
0
0
-1
-1
-2
-2
-3
-3
0
0.5
1
1.5
2
2.5
3
3.5
0
0.5
1
1.5
r
2
2.5
3
3.5
2
2.5
3
3.5
2
2.5
3
3.5
r
a)
b)
4
4
3
3
2
2
1
1
0
0
-1
-1
-2
-2
-3
-3
-4
-4
0
0.5
1
1.5
2
2.5
3
3.5
0
0.5
1
1.5
r
r
c)
d)
4
4
3
3
2
2
1
1
0
0
-1
-1
-2
-2
-3
-3
-4
-4
0
0.5
1
1.5
2
r
2.5
3
3.5
0
0.5
1
1.5
r
e)
f)
Рисунок 3.9 — Орбиты: a) 𝑘 = 0.554 𝜇 = −0.560 𝐹 = 0.000119,
b) 𝑘 = 0.552 𝜇 = −547 𝐹 = 0.000574,
c) 𝑘 = 0.527 𝜇 = −0.421 𝐹 = 0.000141,
d) 𝑘 = 0.527 𝜇 = −0.426 𝐹 = 0.000127,
e) 𝑘 = 0.457 𝜇 = 0.458 𝐹 = 0.000038,
f) 𝑘 = 0.455 𝜇 = 415 𝐹 = 0.000469.
Слева - приближение задачей двух тел, справа - упругим отскоком.
33
3.3
Регуляризация.
Исследование орбит с помощью регуляризации проводилось программой TRIPLE. С одной стороны, поскольку программа написана на языке
FORTRAN, она отличается высокой скоростью работы. С другой стороны,
поскольку для вычисления 𝐹𝑚𝑖𝑛 применялася интерполяция квадратичным
полиномом, точность поиска орбит, близких к периодическим, уменьшается.
При проведении оценочного сканирования были получены следующие результаты, представленные на рис. 3.10. Как можно видеть, отчет-
Рисунок 3.10 — Полная карта области. 𝐹 ≤ 1.
ливо выделяется 3 полумесяца. Самый крупный из них соответствует так
называемой орбите Брука (𝑘 = 0.418, 𝜇 = 0). Логично предположить, что
оставшиеся 2 полумесяца также соответствуют периодическим орбитам.
34
Также именно в центре этих областей наблюдается концентрация точек с
𝐹 < 0.006 (рис. 3.11)
Рисунок 3.11 — Полная карта области. 𝐹 ≤ 0.006.
Первоначально было выполнено сканирование только перспективных
областей. Однако в дальнейшем удалось получить доступ к более мощному
компьютеру и с его помощью просканировать всю область с шагом 0.001
по 𝑘 и 𝜇. Сканирование заняло более недели. Результаты представлены на
рисунке 3.12.
Нужно отметить, что отчетливо проявился еще один полумесяц, как
бы "обнимающий"полумесяц, порожденный орбитой Брука.
Параметр 𝐹𝑚𝑖𝑛 для орбиты Брука равен 5 * 10−4 . Поскольку точность
получения 𝐹𝑚𝑖𝑛 в этой реализации меньше, чем в предыдущей, как близкие
к периодическим были обозначены орбиты с 𝐹𝑚𝑖𝑛 ≤ 10−3 . Их распределение
можно видеть на рисунке 3.13.
Основную долю близких к периодическим орбит составляют орбиты, порожденные орбитой Брука. Однако также присутствуют и орбиты,
35
Рисунок 3.12 — Полная карта области. 𝐹 ≤ 0.01.
как принадлежащие полумесяцам, так и внешне случайно распределенные.
Примеры орбит можно видеть на рисунке 3.14.
В таблице 3.3 представлено найденное число орбит с соответствующими значениями 𝐹𝑚𝑖𝑛 .
Таблица 3.3
Количество найденных орбит в зависимости от 𝐹 .
𝐹 (·10−4 ) 9 - 10 8 - 9 7 - 8 6 - 7
N
921
730 602 474
4-5
203
3-4
119
2-3 1-2
78
33
<1
3
5-6
330
Всего
3493
36
Рисунок 3.13 — Полная карта области. 𝐹 ≤ 0.01.
37
a)
b)
c)
d)
e)
f)
Рисунок 3.14 — Орбиты в координатах 𝑟, 𝑅 (левый столбец) и
относительно центра масс (правый столбец): a), b) 𝑘 = 0.418, 𝜇 = 0.0,
𝐹 = 5 * 10−4 , c), d) 𝑘 = 0.7356, 𝜇 = −0.2968, 𝐹 = 5 * 10−4 , e), f)
𝑘 = 0.7349, 𝜇 = 0.3061, 𝐹 = 9 * 10−4
В правом столбце разные цвета показывают орбиты разных тел.
38
Заключение
В дипломной работе можно выделить следующие пункты:
1. Были рассмотрены три метода обработки двойных столкновений в
плоской равнобедренной задаче трех тел:
– приближение задачей двух тел;
– приближение упругим отскоком;
– регуляризация (с помощью программы TRIPLE).
Все три метода показали хорошо согласующиеся результаты, что
позволяет говорить о корректности исследования задачи с их помощью. При этом границы применимости методов различны, благодаря чему орбиты, которые некорректно обрабатываются одним
методом, можно исследовать с помощью другого.
2. Был предложен метод оценки близости орбиты к периодической с
помощью параметра 𝐹𝑚𝑖𝑛 .
3. Было проведено сканирование области начальных условий, завдаваемой параметрами 𝑘, 𝜇 тремя различными методами. Это позволило выделить четыре области, где, предположительно, находятся
периодические орбиты, а также найти большое количество орбит,
близких к периодическим.
Можно выделить следующие перспективные направления для дальнейшего исследования задачи:
1. Провести более детальное сканирование наиболее перспективных
областей с большей точностью. Для этого потребуются значительные вычислительные мощности.
2. Провести сканирование для более длительного промежутка времени 𝑇 , что, вероятно, позволит обнаружить орбиты с периодом,
большим 100𝜏 .
3. Провести сканирование для других значений 𝜀, рассмотреть эволюцию областей орбит, близких к периодическим, при изменении
соотношения масс тел.
39
4. Для найденных орбит, близких к периодическим, рассмотреть их
эволюцию при отклонении начальных условий от плоской равнобедренной задачи трех тел, что, возможно, позволит найти новые
семейства периодичских орбит в общей задаче трех тел.
40
Список литературы
1. L. Euler. De moto rectilineo trium corporum se mutuo attrahentium. //
Novi Comm. Acad. Sci. Imp. Petrop. — 1767.
2. Orlov V.V., Petrova A.V., Martynova A.I. Classification of orbits in the
plane isosceles three-body problem // Mon. Not. R. Astron. Soc. — 2002.
— Vol. 333. — Pp. 495–500.
3. Аносова Ж.П., Орлов В.В. Динамическая эволюция тройных систем. //
Труды АО ЛГУ. — 1985. — Т. 40. — С. 66–144.
4. Broucke R. On the Isosceles Triangle Configuration in the Planar General
Three-body Problem. // Astron. Astrophys. — 1979. — Vol. 73. — Pp. 303–
313.
5. Орлов В.В., Петрова А.В., Мартынова А.И. Тройные сближения в
плоской равнобедренной задаче трех тел равных масс. // Письма в Астрономический журнал. — 2001. — Т. 27. — С. 795–800.
6. Динамика тройных систем. / А.И. Мартынова, В.В. Орлов, А.В. Рубинов и др. — СПб: Изд-во С.-Петерб. ун-та, 2010. — 216 с.
7. Aarseth S.J., Zare K. // Celestial Mechanics. — 1974. — Vol. 10. — P. 185.
8. MathWorks. — 2016. http://www.mathworks.com/help/matlab/ref/ode113.
html.
Отзывы:
Авторизуйтесь, чтобы оставить отзыв