The study of the efficiency of gauss-everhart and lobbie integrators in the problems of asteroid dynamics.pdf Введение Численное моделирование является преимущественным методом анализа и прогнозирования движения астероидов, спутников, космического мусора и других небесных тел. В большинстве случаев дифференциальные уравнения движения не имеют аналитического решения, поэтому единственным выходом является численное интегрирование. В настоящее время разработано множество методов численного интегрирования, однако все возрастающие требования к точности делают актуальной разработку новых более эффективных интеграторов. Кроме того, несмотря на повсеместное использование многопроцессорных вычислительных систем, важным фактором остается время, затрачиваемое на вычисление. Это особенно критично в случае предсказания возможного столкновения астероидов с Землей, поскольку для точной оценки требуется исследовать большое число тестовых частиц. Работа посвящена изучению эффективности применения интегратора Гаусса - Эверхарта [1], который является улучшенной версией оригинального интегратора Эверхарта [2], и нового коллокационного интегратора Lobbie [3] к решению задач астероидной динамики. В процессе исследования были проведены численные эксперименты по моделированию орбит астероидов с различными эксцентриситетами, рассматривалась эффективность использования интеграторов для решения невозмущенной и возмущенной задачи двух тел, а также смешанной системы уравнений первого и второго порядков с параметром хаотичности MEGNO (Mean Exponential Growth factor of Nearby Orbits) [4]. 1. Дифференциальные уравнения движения Движение небесных тел под действием различных факторов описывается дифференциальными уравнениями (ДУ). Выбор вида уравнений играет важную роль, так как от них зависит эффективность интегрирования. В данной работе использовались уравнения движения возмущенной и невозмущенной задачи двух тел в прямоугольной гелиоцентрической системе координат, отнесенной к экватору J2000.0. Возмущенная задача двух тел описывается системой уравнений вида , (1) где P - возмущающее ускорение от других тел: , (2) x = {x1, x2, x3} - вектор положения; - расстояние от Солнца до астероида; - гравитационный параметр; - постоянная Гаусса; и m - массы Солнца и исследуемого объекта соответственно. Дифференциальное уравнение второго порядка (1) можно представить в виде системы уравнений первого порядка (3) В член P входят возмущения от влияния планет Солнечной системы, Плутона и Луны; Mi - масса возмущающего тела; xi - вектор положения возмущающего тела; Δi - расстояние между возмущающим телом и исследуемым объектом. Координаты, скорости и массы планет, Плутона и Луны считывались из фонда координат больших планет DE431 [5]. При P = 0 получаем невозмущенную задачу двух тел. В работе также исследовалась смешанная система дифференциальных уравнений с параметром хаотичности MEGNO [4]. В обобщенном виде уравнения движения (1) можно записать так: (4) где - вектор состояния исследуемого астероида в фазовом пространстве. Тогда зададим вектор δ как начальное отклонение вектора состояния q, который будет содержать в себе вариации шести параметров в фазовом пространстве . Эволюция вектора δ с точностью до бесконечно малых первого порядка описывается вариационным уравнением , (5) где - матрица Якоби системы дифференциальных уравнений (4). Совместно интегрируя с уравнениями движения (4) и уравнениями в вариациях (5), можно получить вспомогательные величины y и w . (6) С помощью этих величин определяется параметр MEGNO и его усредненная величина как . (7) Эволюция параметра во времени позволяет точно разделить регулярный и хаотический режимы движения. В итоге общая смешанная система ДУ первого и второго порядков имеет вид (8) 2. Интегратор Гаусса - Эверхарта Интегратор Гаусса - Эверхарта [1] является модификацией интегратора Эверхарта [2]. Классический (RA15), а также и модифицированный (RADAU_27) программный код Эверхарта значительно ограничивают возможности интегратора, поэтому у В.А. Авдюшева возникла идея подвергнуть дальнейшей модификации метод Эверхарта [1]. Благодаря программной реализации на языке Fortran программный код стал меньше, читабельнее и проще для понимания. Были исключены константы, связанные с порядком метода, а также исправлен алгоритм выбора переменного шага, теперь его величина определяется порядком интегратора. Кроме того, в новой версии откорректирован выбор начального шага, теперь он осуществляется по оценкам интегрирующей схемы второго порядка с учетом поведения правых частей уравнений, появилась возможность вести интегрирование на шаге до полной сходимости итерационного процесса и др. 3. Интегратор Lobbie Коллокационный интегратор Lobbie, основанный на разбиении Лобатто, предназначен для численного решения смешанных систем первого и второго порядков. Его прототипом стал широко используемый в динамической астрономии интегратор Эверхарта [2]. Примечательной особенностью коллокационных интеграторов является то, что их теоретическая основа, как и программная реализация, универсальна для любого порядка [6]. Практически порядок определяется разбиением на шаге, а именно количеством и спецификой распределения узловых значений, через которые выражаются все остальные константы интегратора. Кроме того, на разбиениях гауссовых квадратур Лежандра и Лобатто коллокационные интеграторы становятся геометрическими [7]: симметричными и орбитально устойчивыми, а на разбиении Лежандра - еще и симплектическими. Следует также отметить, что, в отличие от других интеграторов, коллокационные позволяют на каждом шаге легко конструировать приближенное аналитическое решение, чем удобно пользоваться для частого вывода результатов на плотной временной сетке. Алгоритм интегратора Lobbie подробно описан в [3]. 4. Методика исследования Для тестирования интеграторов были выбраны следующие объекты Солнечной системы с различными значениями эксцентриситетов их орбит: (4964) Kourovka, (3753) Cruithne, (3200) Phaethon. В таблице для перечисленных астероидов представлены кеплеровы элементы: большая полуось a, эксцентриситет e, наклонение i, долгота восходящего узла Ω, аргумент перицентра ω, средняя аномалия M и период обращения T. Параметры орбит взяты с сайта Центра Малых Планет MPC (Minor Planet Center) на момент времени t. Орбитальные характеристики астероидов Объекты Phaethon Cruithne Kourovka a, a.e. 1.271348 0.997716 2.2624391 e 0.889789 0.514934 0.1233911 i, град 22.2571 19.80691 4.91924 Ω, град 265.22 126.22 167.40334 ω, град 322.1801 43.84422 204.864 M, град 143.9764 202.0643 66.34177 T, годы 1.43 1 3.4 t, JD 2459400.5 2459400.5 2459400.5 Из таблицы можно видеть, что Kourovka имеет почти круговую орбиту, является астероидом главного пояса и обладает самым большим периодом из рассматриваемых объектов. Более вытянутая орбита наблюдается у астероидов Cruithne и Phaethon. В ходе эксперимента сравнивалась эффективность интеграторов Гаусса - Эверхарта и Lobbie 10-го порядка на интервале времени 100 орбитальных периодов астероида. Так как интегратор Lobbie способен решать системы первого и второго порядков, то имеются две его реализации: Lobbie(I) и Lobbie(II). Для каждого случая задавался параметр локальной точности, который варьировался от 10-2 до 10-20. Выходными результатами являются параметры быстродействия и точности. На всех графиках, которые будут представлены далее в данной работе, оси абсцисс соответствуют значения параметра быстродействия nf, а оси ординат - точности интегрирования r. Под быстродействием понимается число обращений к процедуре правых частей, которые фактически несут в себе основной физический смысл и на вычисление которых тратится основное время. Так как интегратор Lobbie из-за своих технологических особенностей является симметричным, то использование метода прямого и обратного интегрирования не является достоверным критерием для оценивания точности. Поэтому оценка точности интегрирования определялась путем сравнения с эталонными орбитами, которые были посчитаны интеграторами более высокого порядка: , (9) где и - i-координата для эталонной и расчетной орбиты соответственно. 5. Эффективность интеграторов на примере решения невозмущенной задачи двух тел Интеграторы были протестированы на примере численного моделирования орбит исследуемых объектов. Результаты в виде диаграмм «точность - быстродействие» представлены на рис. 1. Рис. 1. Диаграммы «точность - быстродействие» интеграторов Lobbie(I, II) и Гаусса - Эверхарта для астероида Kourovka (e = 0.12) (a), Cruithne (e = 0.51) (б) и Phaethon (e = 0.89) (в) на примере решения задачи двух тел Видно, что графики для Lobbie(I) и Lobbie(II) практически параллельны, а кривая для Gauss15 начинается около кривой интегратора Lobbie(I) и постепенно сходится к Lobbie(II). При достаточно малом параметре задаваемой точности эффективность интегратора с системой уравнений второго порядка выше в 2 раза, чем первого. В среднем при более высокой локальной погрешности эффективность лучше примерно в 3-4 раза. Наибольшее преимущество наблюдается в самом широком месте между графиками и соответствует а.е., быстродействие при этом улучшается в 6 раз. Однако, если провести вертикальную линию в правом конце графиков, то можно увидеть, что при одинаковом количестве обращений к функции правых частей алгоритм Гаусса - Эверхарта дает на два порядка более точный результат, чем Lobbie(II), и на шесть порядков лучше, чем схема Lobbie(I). Два других объекта имеют более вытянутые орбиты и при численном моделировании требуют больше ресурсов. При прохождении астероида вблизи перицентра влияние центрального тела существенно увеличивается и шаг интегрирования начинает уменьшаться. Поэтому общее количество обращений к функции правых частей возросло. Для астероида Phaethon переломный момент в интеграторе Lobbie наступил гораздо раньше при а.е. Анализируя все три графика, можно сделать вывод, что интегратор Lobbie обладает преимуществом в силу возможности интегрировать уравнения второго порядка. Но при высокой точности он уступает интегратору Гаусса - Эверхарта. Наклон оси графиков интегратора Гаусса - Эверхарта является более крутым, что дает ему большее преимущество при более высокой точности. Но ее использование зачастую не совсем целесообразно, так как ошибка в начальных данных оказывается не меньше, чем 10-9 а.е. Таким образом, в диапазоне значений от 10-2 до 10-15 а.е. явное преимущество у схемы Lobbie(II). Однако использование более высокой точности не является бессмысленным и может помочь при интегрировании на длинные промежутки времени. Важно отметить то, что в данном эксперименте использовались три итерации (NI = 3). При таком значении в интеграторах могут проявляться геометрические свойства. Увидеть это можно при исследовании поведения интегралов, в частности интеграла энергии. Известно, что полная энергия системы есть величина неизменная, тогда замеряя ее начальное значение, можно увидеть характер ее изменения при интегрировании. Эксперимент проводился с постоянным шагом h = T/16 для интеграторов 8-го порядка. Результаты приведены на рис. 2. Из данного графика видно, что при использовании в Lobbie двух итераций идет линейный рост ошибки, а при трех итерациях поведение ошибки устойчиво. Для интегратора Гаусса - Эверхарта при использовании четырех итераций идет резкий линейный рост, а геометрические свойства проявляются при использовании пяти итераций. Это благотворно влияет на поведение решения и помогает удерживать его вблизи орбиты. В данном случае сохранность энергии обеспечивает сохранение периода и большой полуоси [3]. Это наглядно показывает, что итерационная сходимость в схеме Lobbie лучше. Рис. 2. Отклонение в вычисляемой энергии H для астероида Kourovka (e = 0.12) 6. Эффективность интеграторов на примере решения возмущенной задачи двух тел На данном этапе исследовалось поведение интеграторов Lobbie и Гаусса - Эверхарта в рамках возмущенной задачи двух тел. Возмущающими факторами являются планеты Солнечной системы, Луна и Плутон. Результаты исследования представлены на рис. 3. Из результатов исследования движения астероида Kourovka (рис. 3, а) видно, что все три алгоритма имеют приблизительно одинаковую эффективность и находятся вблизи друг друга. Графики стали менее гладкими, угол наклона уменьшился, что тоже говорит о снижении эффективности. Более наглядно преимущество демонстрируют интеграторы Lobbie второго порядка и Гаусса - Эверхарта для астероида Cruithne (рис. 3, б). Однако для более вытянутой орбиты астероида Phaethon схема Lobbie(II) на точности от 10-6 до 10-12 уступает схемам первого порядка (рис. 3, в). Рис. 3. Диаграммы «точность - быстродействие» интеграторов Lobbie(I, II) и Гаусса - Эверхарта для астероида Kourovka (e = 0.12) (а), Cruithne (e = 0.51) (б) и Phaethon (e = 0.89) (в) на примере решения возмущенной задачи двух тел Для всех астероидов три интегратора достигают своего предела точности при 10-18-10-19 а.е. Это обусловлено влиянием ошибок, накапливаемых в возмущающем члене, для вычисления которого необходимы массы и координаты планет. Значения этих величин считывались из фонда координат DE431. Данные в фонде имеют двойную точность, что и налагает ограничения в нашем исследовании. На рис. 3, в для кривой Gauss15 также заметно ухудшение эффективности при достижении точности 10-20 а.е., что обуславливается нарастанием ошибок округления. Таким образом, в ходе исследования возмущенной задачи двух тел не было получено значительного преимущества интегратора Lobbie. 7. Эффективность интеграторов для решения смешанных систем В качестве объекта исследования выступает система смешанных дифференциальных уравнений (8) с учетом возмущений от планет, Плутона и Луны. Интегратор Гаусса - Эверхарта не имеет возможности решать такие системы, поэтому система приводится к системе уравнений первого порядка. Результаты исследования эффективности представлены на рис. 4. Согласно работе [3], преимущество интегратора Lobbie должно сохраняться. Однако наши результаты показывают (рис. 4), что ожидаемого улучшения не удалось получить. Некоторый выигрыш имеется для астероида Cruithne. Причиной пониженной эффективности является вид добавочных уравнений в вариациях. Общее количество переменных увеличилось на 5, и вклад от их вычисления с заданной точностью стал существенным. Поэтому шаг начал сильно уменьшаться, что привело к снижению эффективности. Было принято решение задать точность только по части переменных, а именно по координатам, что дало возможность повысить быстродействие более чем в 2 раза (график Lobbie coor). Кроме того, можно предположить, что эффективность должна повыситься, если использовать при вычислении одну итерацию. При таких условиях вклад от уравнений первого порядка не успевает накапливаться, что дает заметное преимущество интегратору Lobbie. Результаты, представленные на рис. 5, подтверждают сделанное предположение - эффективность использования интегратора Lobbie для всех трех объектов повысилась на 3-4 порядка. Рис. 4. Диаграммы «точность - быстродействие» интеграторов Lobbie и Гаусса - Эверхарта для астероидов Kourovka (e = 0.12) (а), Cruithne (e = 0.51) (б) и Phaethon (e = 0.89) (в) на примере решения смешанных систем Рис. 5. Диаграммы «точность - быстродействие» интеграторов Lobbie и Гаусса - Эверхарта для астероидов Kourovka (e = 0.12) (а), Cruithne (e = 0.51) (б) и Phaethon (e = 0.89) (в) на примере решения смешанных систем с параметром NI = 1 Таким образом, решающим фактором в данной части исследования оказался вид уравнений в вариациях. Стоит заметить, что для уравнения (8) в силу наличия операции деления важна не столько абсолютная погрешность, сколько относительная. Повысить эффективность интегрирования возможно также путем введения в систему стабилизирующих уравнений или регуляризирующих преобразований. Заключение Таким образом, в ходе данного исследования изучена эффективность использования интеграторов Гаусса - Эверхарта и Lobbie в задачах астероидной динамики. Были получены следующие результаты: 1. При решении задачи двух тел интегратор Lobbie имеет значительное преимущество за счет уравнений второго порядка. 2. Оба интегратора обладают геометрическими свойствами, но для интегратора Гаусса - Эверхарта требуется большее число итераций. 3. Не было обнаружено существенного преимущества какого-либо из интеграторов при решении возмущенной задачи двух тел. 4. Интегратор Lobbie показал более высокую эффективность в решении смешанных систем ДУ, чем интегратор Гаусса - Эверхарта. Обобщая, можно сделать вывод, что оба интегратора хорошо подходят для решения задач астероидной динамики. Однако более универсальным и эффективным является интегратор Lobbie.

Скачать электронную версию публикации