Из рис. 2 б видно, что в неподвижных координатах угловое положение перигелия в плоскости его орбиты определяется дугами GD, DA и AB:
р0a = GD +DA+ AB. (22)
Для определения этих дуг рассмотрим сферический треугольник 2DА, в котором известны два угла: 2 = iEa, А = iеа, и сторона 2А = а - 2, лежащая против неизвестного угла D. Его можно определить по «теореме косинусов» (формула (1.1.010) из [8]):
cosD = -cos iEa·cosieа + sin iEa·sinieа·cos(а - 2). (23)
Так как угол наклона орбиты к плоскости неподвижного экватора iа = - D, то тогда из (23) получаем его значение в виде:
ia = - arccos(-cos iEa·cosiеа + sin iEa·sinieа·cos(а - 2)). (24)
Далее, из треугольника 2DА, по теореме синусов (см. формулу (1.1.007),
[8]) находим величину дуги 2D и положение восходящего узла в виде формул:
= 0D =02 + 2D = 02 + arcsin[sin(а - 2) sin iEa /sinia]. (25)
Для определения положения перигелия на орбите относительно неподвижной точки G (см. рис. 2 б), согласно (22), найдём аналогично дугу DA из 2DA :
sin iEa / sin DA = sin ( - iа) / sin 2A.
Из этого вытекает, что дуга
DA = arcsin[sin iEa·sin(а - 2) / siniа]. (26)
Дуга AB, фигурирующая в выражении (22), определяется по формуле
AB = щ = а - а. (27)
Из прямоугольного сферического треугольника 0GD, в котором известны дуга 0D = и угол D = iа, определяем дугу 0G согласно формуле
sin оG = sin оD sin ia; cos GD = cos оD/ cos оG,
и далее, определяем дугу
GD = arcos{cos a/[1 - (sin ·sin ia)2]0.5}. (28)
После подстановки значений дуг DA, AB, и GD в соотношение (22), получаем угловое расстояние перигелия Марса в неподвижных координатах, которое отсчитывается от неподвижной точки G на его орбите, в следующем виде:
p0 = a - а + arcsin[siniEa sin(а - 2)/sin iа] + arccos[cos /(1 - (sin sin iа)2)0.5]. (29)
Таким образом, используя наблюдения, можно вычислить наклон орбиты Марса к неподвижной плоскости экватора iа, долготу восходящего узла a с помощью соотношений (24) и (25), соответственно, и угловое расстояние перигелия p0 от неподвижной точки G на орбите Марса с помощью соотношения (29).
2.3 Сравнение вычисленных вековых возмущений элементов орбиты Марса со значениями, полученными из наблюдений
Для этого мы проинтегрировали дифференциальные уравнения (3) с использованием, приведенных выше, двух вариантов начальных условий на интервалах времени: -2.6 тыс. лет T 2 тыс. лет для первого варианта и -3.4 тыс. лет T 3.6 тыс. лет - для второго. Мы решали эти задачи с шагом интегрирования dt = 1?10-5 года и с «расширенной длиной» чисел ( с четырех кратной точностью), при которой вычисляемые координаты представляются числами с 35-ю десятичными знаками. Вычисления показали, что погрешность в этом случае пропорциональна [8] интервалу времени интегрирования и изменяется со скоростью, примерно равной dM/dt = 1.47·10-21 в столетие. Вычисления с заданными начальными условиями и параметрами на интервале времени в 100 млн. лет, привели бы к относительной погрешности, примерно равной, M = 1.47·10-15, что значительно меньше, приведенной нами ранее, величины M = 8 10-11. Не смотря на то, что погрешность уменьшилась в 54000 раз, результаты вычислений, по сравнению с результатами интегрирования «с двойной точностью» и с шагом интегрирования dt = 1?10-4 года, практически совпали. Это соображение говорит в пользу достоверности проведенных вычислений.
На рис. 3 точками представлена динамика параметров орбиты Марса: эксцентриситета e, долготы восходящего узла , угла наклона плоскости орбиты i и угла положения перигелия р0 в системе неподвижного экватора, вычисленная методам численного интегрирования системы (3), описывающей взаимодействия одиннадцати тел Солнечной Системы. Непрерывные линии представляют те же величины (ea, iа, а и р0а), полученные из наблюдений в промежутке времени, порядка 2500 лет, после некоторой математической обработки. По этой причине их достоверность ограничена этим периодом времени. Из графиков видно, что результаты наблюдений хорошо согласуются с вычисленными из дифференциальных уравнений значениями элементов e, i, и р в пределах нескольких тысяч лет «в прошлое» и «в будущее» от исходной эпохи. Наибольшее отклонение имеет наклон i при Т = -2.6 тыс. лет с относительной погрешностью (i -iа)/ iа = 7.4·10-4.
Рис. 3. Сопоставление рассчитанной (1) эволюции орбиты Марса за -2.6 ч +2 тыс. лет с аппроксимацией С. Ньюкомба данных наблюдений (2): e - эксцентриситет; i -наклон плоскости орбиты к плоскости экватора 1950.0 г. в радианах; - угловое положение восходящего узла орбиты от оси x на эпоху 1950.0 г. в радианах; угловое положение перигелия в плоскости орбиты: p0 - от неподвижной точки G на орбите (1, 2) и p - от восходящего узла в радианах (3). Т - время в тысячах юлианских лет от 30.12.1949 г.; интервал между точками - 200 лет.
Близость результатов «теоретических» вычислений с данными наблюдений, представленных «эмпирическими» соотношениями, свидетельствует в пользу достоверности наших расчетов. Анализ также показывает, что при изучении эволюции орбиты Марса в течение 300 непрерывных оборотов, появляются небольшие короткопериодические колебания в элементах орбиты и, прежде всего, в эксцентриситете орбиты Марса, период которых равен двенадцати годам. Эти колебания показаны на рис. 3.
На рис. 3 представлены результаты вычислений с использованием первого варианта начальных условий (согласно табл.1), отнесенных к эпохе 1950.0 г. Такие же вычисления были нами выполнены и со вторым вариантом начальных условий (согласно табл.2), но отнесенных к эпохе 2000.0 г. Их результаты еще ближе к «эмпирическим» соотношениям. Эти сравнения убеждают нас в справедливости, приведенных выше, количественных и качественных выводов.
3. Эволюция орбиты Марса в интервале времени 3 млн. лет
Эволюция вычисленных элементов орбиты Марса за 3 млн. лет «в прошлое» представлена на рис. 4. Она показывает, что эксцентриситет орбиты е испытывает короткопериодические изменения с амплитудой 0.019, с «главным» периодом, равным Те1 = 95.2 тыс. лет при среднем (за 50 млн. лет) значении еm = 0.066. Долгота восходящего узла цЩ меняется с «средним» периодом ТЩ = 73.1 тыс. лет вокруг среднего значения цЩт = 0.068 радиан. Угол наклона плоскости орбиты і к плоскости экватора (эпоха 1950 г.) испытывает колебания с таким же периодом Ті = 73.1 тыс. лет вокруг среднего значения іm = 0.405 радиан. Долгота перигелия цр практически линейно увеличивается со временем, т.е. перигелий перемещается в направлении обращения Марса вокруг Солнца, совершая в среднем один оборот за Tp = 76.8 тыс. лет. Угловая скорость ?p вращения перигелия колеблется вокруг среднего значения ?pm = 1687 " за столетие, при этом в момент времени t = -1.35 млн. лет она принимает отрицательное значение.
Рис. 4. Эволюция орбиты Марса за 3 млн. лет: обозначения те же, что и на рис. 3; Т - время в миллионах лет; интервал между соседними точками равен 10 тыс. лет; ?p - угловая скорость вращения перигелия в "/столетие на интервале времени в 20 тыс. лет; средняя за 50 млн. лет угловая скорость вращения перигелия ?pm = 1687 "/столетие - средняя за 50 млн. лет.
4. Эволюция орбиты Марса за 100 млн. лет.
На рис. 5а представлено изменение параметров орбиты Марса за 50 млн. лет. На нем отчетливо виден второй период изменения эксцентриситета, равный Те2 = 2.31 млн. лет.
Рис. 5. Эволюция орбиты Марса за 100 млн. лет (обозначения те же, что и на рис. 3 - 4).
При этих колебаниях достигаются «крайние» значения эксцентриситета е = 0.0014 и е = 0.126. Колебания угла наклона i происходят в пределах от 0.288 до 0.521 радиан, их диапазон составляет 13.37є. Долгопериодические колебания восходящего узла цЩ и наклона орбиты i обусловлены колебаниями ее угла нутации ? со средним периодом, равным Т? = 1.15 млн. лет. Следует отметить, что период Т? в два раза меньше Те2. Перигелий со средним за 50 млн. лет периодом Tp перемещается в направлении обращения Марса вокруг Солнца. Угловая скорость ?p вращения перигелия колеблется вокруг среднего значения ?pm. При сопоставлении с графиком для e видно, что возвратные движения перигелия происходят, когда эксцентриситет орбиты приближается к нулю. Из представленных графиков следует, что амплитуды колебаний и периоды стабильны, т.е. орбита Марса совершает «стационарные» движения.
На рис. 5б показано изменение параметров орбиты Марса в диапазоне времени -50 млн. лет T -100 млн. лет. Из сравнения графиков, на рисунках 5а и 5б, видно, что характер эволюции параметров орбиты в течение двух периодов времени по 50 миллионов лет каждый сохраняется. Периоды колебаний эксцентриситета е, восходящего узла орбиты и угла наклона её плоскости i к плоскости неподвижного экватора, а также их амплитуды и период Тp обращения перигелия остаются теми же. Аналогичные результаты получены для орбит других планет на том же промежутке времени в 100 млн. лет. Это свидетельствует об устойчивом характере планетных движений в Солнечной Системе.
5. Сравнение наших результатов с результатами других исследователей
В 1950 г. Д. Брауэр и А. Вурком опубликовали новую аналитическую теорию вековых возмущений элементов орбит планет Солнечной Системы [21]. На ее основе К. Остервинтер и др. [22] вычислили изменения элементов орбит больших планет, от Меркурия до Нептуна, за период времени -5 млн. лет T 5 млн. лет. В этой работе приведены, в частности, углы наклона плоскости орбиты к неподвижной эклиптике 1оА1 (см. рис. 2б). Для сравнения с ними мы определили значения наклонов планетных орбит к плоскости неподвижной эклиптики по значениям рассчитанных нами наклонов к плоскости неподвижного экватора. С этой целью мы рассмотрели треугольник 0DA1, в котором известны два угла = 0 и DA1 = - і, сторона между ними 0D = . И на основе сферической тригонометрии вычислили другие его параметры Угол наклона орбиты іе0 = A1D к неподвижной эклиптике может быть вычислен по «теореме косинусов» [8] и будет равен
іе0 = arccos[-cos0 cos( - і) + sin0 sin( - і) cos.]. (30)
Из этого сравнения видно, что значения эксцентриситетов е в начале (при T = 0) совпадают количественно, затем («после» -1 млн. лет) начинают различаться качественно из-за того, что второй период изменения эксцентриситета Те2 = 2 млн. лет меньше найденного нами периода (Те2 = 2.31 млн. лет), который, как будет показано на рис. 7, совпадает с периодом колебаний Ж.Ляскара и др. [4]. Это говорит о том, что в работе [22]
Рис. 6. Сравнение численных решений (1) эксцентриситета и наклона орбиты Марса в промежутке времени, равном 5-ти миллионам лет (от -5 млн. лет до «сегодняшнего дня», Т=0) с расчетами (2) К. Остервинтера и др. [22], полученных ими на основании теории Брауэра и Вуркома [21].
долгопериодические изменения эксцентриситета вычислены, по-видимому, с большой погрешностью. Что касается углов наклона, наши результаты и результаты, приведенные в работе [22], фактически совпадают. В статье [25] период больших колебаний наклона ie0 равен, Тie2 = 1.166 млн. лет, т.е. близок к «нашему периоду колебаний», Т? = 1.15 млн. лет, для угла нутации ?. Однако периоды коротких колебаний, приведенные в статье [22], в два с лишним раза больше вычисленных нами периодов колебаний и периодов колебаний, приведенных в работе [4] (см. рис. 7). Короткопериодические колебания углов наклона, приведенные в работе [22], по-видимому, вычислены с большей погрешностью. Однако, в отличие от работы [4], эволюция колебаний плоскости орбиты Марса, описанная в статье [22] так же, как и у нас, имеет устойчивый характер.
В работе [4] долгота перигелия L отсчитывается от подвижного восходящего узла (moving equinox), а угол наклона IMr плоскости орбиты (obliquity) отнесен к подвижной плоскости экватора Марса. Для дальнейшего необходимо привести эти результаты к неподвижному земному экватору и с этой целью выполним некоторые геометрические построения.
Восходящий узел находится на пересечении кругов обиты Марса DАВ и его экватора АMr А?Mr (см. рис. 2б). Восходящий узел Mra расположен на другом конце диаметра от нисходящего узла Mrd, показанного на этом рисунке. Поэтому долгота перигелия L определяется дугой MraDMrdB, которую можно записать в виде суммы дуг: