В. В. Абрамов
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ДВИЖЕНИЯ МАЛЫХ ТЕЛ СОЛНЕЧНОЙ СИСТЕМЫ НА ОСНОВЕ МЕТОДОВ АДАМСА
Рассмотрены методы Адамса численного интегрирования дифференциальных уравнений. Проведено исследование эволюции орбиты астероида 99942 Apophis на интервале времени с 1800 по 2206 гг., определены моменты его сближений с другими небесными телами. Результаты вычислений хорошо согласуются с исследованиями, проводимыми в нашей стране и за рубежом.
В данной работе проводилось исследование многошаговых методов Адамса, на основе которых осуществлялось математическое моделирование движения небесных тел Солнечной системы.
Методы Адамса относятся к числу наиболее экономичных методов, так как при одинаковой точности на одном шаге численного интегрирования требуется меньшее количество вычислений правых частей дифференциальных уравнений по сравнению с одношаговыми методами. Также несомненным достоинством методов Адамса является их универсальность, то есть они могут быть применены для решения широкого класса задач, описываемых дифференциальными уравнениями. Следовательно, усложнение математической модели в результате учёта дополнительных факторов, таких как релятивистские эффекты и несферичность тел, приводящие к изменению дифференциальных уравнений, не потребует изменения самого алгоритма численного интегрирования.
Общая форма записи многошаговых методов имеет вид [1]
Если ак = 1, а£-1 = -1, а£-2 = • = а0 = 0, то выражение (1) представляет собой общую форму записи методов Адамса.
Из формулы (1) следует, что для вычисления значений уп необходимо сначала получить к начальных значений Уо, У1, • , Уk-l. При = 0 метод (1) называется явным многошаговым методом. При ^ 0 правая часть (1) содержит /(хп +к,уп+£), и необходимо решать нелинейное уравнение относительно уп : метод (1) в этом случае называется неявным многошаговым методом.
Например, выражение
является формулой метода Адамса-Мултона 6-го порядка.
Явные методы Адамса-Бэшфорта значительно проще неявных методов Адамса-Мултона, однако неявные методы являются более точными, поэтому они являются более предпочтительными для решения практических задач. При этом требуется решать нелинейное уравнение, например с помощью итерационного процесса - метода прогноза и коррекции:
п = 0,1,2, •
(1)
Уп +1 = Уп + 60480 (19087/ +1 + 65112Гп - 46461/п_1 + +37504/п-2 - 20211/п-з + 6312/п_А - 863,/п_ъ)
Уп +1 = Уп +
к
Р:
Е:
С:
с : у„+1 = у„ + кМо/П+1 + М,/(х„+1_г, у„+1_г) = у„+1;
г=1
Е : /„+1 = / (х„+1, ул+1). (2)
В данном методе в качестве начального приближения используется значение, полученное по формуле метода Адамса-Бэшфорта, а дальнейшие уточнения производятся с помощью метода Адамса-Мултона [1]. Итерационный процесс Р(ЕС)пЕ в методе (2) выполняется либо фиксированное число раз, либо до тех пор, пока не будет достигнута сходимость.
В настоящей работе исследование движения небесных тел проводилось в экваториальной барицентрической системе координат. Дифференциальные уравнения движения с учётом релятивистских эффектов имеют вид [2]
,_ 2Ф+) ^_ 26_1 ^ + ТЦ. + (1 + т)4_ 2СЪ+Г) _
С к ф1г1к с к Ф]г]к с с с
2с2
-|2
(Г Т г,) • Г
Г
,
+А<- Г) • Г,}+Й{ ТТ Т Т)]>[(2+2у)г - (1+пу,] }>
/■> 2 ' У '' У 2
2с J с і ^гГї
ХГ (3)
2с І*і Г,
Здесь г, Г, Г — координаты, скорость, ускорение в барицентрической системе координат /—того возмущаемого тела; г, , г,, Г, — координаты, скорость, ускорение в барицентрической
22
системе координат , -го возмущающего тела; |т, = к т,, где к — гравитационная постоянная и т, — массатого тела; г, = \г, — г|, Ь и у — релятивистские параметры; Ь = у = 1; и = | Гг| и с — скорость света.
Начальные данные координат и скоростей больших планет, Луны и Солнца были взяты из банка данных координат и скоростей этих объектов, описанного в работе [3] и согласованного с базой данных координат и скоростей в ББ 405 [4]. Закладка начальных данных координат и скоростей больших планет проводилась с помощью одношагового метода Эверхарта [5] по программе ЯАБА-27 [6].
Т а б л и ц а 1
Сближения астероида 99942 АрорЫ* с Землёй и Луной
Объект Дата Время а. е. км
Земля 13.04.1823 23:50:12 0,007648 1144158
Луна 14.04.1823 07:09:07 0,008391 1255234
Луна 13.04.1866 18:10:56 0,014404 2154828
Земля 13.04.1866 19:26:58 0,015762 2358005
Земля 13.04.2029 21:46:22 0,000258 38618
Луна 14.04.2029 14:35:31 0,000631 94437
Луна 13.04.2066 15:57:53 0,001120 167605
Земля 14.04.2066 01:35:36 0,002995 448113
Луна 12.04.2085 22:43:58 0,007848 1174026
Земля 13.04.2085 04:42:14 0,009351 1398892
Луна 28.09.2117 07:35:36 0,012682 1897199
Земля 28.09.2117 19:10:50 0,012474 1866059
Луна 29.09.2127 00:24:46 0,012212 1826944
Земля 29.09.2127 01:38:29 0,011941 1786292
Земля 28.09.2137 17:56:32 0,011831 1769913
Луна 29.09.2137 11:49:37 0,011214 1677585
Земля 16.04.2156 07:54:37 0,012148 1817373
Луна 16.04.2156 10:14:35 0,010904 1631192
С помощью метода Адамса-Мултона было проведено исследование эволюции орбиты астероида 99942 АрорЬІ8 на интервале времени с 1800 по 2206 гг. В табл. 1 приведены все сближения данного астероида с Землёй и Луной на расстояние менее 0,02 а. е. на исследуемом интервале времени. Расстояние от астероида до центра Земли или Луны в момент сближения указано в астрономических единицах и километрах.
Как видно из табл. 1, минимальное сближение с Землёй произойдёт 13 апреля 2029 года в 21:46 по всемирному времени на геоцентрическом расстоянии 0,000258 а. е., что составляет примерно 38 тыс. км или 5 земных радиусов от поверхности нашей планеты. Спустя 17 часов Апофис пройдёт на минимальном расстоянии от центра Луны, которое составит 0,000631 а. е. или примерно 94 тыс. км. Данные результаты очень хорошо согласуются с исследованиями, проводимыми как в нашей стране, так и за рубежом.
В табл. 2 представлены все элементы орбиты астероида Апофис кроме средней аномалии М для некоторых моментов времени. Таким образом, табл. 2 показывает эволюцию орбиты этого малого тела на исследуемом временном интервале.
Т а б л и ц а 2
Элементы орбиты астероида 99942 Apophis
Дата a e w п i
05.01.1800 0,934664 0,188071 116,5963 211,4251 3,118633
11.05.1823 0,929381 0,189654 119,0503 210,4626 3,266028
09.05.1866 0,926224 0,190483 121,2233 209,0425 3,293669
16.03.2029 0,922299 0,191267 126,7488 203,8120 3,348612
11.05.2029 1,100632 0,188382 72,10625 203,5068 2,252664
17.03.2066 1,101449 0,188610 72,83814 202,7343 2,273405
11.05.2066 1,077099 0,183010 79,15851 202,5629 2,183697
09.05.2085 1,069004 0,181526 81,41254 202,2513 2,183262
25.10.2117 1,072704 0,181895 83,22313 200,7233 2,302441
26.10.2127 1,072770 0,181590 84,41660 199,5995 2,437007
26.10.2137 1,073652 0,181248 85,69574 198,6057 2,571903
13.05.2156 1,068532 0,180517 87,90469 198,2243 2,524838
05.08.2206 1,067732 0,180361 89,06481 197,5835 2,526403
Из табл. 2 видно, что до максимального сближения с Землёй в 2029 году астероид Апофис принадлежит к группе Атона, так как его большая полуось a < 1, а после этого сближения орбита астероида резко изменит свою форму, и он перейдёт в группу Аполлона (a > 1). Остальные сближения не приводят к кардинальным изменениям элементов орбиты.
Отклонения координат небесных тел, полученных в результате численного интегрирования уравнений движения (3) с помощью метода Адамса- Мултона 11 порядка с шагом 0,25 дня и модифицированного метода Эверхарта 27 порядка с шагом 1 день, очень малы и находятся в пределах точности оптических наблюдений. При этом вычисления с помощью метода прогноза и коррекции производятся примерно в четыре раза быстрее, чем с помощью метода Эверхарта. Данные результаты указывают на высокую эффективность методов Адамса при исследовании движения небесных тел на длительных интервалах времени.
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Холл Дж., Уатт Дж. Современные численные методы решения обыкновенных дифференциальных уравнений. — М.: Мир, 1979. — 312 с.
2. NewhallX. X., Standish E. M. Jr., Williams J. G. DE 102: a numerically integrated ephemeris of the Moon and planets spanning forty-four centuries// Astron. and Astrophys., 1983. — No. 125. — P. 150-167.
3. Заусаев А. Ф., Заусаев А. А. Каталог орбитальной эволюции короткопериодических комет с 1900 по 2100 гг. — М.: Машиностроение-1, 2005. — 346 с.
4. Standish E. M. JPL Planetary and Lunar Ephemerides, DE 405 / LE 405// Jet Prop Lab Technical Report. IOM 312. F-048. 1998. — P. 1-7.
5. EverhartE. Implicit single methods for integrating orbits// Celestial Mech., 1974. —No. 10. — P. 35-55.
6. Заусаев А. Ф., Заусаев А. А., Ольхин А. Г. Применение метода Эверхарта 31 порядка для решения уравнений движения больших планет// ВАК-2004. Горизонты Вселенной: Тез. докл. на Всерос. астрономич. конф. — М.: МГУ, ГАИШ, 2004. — С. 209.
Работа выполнена при финансовой поддержке Федерального агентства по образованию (проект РНП. 2.1.1.1689).
Поступила 19.09.2006 г.