На правах рукописи

Заусаев Дмитрий Анатольевич

МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ДВИЖЕНИЯ

НЕБЕСНЫХ ТЕЛ С ИСПОЛЬЗОВАНИЕМ БАНКА ДАННЫХ

КООРДИНАТ БОЛЬШИХ ПЛАНЕТ

05.13.18 – Математическое моделирование, численные методы

и комплексы программ

Автореферат диссертации на соискание ученой степени кандидата физико-математических наук

Самара – 2013

Работа выполнена на кафедре «Прикладная математика и информатика»

в федеральном государственном бюджетном образовательном учреждении высшего профессионального образования «Самарский государственный технический университет»

Научный руководитель: Радченко Владимир Павлович, доктор физико-математических наук, профессор.

Официальные оппоненты: Курушина Светлана Евгеньевна, доктор физико-математических наук, доцент, ФГБОУ ВПО «Самарский государственный аэрокосмический университет имени академика С.П. Королева (национальный исследовательский университет)», кафедра физики, профессор Обрубов Юрий Викторович доктор физико-математических наук, профессор, Калужский филиал ФГБОУ ВПО «Московский государственный технический университет им. Н.Э. Баумана», кафедры «Высшая математика», профессор

Ведущая организация: федеральное государственное бюджетное учреждение науки «Институт астрономии Российской академии наук» (ИНАСАН)

Защита диссертации состоится 6 декабря 2013 года в 12 часов на заседании диссертационного совета Д 212.215.05, созданного на базе ФГБОУ ВПО «Самарский государственный аэрокосмический университет имени академика С.П. Королева (национальный исследовательский университет)» по адресу: 443086 Россия, г. Самара, Московское шоссе,

С диссертацией можно ознакомиться в библиотеке ФГБОУ ВПО «Самарский государственный аэрокосмический университет имени академика С.П. Королева (национальный исследовательский университет)».

Автореферат разослан _ октября 2013 года.

Ученый секретарь диссертационного совета д.т.н., доцент Востокин С.В.

ОБЩАЯ ХАРАКТЕРИСТИКА РАБОТЫ

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

Проблеме астероидной опасности посвящен ряд работ1, 2, 3, 4, 5, в которых дается оценка степени сближения и вероятности столкновения потенциально опасных астероидов с Землей. Помимо этого при решении задач, связанных с этой проблемой, требуется проводить регулярные исследования движения свыше 10000 астероидов групп Аполлона, Амура, Атона и более 250 короткопериодических комет на интервале времени порядка нескольких столетий. Даже при наличии современных средств вычислительной техники решение данной задачи вызывает определенные трудности. Создание высокоэффективных алгоритмов и программ для численного интегрирования уравнений движения небесных тел является необходимым условием для своевременного получения результатов проведенных исследований. Повышения эффективности работы программы при решении уравнений движения небесных тел можно достичь двумя путями: либо путем создания банка данных координат больших планет, либо путем модификации метода численного интегрирования.

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

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

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

Создание банка данных координат больших планет способствует понижению порядка системы дифференциальных уравнений с 72 до 6, что более чем на порядок сокращает расчетное время практически без потери точности решения. ОднаАстероидно-кометная опасность: вчера, сегодня, завтра под редакцией Б.М. Шустова, Л.В. Рыхловой. М.:

Физматлит, 2010. 384 с.

Башаков А.А., Питьев Н.П., Соколов Л.Л. Особенности движения астероида 99942 Апофис. – 2008. – Т. 42, № 1. – С. 20-29.

Виноградова Т.А., Железнов Н.Б., Кузнецов В.Б. Каталог потенциально опасных астероидов и комет // Тр. ИПА.

– 2003. – Т.9 – С. 11-218.

Заботин А.С.,Кочеткова О.М, Шор В.А. Сближение малой планеты (99942) Apophis 2004 MT4 с Землей в 2029 г.

//Всероссийская конференция «Астероиднокометная опасность – 2005. – С. 134-137»

Yeomans. D.K., Chesley S.R., Chodas P.W. NASA’s Near-Earth Object Program Office. Proceedings of the International Conference “Asteroid-Comet Hazard – 2009” Saint Petersburg “Nauka”, 2010. – P. 244-254.

Беляев Н.А. Эволюция орбиты кометы Даниэля 1909IV за 400 лет (1660-2060 гг.) // Бюлл.ИТА. 1966. Т.10. №10.

С.696-710.

Standish E.M. DE403/LE403: Annooncement // Jet Prop Lab Technical Report. IOM 314. 1995. P.10-124.

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

Вышеизложенное и определяет актуальность тематики диссертационного исследования.

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

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

1) разработка модифицированной математической модели для описания движения малых тел Солнечной системы, сближающихся с Землей, с использованием банка данных координат больших планет;

2) создание банка данных для барицентрических координат и скоростей больших планет, Луны и Солнца;

3) разработка алгоритмов и комплекса программ численного интегрирования уравнений движения небесных тел на основе метода Эверхарта с использованием банка данных координат больших планет;

4) исследование сходимости и устойчивости как для систем дифференциальных уравнений задачи Коши, так и для численного метода;

5) оценка погрешности математической модели и численных методов;

6) проведение исследования эволюции орбит небесных тел на интервале времени с 1800 по 2200 гг.

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

Научная новизна.

1. Разработана модифицированная математическая модель для описания движения малых тел Солнечной системы с использованием банка данных координат планет в форме оскулирующих элементов, коэффициентов полиномов Чебышева и Эверхарта, что, в отличие от существующих подходов, позволило понизить порядок системы дифференциальных уравнений с 72 до 6 и более чем на порядок сократить расчетное время без существенной потери точности решения.

2. Созданы алгоритмы и программное обеспечение с использованием банка данных координат планет для исследования эволюции орбит небесных тел.

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

4. Проведено исследование устойчивости решений системы дифференциальных уравнений модифицированным методом Эверхарта для различных небесных тел, представляющих потенциальную угрозу для Земли.

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

Теоретическая и практическая значимость работы.

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

2. Созданный программный комплекс для исследования движения небесных тел может быть использован для выявления объектов, представляющих опасность столкновения с Землей.

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

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

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

сопоставлением с результатами других исследователей.

На защиту выносятся следующие положения.

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

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

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

4. Результаты исследования эволюции орбит малых тел Солнечной системы, представляющих потенциальную угрозу для Земли.

Связь диссертационной работы с планами научных исследований.

Работа выполнялась в рамках плана НИР СамГТУ (тема «Разработка методов математического моделирования динамики и деградации процессов в механике сплошных сред, технических, экономических, биологических и социальных системах и методов решения неклассических краевых задач и их приложений».);

проекта Федерального агентства по образованию РФ (проект РНП 2.1.1.1689):

«Создание информационной среды на базе современных математических моделей и методов для исследования эволюции малых тел в Солнечной системе» в рамках аналитической ведомственной целевой программы: «Развитие научного потенциала высшей школы (2006–2008 гг.)»; проекта Министерства образования и науки РФ (проект РНП 2.1.1.745): «Создание научно-информационной базы данных эволюции орбит малых тел Солнечной системы, представляющих потенциальную опасность для Земли» в рамках аналитической ведомственной целевой программы: «Развитие научного потенциала высшей школы (2009–2010 гг.)»; проекта Министерства образования и науки РФ (проект РНП 2.534.2011) «Разработка математического и программного обеспечения для исследования эволюции орбит главных метеорных потоков (2011-2012 гг.)».

Апробация работы. Основные результаты диссертации докладывались на десятой Международной конференции «Актуальные проблемы современной науки» (Самара, 2010 г.); на Всероссийских научных конференциях с международным участием «Математическое моделирование и краевые задачи» (Самара, 2009, 2010, 2011, 2013); XIX Международной конференции «Решетневские чтения» (г.

Красноярск, 2010 г.); Международной молодежной научной конференции по естественнонаучным и техническим дисциплинам «Научному прогрессу – творчество молодых» (г. Йошкар-Ола, 2010 г.); на Седьмой Международной конференции «Математическое моделирование физических, экономических, технических, социальных систем и процессов» (г. Ульяновск, 2009 г.); на научном семинаре «Механика и прикладная математика» Самарского государственного технического университета (рук. проф. Радченко В.П., 2011-2013 гг.).

Публикации. По теме диссертации опубликовано 16 печатных работ, из которых 6 входят в список изданий, рекомендованных ВАК. Список публикаций приведен в конце автореферата.

Личный вклад автора. Работы [6,9,10,11,13,14] выполнены самостоятельно, в работах [1-5, 12, 15] диссертанту принадлежит совместная постановка задачи, лично соискателем построены решения задач, разработано алгоритмическое и программное обеспечение, выполнены расчеты и анализ результатов.

Структура и объем диссертации. Диссертация состоит из введения, пяти глав, общих выводов, списка литературы и двух приложений, в которых приведены выдержки из печатного варианта каталога орбитальной эволюции короткопериодических комет и листинги разработанных программ. Общий объем диссертации составляет 164 страницы, включая 4 рисунка и 55 таблиц. Библиографический список включает 137 наименований.

СОДЕРЖАНИЕ РАБОТЫ

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

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

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

Рассматривается построение вычислительных алгоритмов метода Эверхарта. Подробно описывается модификация метода Эверхарта, позволяющая расширить его практическое использование для численного интегрирования уравнений движения небесных тел до 33 порядка включительно. Доказывается связь неявного метода Эверхарта с неявными методами Рунге-Кутты. Рассматриваются вопросы, связанные с устойчивостью метода Эверхарта. Проводится сопоставление метода Эверхарта с неявным многошаговым методом Адамса-Мултона. Отмечается преимущество метода Эверхарта перед многошаговым методом Адамса-Мултона по точности и устойчивости. Однако по быстродействию метод Эверхарта уступает многошаговому методу Адамса-Мултона.

В результате проведенных исследований метода Эверхарта были доказаны следующие теоремы.

Теорема 1. Неявный метод Эверхарта принадлежит классу неявных методов Рунге-Кутты.

Теорема 2. Метод Эверхатра, использующий квадратурную формулу Лобатто, А – устойчив.

Третья глава посвящена созданию различных банков данных для получения координат и компонент скоростей возмущающих планет (Меркурий-Плутон) и Луны на основе современной численной теории DE-405. Рассмотрены три формы хранения информации, с помощью которой с высокой степенью точности вычисляются координаты и скорости возмущающих планет на любой момент времени в интервале с 1600 по 2200 гг. Показано, что создание банка данных координат больших планет и использование их в программах численного интегрирования позволяет понизить порядок системы дифференциальных уравнений с 72 до 6, что более чем на порядок сокращает расчетное время. Однако эффективность работы программ существенным образом зависит от формы хранения координат возмущающих тел.

Банк данных в форме оскулирующих элементов включает в себя для каждой планеты систему оскулирующих элементов. При этом внутри каждого промежутка координаты планеты вычисляются по формулам невозмущенного движения, а при переходе от одного промежутка к другому меняются системы оскулирующих элементов планеты, поэтому точность координат существенным образом зависит от длины промежутка, в пределах которых движение планеты принято считать невозмущенным. Рассмотрено два способа хранения положений больших планет в форме оскулирующих элементов орбит с шагом один и десять дней. Показано, что при способе хранения оскулирующих элементов с шагом 1 день (вместо 10 дней) уменьшается погрешность в координатах планет и Луны более чем на порядок.

Максимальные различия в координатах Земли и Луны, полученных с помощью банка данных DE405 и по оскулирующим элементам, составляет 8,6 10 а.е.

(астрономических единиц) и 7,0 10 5 а.е.

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

Самые точные координаты и скорости планет получены с помощью банка данных коэффициентов полиномов Чебышева с шагом 10 дней, созданного на основе банка DE405. На интервале времени с 1660 по 2200 гг. был создан банк данных, представляющий собой последовательность групп коэффициентов разложений координат и скоростей больших планет, Луны и Солнца в ряды по полиномам Чебышева. В каждой группе объединены коэффициенты разложений, соответствующие одному 10-дневному интервалу. Хронологически упорядоченная последовательность этих групп записана на жестком диске в двоичном формате в виде файла прямого доступа. Максимальные различия в координатах Земли и Луны при сопоставлении данного банка с банком данных DE405 составляют 2 1012 а.е. и 9 1012 а.е. соответственно.

На точность прогнозирования движения небесных тел влияет множество факторов. Основными из них являются: степень точности начальных данных; учет физических сил в математической модели; применяемый метод численного интегрирования уравнений движения исследуемого объекта; проблема, связанная с устойчивостью как самой задачи Коши, так и применяемого метода.

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

Брумберг В.А. Релятивистская небесная механика. М.: Наука, 1972. 382 с.

где f – постоянная тяготения; m0, mi – массы возмущаемого тела и возмущающих планет; x, y, z – координаты возмущаемого тела; xi, yi, zi – координаты возмущающих планет; x, y, z – компоненты скорости возмущаемого тела;

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

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

Разработаны алгоритмы и программы численного интегрирования движения небесных тел методом Эверхарта с использованием банков данных координат планет и Луны. Для определения области применимости математического и программного обеспечения исследовалась эволюция орбит 10 короткопериодических комет и 10 астероидов, принадлежащих к группам Аполлона, и Атона. Элементы орбит комет и астероидов взяты из каталогов9, 10, причем все 10 комет взяты из начала списка каталога. Это кометы: P/Halley, P/Encke, D/Biela, P/Faye, D/Brorsen, Заусаев А.Ф., Заусаев А.А. Каталог орбитальной эволюции короткопериодических комет с 1800 по 2104 гг. М.:

Машиностроение - 1. 2005. 410 с.

Заусаев А.Ф., Абрамов В.В., Денисов С.С. Каталог орбитальной эволюции астероидов, сближающихся с Землей с 1800 по 2204 гг.. М.: Машиностроение - 1, 2007. 608 с P/d'Arrest, P/Pons-Winnecke, P/Tuttle, P/Tempel 1, P/Tempel 2, в то время как 10 астероидов находятся в конце каталога. Это астероиды: Aten/2001 CP36, Apollo / 2001 BF10, Aten / 2002 VX91, Apollo / 2004 BN41, Apollo / 2005 WM3, Apollo / 2003 HP32, Aten / 2003 HT42, Apollo / 2005 VG7, Aten / 1998 SZ27, Aten / 1994 GL.

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

Численное интегрирование для каждого объекта было проведено на интервале времени с 1800 по 2200 гг. тремя методами с использованием банков данных координат планет и Луны. Результаты вычислений сравнивались с данными сайта smallbodies.ru и между собой. Результаты сопоставления с данными сайта smallbodies.ru показали, что максимальные расхождения в элементах кометных орбит для всех 10 комет на конце интервала интегрирования незначительны. При этом в угловых элементах максимальное расхождение наблюдается в средней аномалии у кометы D/Biela – 0.0043 градуса, а отличия в остальных элементах орбит не превышают точности оптических наблюдений.

Следует отметить, что результаты численного интегрирования уравнений движения для этих же 10 комет с использованием полиномов Чебышева, полиномов Эверхарта и оскулирующих элементов практически не отличаются друг от друга. Так как все четыре метода дают незначительные расхождения в элементах орбит, это позволяет считать, что их можно применять для исследования эволюции орбит короткопериодических комет на интервале времени с 1800 по 2190 гг.

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

Результаты численного интегрирования уравнений движения этих астероидов на интервале времени с 2006 по 2190 гг. получены четырьмя методами: путем совместного интегрирования уравнений движения астероида с возмущающими планетами, с использованием оскулирующих элементов, полиномов Чебышева, а также полиномов Эверхарта. Каждый из рассматриваемых астероидов на интервале времени с 1800-2200 гг. сближается с Землей на расстояние менее 0,01 а.е.

Как следует из результатов вычислений, оскулирующие элементы орбит, полученные различными методами, отличаются друг от друга. Однако степень различия оскулирующих элементов орбит для различных астероидов неодинакова. У двух из этих астероидов – у Apollo /2001 BF10 и Aten /2002 VX91 – расхождение оскулирующих элементов орбит, полученные различными методами, наибольшие.

Особенно большие расхождения имеют место в средней аномалии. У астероида Apollo /2001 BF10 разброс в средней аномалии, полученной различными методами, составляет около 20 градусов, а у астероида Aten /2002 VX91 – свыше 50 градусов.

Основная причина этих расхождений связана с орбитальной устойчивостью.

Существенное нарушение орбитальной устойчивости происходит при умеренных и тесных сближениях астероидов с большими планетами. Тестовые расчеты показали, что изменение начальных данных большой полуоси у астероида Aten / VX91 на величину a 10 7 а. е. приводит к расхождению средней аномалии в конце интервала интегрирования почти на 60 градусов. Вследствие того, что орбиты астероидов Apollo /2001 BF10 и Aten / 2002 VX91 лежат вблизи плоскости эклиптики, а большие полуоси близки к большой полуоси Земли, возмущающее действие Земли постоянно оказывает существенное влияние на движение этих астероидов, являясь причиной орбитальной неустойчивости Следует отметить, что наряду с учетом действующих сил на точность полученных результатов существенное влияние оказывает метод решения дифференциальных уравнений. В методе Коуэлла уравнения движения в прямоугольных координатах интегрируются непосредственно. В методе же Энке координаты не получаются непосредственно, а вместо этого интегрирование дает разности между действительными координатами и координатами в оскулирующей орбите, т.е.

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

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

Данные численного интегрирования уравнений Энке сопоставлены с данными, размещенными на сайте http://smallbodies.ru, которые получены путем совместного интегрированием. Показано, что максимальные расхождения в элементах кометных орбит для всех 10 комет в обоих методах на конце интервала интегрирования незначительны. При этом в угловых элементах максимальное расхождение наблюдается в средней аномалии у комет P/Encke и D/Biela – 0.0030 и 0. градуса соответственно, что составляет около 11 и 15 угловых секунд. Отличия в остальных элементах орбит не превышают точности оптических наблюдений.

Так как оба метода дают незначительные расхождения в элементах орбит, то это позволяет считать, что метод Энке можно использовать для исследования эволюции орбит короткопериодических комет на интервале времени около 400 лет.

Аналогичные исследования были проведены для 10 астероидов, принадлежащих к группам Аполлона и Атона.

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

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

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

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

Если в уравнениях движения (1) использовать фиктивное время, связанное с физическим соотношением d t / d s = r, то в векторной форме уравнения (1) примут вид где штрих указывает на дифференцирование по s, t r ; r – вектор с координатаxx j yy j zz j Однако уравнения (2) все еще являются сингулярными.

Второй этап регуляризации заключается в устранении сингулярности относительно центрального тела в самих дифференциальных уравнениях. Это достигается путем преобразования вектора физического пространства в четырехмерное параметрическое пространство с помощью матрицы преобразования Кустанхеймо Штифель Е., Шейфеле Г. Линейная и регулярная небесная механика. М. Наука, 1975. 304 с.

Everhart E. Implicit single methods for integrating orbits // Celestial mechanics. 1974. №.10. Р.35-55.

где u – четырехмерный вектор с параметрами u1, u 2, u3, u 4, которые выражаются через трехмерные гелиоцентрические координаты по формулам:

Формулы связи для скорости имеют следующий вид:

Переходя в уравнении (2) от старых переменных x, y, z к новым переменным u1, u2, u3, u4, получим следующую систему дифференциальных уравнений:

Таким образом, система (6) является регулярной при соударении ( u 0 ) с центральным телом при условии, что возмущающие силы имеют конечные значения.

Гелиоцентрические координаты по параметрам 4-мерного вектора вычисляются по формулам:

Гелиоцентрические скорости по производным параметров по фиктивному времени находятся из соотношений:

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

Тогда построение решения уравнений движения возмущаемого тела по формулам (2)-(9) можно проводить в следующей последовательности.

1. Зная координаты и скорости возмущаемого тела на заданный момент времени, по формулам (4), (5) находим компоненты 4-мерного вектора.

2. Значения величин, входящих в уравнение (6), находим с помощью вышеуказанных соотношений. При этом координаты больших планет выбираются из банка данных на известные моменты физического времени.

3. Гелиоцентрические координаты и скорости возмущаемого тела по компонентам координат и скоростей 4-мерного вектора вычисляем по формулам (7), (8).

4. Для решения систем дифференциальных уравнений (6) рекомендуется использовать высокоточные численные методы.

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

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

На основании проведенных исследований относительно эволюции орбиты астероида Апофис можно сделать следующие выводы:

а) использование начальных данных элементов орбит астероида Апофис, заданных на различные моменты времени, не позволяет получить удовлетворительный результат в определении величины минимального сближения астероида с Землей;

б) эволюция элементов орбит Апофиса после его сближения с Землей существенным образом зависит от условий и величины сближения;

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

В диссертации разработано 6 вариантов алгоритмов и программ для проведения численного интегрирования уравнений движения небесных тел с использованием численного метода интегрирования Эверхарта. Основное отличие рассматриваемых алгоритмов заключается как в использовании различных математических моделей, описывающих движение объекта, так и в различных формах хранения координат возмущающих тел. Сравнение эффективности различных алгоритмов проводилось на примере интегрирования уравнений движения 5 короткопериодических комет. В их число вошли кометы: Галлея с относительно большим периодом обращения; комета Энке с малым периодом обращения и три кометы: Брорзена, Понса-Виннике и Брукса 2, имеющие тесные сближения с Юпитером. Кометы являются сложными объектами для исследования, так как без учета негравитационных эффектов не представляется возможным исследовать их движение с высокой степенью точности. Однако элементы их орбит в момент прохождения через перигелий могут быть использованы при сравнении методов численного интегрирования. Средняя аномалия в момент прохождения кометы через перигелий равна нулю или 360 градусам, поэтому по близости к этой величине можно судить об адекватности работы численного метода.

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

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

Проведенные в диссертации исследования показали, что при отсутствии тесных сближений с большими планетами не возникает трудностей при численном интегрировании уравнений движения астероидов на интервале времени 400 лет (1800-2200 гг.). Однако после сближения с Землей возникают проблемы, связанные с устойчивостью решений. Смещение начальных данных даже в пределах среднеквадратических погрешностей приводит к значительным ошибкам в результатах вычисления.

На основании проведенных исследований можно сделать следующий вывод:

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

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

Начальные данные для каждого астероида брались на эпоху 08.01.2013, результаты вычислений сопоставлялись на дату после сближения астероида с планетой.

Сопоставление элементов орбит проводилось для астероидов: 2005 WY55, YR3, 2005 YU55, 2004 MN4 Apophis, 2001 WN5, 1999 AN10, 2002 CU11, RQ252 через 10 лет после их сближений с Землей.

Результаты сопоставления показали, что наибольшие различия в элементах орбит с данными NASA получились у астероидов 2004 MN4 Apophis и 2005 YR3, т.е. у астероидов, проходящих через сферу действия Земли. Расхождения в значениях средних аномалий составили около 4 градусов. Причиной этих расхождений являются проблемы, связанные с орбитальной устойчивостью. Различие элементов орбит у других 6 астероидов с данными NASA является вполне удовлетворительным, так как расхождения средних аномалий составляют не более 0.3 градусов. Различия в остальных элементах орбит находятся в пределах среднеквадратических погрешностей оптических наблюдений.

В пятой главе дается описание разработанного программного обеспечения.

Исследование движения астероидов и короткопериодических комет сопряжено с большим объемом вычислений. Успешное выполнение проводимых исследований существенным образом зависит от эффективности алгоритмов и программ, используемых при решении поставленных задач. Значительную роль в реализации данной задачи играет выбор языка и среды программирования. В качестве языка программирования был выбран Си++, который поддерживает парадигму объектно-ориентированного программирования. Данный язык широко используется для создания высокоэффективных комплексов программ. В качестве среды для разработки программы была выбрана среда Microsoft Visual Studio Express Edition.

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

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

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

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

Банк данных координат больших планет в форме коэффициентов полиномов Чебышева позволяет определять координаты и скорости планет приблизительно на порядок точнее по сравнению с банком данных полиномов Эверхарта.

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

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

Элементы новизны разработанного программного комплекса определяются следующим:

а) в банке данных предложен новый способ хранения координат планет в форме оскулирующих элементов орбит, полиномов Эверхарта и коэффициентов полиномов Чебышева;

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

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

В ходе выполнения диссертационной работы получены следующие результаты:

1. Разработана модифицированная математическая модель для описания движения малых тел Солнечной системы с использованием банка данных координат планет в форме оскулирующих элементов, коэффициентов полиномов Чебышева и Эверхарта, что позволило понизить порядок системы дифференциальных уравнений с 72 до 6 и сократить расчетное время в 12 – 18 раз.

2. Созданы алгоритмы и программное обеспечение с использованием банка данных координат планет для исследования эволюции орбит небесных тел.

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

4. Проведено исследование сходимости и устойчивости решений дифференциальных уравнений модифицированным методом Эверхарта для различных небесных тел, представляющих потенциальную угрозу для Земли.

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

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

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

СПИСОК ОСНОВНЫХ ПУБЛИКАЦИЙ В РЕЦЕНЗИРУЕМЫХ

ЖУРНАЛАХ ИЗ ПЕРЕЧНЯ ВАК:

1. Заусаев, А.Ф. Методы интерполяции, используемые для получения координат и элементов орбит больших планет и малых тел Cолнечной системы [Текст] / А.Ф. Заусаев, Д.А. Заусаев // Вестник Самар. гос. техн. ун-та. Серия: Физ.мат. науки. 2008. № 2 (17). С. 231-238.

2. Заусаев, А.Ф. Интегрирование уравнений движения малых тел Солнечной системы методом оскулирующих элементов [Текст] / Заусаев А.Ф., Заусаев Д.А.

// Вестник Самар. гос. техн. ун-та. Серия: Физ.-мат. науки. 2009. № 1 (18).

С. 228-238.

3. Заусаев, А.Ф. Численное интегрирование уравнение движения малых тел Солнечной системы с использованием оскулирующих элементов больших планет [Текст] / А.Ф. Заусаев, Д.А. Заусаев // Вестник Самар. гос. техн. ун-та. Серия:

Физ.-мат. науки. 2009. № 2 (19). С. 231 - 239.

4. Заусаев, Д.А. Численное интегрирование уравнение движения небесных тел с учетом регуляризации и использованием оскулирующих элементов больших планет [Текст] / Д.А. Заусаев, Л.А. Соловьев // Вестник Самар. гос. техн. унта. Серия: Физ.-мат. науки. 2010. № 5 (21). С. 305-308.

5. Радченко, В.П. Использование банка данных координат больших планет для численного интегрирования уравнений движения небесных тел. [Текст] / В.П., Радченко, Д.А. Заусаев, // Вестник Самар. гос. техн. ун-та. Серия: Физ.мат. науки. 2011. № 3 (24). С. 202-207.

6. Заусаев, Д.А. Использование дифференциальных уравнений в форме Энке для исследования движения малых тел Солнечной системы [Текст] / Д.А. Заусаев // Вестник Самар. гос. аэрокосмического ун-та имени академика С.П. Королева (национальный исследовательский университет). 2013. №5(36). Часть 2. С.

376-385.

Патенты и авторские свидетельства 7. Заусаев, Д.А. Программа для численного интегрирования уравнений движения небесных тел методом Энке. Роспатент. Свидетельство о государственной регистрации программы для ЭВМ №2013616434 от 08.06.2013.

Публикации в прочих изданиях 8. Заусаев, А.Ф. Численное интегрирование уравнений движения малых тел Солнечной системы с использованием оскулирующих элементов больших планет [Текст] / А.Ф. Заусаев, Д.А. Заусаев // Математическое моделирование и краевые задачи. Тр. Шестой Всерос. науч. конф. с междун. участием. Ч. 3. Самара:

СамГТУ. 2009. С. 125–130.

9. Заусаев, А.Ф. Метод оскулирующих элементов, используемый для исследования эволюции орбит малых тел Солнечной системы [Текст] / А.Ф. Заусаев, Д.А. Заусаев // Тр. Седьмой Междун.конф. «Мат. моделирование физических, экономических, технических, социальных систем и процессов». Ульяновск:

УлГУ. 2009. С. 109–110.

10. Заусаев, Д.А. Численный метод оскулирующих элементов, используемый для решения уравнений движения малых тел Солнечной системы [Текст] / Д.А.

Заусаев // Тезисы докладов Международной научной студенческой конференции по естественным и технических дисциплинам. Йошкар-Ола: МарГТУ.

11. Заусаев, Д.А. Применение банка данных оскулирующих элементов больших планет к исследованию эволюции астероидов, сближающихся с Землей за период с 1600-2200 гг [Текст] / Д.А. Заусаев //Актуальные проблемы современной науки. Естественные науки. Ч. 1-3 Математика, математическое моделирование, механика. Труды 5-го Международного форума (10-й Международной конференции) Самара: СамГТУ. 2010. С. 90-95.

12. Заусаев, Д.А. Банк данных оскулирующих элементов больших планет и Луны за период 1600-2200 гг [Текст] / Д.А. Заусаев // Математическое моделирование и краевые задачи. Тр. Седьмой Всерос. науч. конф. с междун. участием. Ч.3.

Самара: СамГТУ. 2010. С. 119-123.

13. Заусаев, А.Ф. Численное интегрирование уравнение движения небесных тел с использованием оскулирующих элементов больших планет [Текст] / А.Ф.

Заусаев, Д.А. Заусаев // Решетневские чтения. Материалы XIV международной конференции, посвященной памяти генерального конструктора ракетнокосмических систем академика М. Ф. Решетнева. Ч. 2. Красноярск. 2010. С 675Заусаев, Д.А. Использование банка данных координат больших планет для численного интегрирования уравнений движения астероидов групп Аполлона и Атона [Текст] / Д.А. Заусаев // Математическое моделирование и краевые задачи. Тр. Седьмой Всерос. науч. конф. с междун. участием. Ч. 3. Самара:

СамГТУ. 2011. С. 77-82.

15. Заусаев, Д.А. Использование банка данных координат больших планет для численного интегрирования уравнений движения астероидов групп Аполлона и Атона [Текст] / Д.А. Заусаев // Математическое моделирование и краевые задачи. Тр. Девятой Всерос. науч. конф. с междун. участием. Ч. 3. Самара: СамГТУ. 2013. С. 25-28.

16. Заусаев, А.Ф. Исследование устойчивости решений уравнений движения астероида 2004 MN4 Apophis после его сближения с Землей [Текст] / А.Ф. Заусаев, Д.А. Заусаев // Математическое моделирование и краевые задачи. Тр. Девятой Всерос. науч. конф. с междун. участием. Ч. 3. Самара: СамГТУ. 2013. С.

Подписано в печать 24.10.2013. Формат 60 84 1/16.

ФГБОУ ВПО «Самарский государственный технический университет»

443100, г. Самара, ул. Молодогвардейская, 244.