Вычислительные технологии
Том 17, № 4, 2012
Символьно-численные расчёты нормальной формы,
интеграла движения и сечений Пуанкаре для консервативных двумерных гамильтоновых систем
В.Е. БоглчЁв, Н.А. Чеканов Белгородский государственный национальный исследовательский университет, Россия e-mail: [email protected], [email protected]
Представлено краткое описание разработанных в среде Maple 8 алгоритмов для символьно-численных расчётов нормальных форм и формальных (приближённых) интегралов движения по методу Биркгофа — Густавсона для гамильтоновых систем с произвольным числом степеней свободы. Приведены результаты расчётов сечений Пуанкаре, полученные при данном подходе с помощью приближённого интеграла движения, и их сравнение с точными вычислениями в случае гамильтоновых систем с двумя степенями свободы. Показано, что имеется хорошее согласие результатов расчётов в области энергий, где классическое движение является регулярным.
Ключевые слова: классическая механика, функция Гамильтона, нормализация, сечения Пуанкаре, динамический хаос, интегралы движения, математическое моделирование.
Введение
Поскольку большинство динамических задач классической механики не могут быть представлены в явном виде (см., например, [1-3]), то продолжают разрабатываться как различные приближённые аналитические методы [4-8], так и методы непосредственных прямых численных расчётов [9-11]. На сегодня перспективным направлением является создание гибридных, или комбинированных, методов, в которых вначале производятся аналитические вычисления, затем — численные расчёты с использованием современных математических пакетов прикладных программ, например MAPLE, REDUCE, MATHEMATICA и др. [12-14].
В настоящей работе представлены результаты символьно-численных расчётов нормальной формы, интеграла движения и сечений Пуанкаре для гамильтоновых систем с двумя степенями свободы. Проведено сравнение сечений Пуанкаре, полученных с помощью интегралов движения по методу Биркгофа — Густавсона, с их прямыми численными расчётами. Сечения Пуанкаре для неинтегрируемых гамильтоновых систем, вычисленные с помощью полученных интегралов движения, достаточно хорошо воспроизводят результаты прямых численных расчётов для энергий, при которых доля фазового пространства, занятого хаотическими траекториями, небольшая. Для энергий выше критических при сильно развитом динамическом хаосе согласие не имеет место.
1. Метод Биркгофа — Густавсона
Пусть дана классическая система с п степенями свободы, гамильтониан которой можно записать в виде разложения
п те
н = £ у (рк + о2) + £ £ С 1,...,1п,т 1,...,тп р1 р22 .. .рП • о?1 о??2 ...от-, (1)
к=1 5=3 11+...+1- +
+т 1 +...+тп=Я
где р = (р1,р2,... ,рп), о = (?1, ?2,... , 0п) — канонически сопряжённые переменные, С11,..,1п,т1,..,т„ — известные коэффициенты. Частоты Шк могут быть как несоизмеримыми (нерезонансными), так и соизмеримыми (резонансными), для которых имеется г независимых соотношений
к=п
^ аши = 0, г = 1, 2, 3 ..., г, (2)
к=1
где матричные элементы а^ — целые числа. Для приведения исходного гамильтониана (1) к нормальной форме необходимо выполнить весьма трудоёмкие канонические преобразования (о,'р) ^ (С,п). Как известно [15, 16], функция Гамильтона Н(о,'р) находится в нормальной форме, если выполняется условие
БСа,п) = 0, (За)
п ( д д \ О = £ шк ^- й дп)- РЧ
О — так называемый дифференциальный оператор нормальной формы.
Как показано в работах [15, 16], для выполнения процедуры нормализации необходимо решить основное уравнение
О(о,п) Ж(д,п) = -Н(5)(д,п) + Г(5)(д,п), (4)
где Ж(5) — однородный полином степени Б в производящей функции
Р (д,п) = 0 • П + Ж (3)(д,п) + ••• + Ж (5)(д,п) + ••• , (5)
которая приводит к нормальной форме полиномиальные члены исходного гамильтониана (1) степени Б. Слагаемые функции Гамильтона более высокой чем Б степени вычисляются по формуле
г<'>^ = нО(,,) + £1 - (^)). (6)
|к| ^ /
где г = Б + 1, Б + 2,... и I - |к| + |к| (Б - 1) = г, 1 < |к| < I < г, I > 2, в > 3,
|к| = N + к| + ... + |кп|, к! = к1! • к^. • ... • кп!.
В работе [15] доказано, что существует (п - г) приближённых независимых интегралов движения вида
п
I (е,п) = £ | (£ + п2), (7)
г=1
где ^ = ... есть любой п-мерный вещественный вектор такой, что выпол-
няются уравнения
п
^«¿к^к = 0,г = 1, 2, 3,..., г. (8)
к=1
Выполняя обратное преобразование (£,п) ^ (?,р), можно получить интегралы движения в исходных канонически сопряжённых переменных р = (р1,р2,... ,рп),
2. Краткое описание алгоритма
Input:
n — число степеней свободы;
SMAX — желаемый максимальный порядок нормальной формы и интегралов движения;
jmax — максимальная степень исходного гамильтониана, квадратичная часть которого задается в нормальной форме, а v[2] + v[3] + ... +v[j] + ... + v[jmax] — его потенциальная часть;
если gip = 0, то нормальная форма вычисляется в переменных действие — угол; если intg2 = 0, то вычисляются интегралы движения, независимые от интеграла энергии;
если intg2 = 0 и test1 = 0, то скобка Пуассона pb = {H, intg2} не должна равняться нулю;
если intg2 = 0 и test2 = 0, то Htest2 должен быть равен исходному гамильтониану. Output:
gSMAX = g(2) + g(3) + ... + g(SMAX) — нормальная форма;
wSMAX = w(3) + w(4) + ... + w(SMAX) — производящая функция, где — однородные полиномы степени S вида pi1 pl22 ...Pn ■ oT1 oT2 ... On", S = li + l2 + ... + ln + mi + m2 + ... + mn;
intg[i] — дополнительные интегралы движения в переменных (pi1 ,p22, ...,рП),
(oT1 ,oT2,..., oTn).
3. Результаты расчётов для некоторых гамильтоновых систем
Ниже представлены нормальная форма С^млх и приближённые интегралы движения I, вычисленные с помощью разработанной программы, и сечения Пуанкаре, полученные как с помощью приближённого интеграла движения, так и прямыми численными расчётами для трёх консервативных гамильтоновых систем с двумя степенями свободы. 1. Неинтегрируемая система с функцией Гамильтона
н = 2(р2 + р2) + ^2 + ) + + с(?2 + ?2)2, (9)
для которой вычислены нормальная форма и приближённый интеграл движения до степени Б МАХ = 6. Нормальная форма здесь и ниже представлена также в переменных
действие — угол (11,/2, ф 1, ф2):
С(/1, /2, ф2) = 2/1 + (6 + 5с)/? - (^ /22 - (^ + 21Г + 25с2 ) 13+
+ Ч№ + 21Ъс - 9 Л /1 /2 + /6 + (1 - 26/1 - 9с/1)(4/2 - /2)008(2^), (10)
А 8 2 У 1 2 V 8
С6 = С(2) + С(3) + С(4) + С(5) + С(6), (11)
С(2) = 1(п2 + е2) + 2(п2 + £), (12)
С(4) = 3С ((£ + ^ + (п2 + П2)2) + ((6^2 - С2П1)2 + 36 && + П1П2)2 +
3с
+:4(6п1 + 6П2)2, (13)
С(6) = ^ ((£ + е2)3 + (п2 + п2)3) - ^ + С2 + п2 + п2) ((6 + 326с + 60с2)х
x(^1П2 - 6ш)2 + 176(6 + 8с)(66 + П1П2)2 + 204с2(С1П1 + 6П2)2) , (14)
С(3) = 0, С(5) = 0, (15)
3 ( )
/ = 3 с (2р1о2 + р1 + р2 + 2р2о2 + «4 + о4) +
с2 ( )
+ 32 (-35р6 + 13о6 + 1306 - 35р2 - 105р4о2 - 57р2о4 - 105р4о2 - 57р2«4) +
^ 1 » ' ,2 2, о 2 2 , „2 2, о „22
+ (86 + 1 с) (о2р2 + 3«2«2 + р202 + 3р1Р2 + 4Р101Р202) +
+ (- 256С - Тс2 - 8(Р301Р202 + Р1 Р20201) + + (- 176с - 862 - 15(р103р202 + р1 Р2?201) +
+(-5 6с +15862 - 1 с2) (р1«22+«эд+(162+| *+32 с2) ш++ +(-И62 -1 "с -1 с2) (^+р?р4)+(-12862 - 8 "с+32 с2) («м+р?о2)+
+(- 6462 - 76с2 - т6 6с) (р2,2р2+р2р2,2)+
+ (-6562 - 26'с - 26с2) (Р?«2«2 + А&2) • (16)
Как видно из рис. 1, а сечения Пуанкаре, вычисленные с помощью приближённого интеграла (3), достаточно хорошо воспроизводят структуру фазового пространства для энергий, при которых в системе реализуется регулярный характер движения. По мере
+ + + нк
/ + + 1 ++ +
V ................/
а
Рис. 1. Сечения Пуанкаре, полученные прямым численным решением уравнений движения (слева) и при помощи интеграла движения (справа); а — полная энергия Е = 0.8, начальные данные: Ь = 0.5, с = 0.01, (р ,р2,91,92) = (0, 0.1, 0, 0), (0, 0.3, 0, 0), (0, 0.4, 0, 0), (0, 0.5,0,0), (0,0.7, 0,0), (0,0,0, -0.812), значение интеграла I = 0.01,0.02,0.03,0.0396,0.045; б — полная энергия Е = 2.5, начальные данные: Ь = 0.5, с = 0.01, (р1 ,р2,91,92) = (0, 0.2,0,0), (0, 0.5,0,0), (0,0.6,0, 0), (0, -0.45,0,-1.85), (0,1,0, 0), (0,1.2, 0, 0), (0,0,0,-1.9), (0,1.38, 0, 0), значение интеграла I = -1.1, -0.5, -0.2, 0.0123, 0.16, 0.3, 0.4; в — полная энергия Е = 2.5, начальные данные: Ь = 0.5, с = 0.01 (р1 ,р2, 91,92) = (0, 0.6, 0, 0.3), (0, 0, 0, -4), (0, 3.6,0, -4.5), (0,1.6, 0, 0), (0, -2.6,0, -4.9), (0, -3.4,0, 2.55), (0, -2.1,0, -3), значение интеграла I = -2200, -1500, -1000, -500, -125, -6
увеличения энергии в системе развивается динамический хаос, что приводит к разрушению приближённого интеграла движения. В этом случае формальный интеграл движения (3) не приводит к истинным сечениям Пуанкаре (рис. 1, б, в).
2. Система с функцией Гамильтона
Н = 2(р?+р2) + V («1 ,«2),
(17)
где V(«1 ,«2) = 2(«1 + «2) + 6 «2 + 3 «23^ +4 с«12«22 + с («12 + «22)2, является интегрируемой. В самом деле, замена переменных х = «1 - «2, у = «1 + «2, как легко проверить, приводит к разделению переменных. Поверхность потенциальной энергии V(«1, «2) при значениях параметров 6 = 1/4, с = 1/288 имеет пять локальных минимумов и четыре седловых точки, изолинии которой показаны на рис. 2. Были вычислены нормальные формы, в том числе и в переменных угол — действие до степени Б МАХ = 8:
ОД, /2, ф2) = 2/1 + ( - 562 + 9с) /2 + (2462 - 4с) /2 + (13564 + |с2 - Т5с62) /2/! +
1175,4 ,о Л т3 /1125 3
+ 1 —ттгтт64 + 125с62 - 85с2 /3 + ( ^гс3 +
108
V 128
12925 - 6^ с2 6" 1 Л+
1536
82944
256
+
/192925 V 10368
66-
5625
16
с3 +
41575 с262 - 64625с6М /22/? +
32
192
+
291025 2,2 452375 ,4 13125 3 1350475 Л 4 / 5 Л . ,
«-с2 62 + ^с64 + —с3 - 1шт 67/4 + Г - 5 008 (2ф2 >/2+
+ (й62 - 3с) 008 (2ф2) /22 + (-51с2 + 75с62 - 2365ьА 008 (2ф2) /2+
Рис. 2. Изолинии поверхности потенциальной энергии (17) при Ь = 1/4, с = 1/288
4144Ь4 + ТС2 - ТС&2) С°8(2Ф2) 1211 +
(8315 2,2 38585 6 12925 ,4 „ 3) , т2
с2ь2 + з8^б6 - — сб4 - 375С^ С°8 (ЗД 121!2+
( 375 3 8315 2г2 12925 г4 38585 Л + (-16"с3 + "96Т^ - -ш СЬ4 + ЙШЬ ) С°8 (4ф2) 1212+
+ (-58205 С2Ь2 + 90175 СЬ4 - 27°095 ьб + 2625 ¿Л _ ш 14+
V 12 72 3888 2 ) у У2> 1
( 8315 2г2 38585,6 12925 г4 375 Л , г4
+--с2Ь2--Ь6 +--сЬ4 +--с3 с°в (4ф2) 14+
V 48 15552 288 8 ) У Г2> 1
(12925 г4 38585,6 8315 2,2 375 3) , г4
+ -сЬ4--Ь6--с2Ь2 +--с3 с°в (2ф2) 124+
V 1152 62208" 192 " " ' 32 ~ ' 4 Т2у 2
I 38585 ,6 12925 ,4 8315 2,2 375 3 ^ ,
+--Ь6 +-сЬ4--с2Ь2 +-с3 с°в (4ф2) 124, (18)
1 248832 4608 768 128 ' У Г2) 2' ^ '
С(4)
С8 = С(2) + С(3) + С(4) + С(5) + С(6) + С(7) + С(8), (19)
С(2) = 1(п2 + б2) + 2(П2 + £), (20)
(8 с - 48 Ь^ (^ + ^22)2 + (П12 + П22)2 + 4 (66 + П1П2)2 +
+2 (б2+б2 )(П12 + П22)) , (21)
^ = (! Ь2с - 96 с2 - Ц *) (^ + ^22)3 + (П12 + П22)3 +
+3 (б2 + б2 + П12 + П22) (4 (66 + П1П2)2 + (б2 + е22) (П12 + П22))) , (22)
Г(8) =( 38585 Ь6 8315 с2Ь2 + 12925 сЬ4 + 375 Л ((е 2 + е 2)4 + ( 2 + „ 2)4 +
с = Г 497664 Ь - ^ сЬ +^216 сЬ +256 сД ^ + 6 ^ + ^ + П2 ^ + + 16 + П1П2)4 + 2 (б2 + б2) (П12 + П22) (3 (б2 + б2) (П12 + П22) + 2 (б2 + б2)2 + +2 (т2 + П22)2) + 24 (66 + П1П2)2 (С12 + б2 + П12 + П22)2) , (23
С(3) = 0, С(5) = 0, С(7) = 0. (24)
Ниже приведён интеграл движения до степени БМАХ = 6 (интеграл до БМАХ = 8 не приводится из-за громоздкости, хотя на его основе были проведены расчёты сечений Пуанкаре) (рис. 3):
35 ( ) (
I = - 3456 (ИЬ4 - 108Ь2с + 108с2) (р6 + 15р2р4 + 15р1р2 + р? + 3?|р4 + 18д|р?р2 + 3^4+
+24?152Р1Р3 + 24^2р1р2 + 3?2р2 + 18?2р2р2 + 3?2р1) - т1^ (2052с2 + 385Ь4 - 3300Ь2с) +
1152
+ (^Р2 + + 8?1?3Р1Р2 + 6?252Р2 + б^р! + 853?2Р1Р2 + ?4р2 + ^Р2) -
■0.4 -0.2 0 0.2 0.4 Я2 -0.4 -0.2 0 0.2 0.4 д.
Рис. 3. Сечения Пуанкаре, полученные прямым численным решением уравнений движения (слева) и при помощи интеграла движения (справа). Полная энергия Е = 0.1, начальные данные: Ь = 1/4, с = 1/288, (р1 ,р2, 91,92) = (0,0.07,0,0), (0, 0.12, 0, 0), (0, 0.16, 0, 0), (0, 0.2, 0, 0), (0, 0.26, 0, 0), (0, 0.3, 0, 0), значение интеграла I = -0.00042, -0.0004, -0.00035, -0.0003, -0.000265, -0.00023
1 (-2916&2с + 545Ь4 - 1404с2) (д* + д?) (д24 + 14д?д| + д4) - -1- (5^ - 18с) +
3456 ч ------■ - — -----, V*2 ■ ■ 144
+ (3р2 + 18р2 р2 + 3р4 + 6д2 р2 + 6д| + З94 + 24д192Р1Р2 + бд2^2 + бд2^2 + 18д2 д2+ +3д 4 + 4Ьд3р2 + 4Ьд3 Р? + 4Ьд5 + 24Ьд1 д^^ + 12Ьд 2 д2 р2 + 12Ьд 2 ^2 р2 + 40Ьд 2 д3+
+8Ьд?р1Р2 + 20Ьд4д2) . (25)
При сложной топологии поверхности потенциальной энергии, например многих минимумов и седловых точек, согласие с точными сечениями Пуанкаре и вычислениями, полученными с помощью интеграла движения, сильно ухудшается. Это показывает, что в таких системах метод нормальных форм даёт удовлетворительное согласие с истинным движением систем только для относительно небольших энергий по сравнению со значением энергии в седловой точке.
3. Функция Гамильтона с сингулярностью в начале координат.
Метод нормальных форм Биркгофа — Густавсона вместе с вычислением приближённого интеграла движения можно применить к динамическим системам, которые содержат сингулярность в начале координат.
К примеру, рассмотрим так называемую плоскую ограниченную задачу трёх тел. В этой задаче движение трёх тел происходит таким образом, что массой одного из них и его влиянием на движение двух других пренебрегается, в то время как массивные тела взаимодействуют друг с другом и с этим телом. В подвижной системе координат, связанной с телом наибольшей массы, гамильтониан в безразмерных переменных имеет вид [17]
H = й + ^ + - хру - X2 + У— - , 1 , (26)
(здесь рх,ру и х, у — канонически сопряжённые импульсы и координаты), содержащий особенность в начале системы координат (х = 0,у = 0), а величина Н(х,у,рх,ру) = Е есть интеграл движения.
Для устранения сингулярности произведём каноническое преобразование с помощью формул
х = «2 - «2 , У = 2«1«2,
_ «1р1 - «2р2 = «2р1 - «1р2 (27)
рх 2(«2 + «2), ру 2(«2 + «2), (27)
а также вместо независимой переменной времени £ введём новый независимый параметр т согласно формуле
= 4 («2 + «2) ¿т. (28)
В результате получим соотношение
Н - е = 2(р1 + р2) - 2(«2р1 - «1р2)(«2 + «2)--(«2 - «2)2 («2 + «2) - 4Е («2 + «2) - 4 = 0. (29)
Таким образом, финитное, т. е. ограниченное, движение происходит в так называемых областях Хилла, где величина Е отрицательна.
Рассматривая только ограниченное движение (Е < 0), введём новую положительную величину ш = V-8Е, а также новый гамильтониан К исходя из выражения (3): 1 ш2
к = 2(р2 + р2) + -у («2 + «2) + 2(«2р1 - «1 р2)(«2 + «2) + 4(«2 - «2 )2(«2 + «2). (30)
При этом движение системы с новым гамильтонианом происходит с условием, что его значение равно постоянной величине К = 4. Произведя каноническую замену «к ^ «к/л/ш, рк ^ ш • рк (к =1, 2), вводим гамильтониан
1 2 4 ' ° ° ° °ч / \/ 2 2\
ovt-1 ■ х-2 ■ ^ ■ ■ "а (?2Р1 - qiP2)(qi +
2 Ш Ш
k = -(p + р2 + q2 + q2) + (q2Pi - qiP2)(q2 + q^ + -т(q2 - q^^2 + q^, (31)
значение которого равно 4/ш.
Ниже представлены нормальная форма, в том числе и в переменных угол — действие, и интеграл движения до степени SMAX = 6:
4 11 4 I-
G(/l, /2, Ф2) = /1 + -7/3 - -Гh/2 + -2J (4/i2 - h^ Ii sin (^2) +
Ш Ш Ш *
o
+ — (4/i2 - /22) /i cos (2 Ф2) , (32)
Ш4 4 '
G6 = G(2) + G(3) + G(4) + G(5) + G(6), (33)
g(2) = -(n2 + £ ) + 2(n2 + C2), (34)
G(4) = -2 (n2 + n2 + C2 + £2)(ni 6 - П2 Ы, (35)
2
Рис. 4. Сечения Пуанкаре, полученные прямым численным решением уравнений движения (слева) и при помощи интеграла движения (справа). Полная энергия Е = -20, начальные данные: (рх,р2,91,92) = (0,0,0,0), (0,0,0,0.2), (0,0,0,0.3), (0,0,0,0.4), (0,0,0,0.5), (0,0.6,0,0.1), значение интеграла I = -0.00001,0.0002,0.0004,0.0006,0.0008,0.001, 0.0012
G(6) =
-5
4W4
12 + П12)3 + (б2 + П22Г - 3 (П22 + П12 + б2 + б2) 6 + П1П2)2-
- (6П2 - 6П1)2)) ,
2 )3
G(3) = 0, G(5) = 0,
(36)
(37)
1
I = -о (<?2 Р1 + Р2 P192 + 91 92P1 + P1 92 - 91 P2 - P2 91 - P1 91P2 - 92 91P2) +
w2
1
/ ч /I V /1 ) V Л ) Л ) А )
Ь—-г(48Р2 Р19291 + 7 91 Р2 +792 р + 16 92 91Р2Р1 - 9Р1 92 + 15 91 92 + 15р Р1 -4 ш4
42 42 3 3 42 222 222
9 Р2 91 + 15 Р1 Р2 + 48 р1 9291Р2 + 16 91 92Р2Р1 - 15 Р2 92 +6 р2 Р1 92 + 22 92 91 р2 + +6Р22Р12912 + 22 92291V - 5Р16 - 5 916 - 5 Р26 - 5 926 + 15 912924 - 15 92 V-
-15Р14912 - 15 914Р12). (38)
На рис. 4 приведены сечения Пуанкаре, полученные в обоих подходах.
Заключение
Созданные программы в символьно-численном виде выполняют нормализацию и вычисления интегралов движения гамильтоновых систем, в общем, с произвольным числом степеней свободы с применением математического пакета MAPLE 8. В случае га-мильтоновых систем с двумя степенями свободы разработанные программы позволяют строить сечения Пуанкаре. При этом структура сечений Пуанкаре хорошо согласуется с их непосредственным прямым численным расчётом только до энергий, при которых классическое движение является регулярным. Предложенные программы могут быть применены для приближённого решения уравнений движения Гамильтона и использованы при изучение нелинейных гамильтоновых систем.
Список литературы
[1] Арнольд В.И. Математические методы классической механики. М.: Наука, 1974. 432 с.
[2] Переломов А.М. Интегрируемые системы классической механики и алгебры Ли. Ижевск: НИЦ "Регулярная и хаотическая динамика", 2002. 238 с.
[3] Гориэли А. Интегрируемость и сингулярность. М.; Ижевск: НИЦ "Регулярная и хаотическая динамика"; Институт компьютерных исследований, 2006. 316 с.
[4] Джакалья Г.Е.О. Методы теории возмущений для нелинейных систем. М.: Наука, 1979. 320 с.
[5] Боголюбов Н.Н., Митропольский Н.Н. Асимптотические методы в теории нелинейных колебаний. М.: Наука, 1974. 504 с..
[6] Найфэ А. Методы возмущений. М.: Мир, 1976. 456 с.
[7] Найфэ А. Введение в методы возмущений. М.: Мир, 1984. 535 с.
[8] BASIOS V., Chekanov N.A., MARKOVSKI B.L. et al. GITA: A REDUCE program for the normalization of polynomial Hamiltonians // Comp. Phys. Commun. 1995. Vol. 90. P. 355-368.
[9] Березин И.С., Жидков Н.П. Методы вычислений. М.: Гос. изд-во физ.-мат. лит-ры, 1962. 640 с.
[10] Форсайт Дж., Мальколм М., Моулер К. Машинные методы математических вычислений: Пер. с англ. М.: Мир, 1980. 277 с.
[11] Уилкинсон Дж., Райнш К. Справочник алгоритмов на языке АЛГОЛ. Линейная алгебра. М.: Машиностроение, 1976. 392 с.
[12] Дьяконов В.П. Maple 8 в математике, физике и образовании. М.: СОЛОН-Пресс, 2003. 656 с.
[13] Hearn A.C. REDUCE User's Manual. Santa Monrca, CA, USA, 2004. 215 p.
[14] Акритас А. Основы компьютерной алгебры с приложениями. М : Мир, 1994. 544 с.
[15] Gustavson F.G. On TOns^u^ion formal integral of a Hamiltonian system near an equilibrium point // Astronom. J. 1966. Vol. 71, No. 8. P. 670-686.
[16] Мозер Ю. Лекции о гамильтоновых системах: Пер. с англ. М : Мир, 1973. 169 с.
[17] Себехей В. Теория орбит: Ограниченная задача трёх тел. М.: Наука, 1982. 656 с.
Поступила в 'редакцию 13 февраля 2012 г., с доработки — 25 июня 2012 г.