English version of the page

Численное моделирование движения малых тел Солнечной системы

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

Прогнозирование движения малого тела
Улучшение орбиты малого тела
Преобразование координат малого тела
Преобразование даты

Особенностью сайта является возможность использования различных планетных эфемерид при учете возмущений, а именно, эфемерид DE405, DE408, DE414, DE421, DE422, DE423, DE424, DE425, DE430, DE431, DE432, DE433, DE434, DE435, DE436, DE437, DE438, DE440, DE441, [1] EПM2011, ЕПМ2015 [2]. Кроме того, сайт позволяет выполнять преобразование параметров движения малого тела от одного вида к другому (например, прямоугольные координаты и скорости в оскулирующие элементы орбиты), перевычислять эти параметры на другую эпоху по формулам невозмущенного движения, оценивать точность численного интегрирования уравнений движения путем сравнения результатов прямого и обратного прогнозирования, выполнять преобразование даты из одной шкалы времени в другую. Кроме того, для удобства пользователей на сайте предусмотрена возможность выбора как способа выдачи, так и формата для всех рассчитываемых величин.

Ниже приводится описание веб-интерфейса к основным программам, расположенным на сайте

Прогнозирование движения малого тела

Основное назначение данной программы состоит в вычислении различных параметров движения малого тела на заданные моменты времени по известным начальным элементам орбиты. Для использования программы следует задать начальные параметры движения в поле ввода c таким же названием. Формат задания произволен - можно задать все параметры в одной строке, можно в отдельных строках или как-то еще. Необходимо только оставлять между ними хотя бы один пробел, причем символ перехода на новую строку пробелом не считается, поэтому если элементы находятся в разных строках, то впереди или сзади каждого параметра должен находиться хотя бы один пробел. Кроме того, набор кеплеровых элементов орбиты и порядок их следования должен строго соответствовать описанию в колонке "Способ задания". Если элементы орбиты, которые вы хотите использовать, имеют другой состав, то воспользуйтесь вспомогательной программой "Преобразование координат малого тела". Эта программа позволяет преобразовывать произвольный набор кеплеровых элементов орбиты в любой другой, а также делать преобразования элементов в прямоугольные координаты и компоненты скорости и наоборот, а также перевычислять их на другую эпоху по формулам невозмущенного движения. В качестве задаваемых начальных параметров могут быть использованы также прямоугольные координаты и компоненты скорости малого тела. Следует указать с помощью переключателя вверху какой набор начальных данных вы хотите использовать.

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

Ниже начальных параметров задается масса рассматриваемого тела (в единицах массы Солнца, равной 1.989 · 1030 кг). Если тело небольшое, то можно задать для него массу, равную нулю.

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

Здесь же расположена кнопка "Прямое и обратное интегрирование", которая запускает процесс интегрирования сначала в прямом направлении - от начального момента эфемериды до конечного (или наоборот, если начальная эпоха находится после интервала прогнозирования), а затем в обратном направлении с запоминанием при счете прямоугольных координат и скоростей на моменты выдачи. При обратном интегрировании в моменты выдачи выполняется вычисление модулей разности вектора положения (dr, а.е.) и скорости (dv, а.е./сут), полученных при прямом и обратном интегрировании, которые и выдаются в новое окно. Если начальная эпоха расположена между начальным и конечным моментами эфемериды, то описанный процесс выполняется в два этапа: сначала вперед до конечного момента эфемериды и обратно, а затем назад до начального момента эфемериды и обратно.
Моменты времени для величин dr и dv выдаются тем же способом и в том же формате, что и при вычислении эфемериды. Т.е., чтобы задать способ и формат выдачи моментов времени, надо прокрутить страницу вниз и сделать это в последней строке таблицы "Что выдавать в эфемериде".

Далее указывается с помощью переключателя основная плоскость (экватор 2000.0 или эклиптика 2000.0) системы координат, в которой заданы начальные параметры.

Следующие два выключателя позволяют учитывать в модели релятивистские и негравитационные эффекты. Релятивистский эффект вносит заметный вклад в движение малых тел, если они имеют достаточно малые перигелийные расстояния (< 0.5 а.е.) Негравитационные силы учитываются только для комет. Для учета негравитационных возмущений следует включить соответствующий выключатель и ввести в графы значения параметров негравитационных сил A1 и A2.

Далее следует раздел учета возмущений от планет и Луны. С помощью включения и выключения учета возмущений можно исследовать их структуру, т.е. вклад притяжения какой-либо планеты или ряда планет в изменение элементов орбиты малого тела. При учете этих возмущений координаты планет и Луны интерполируются с помощью полиномов Чебышева, коэффициенты при которых извлекаются из эфемерид, выбираемых пользователем (DE405,..., DE441, ЕПМ2011, ЕПМ2015). В скобках рядом с названием каждых эфемерид указан охватываемый ими интервал времени в виде пары юлианских дат.

Здесь также расположен выключатель, который позволяет использовать вариант эфемериды, сглаженный на "стыках" смежных интервалов интерполирования таким образом, что на этих стыках сохраняется непрерывность интеполируемых координат планет и их первых производных с точностью до 33-34 десятичных знаков. Алгоритм такого сглаживания и результаты его применения описаны, например, в работах [3, 4]. Сглаженные варианты эфемерид имеет смысл применять лишь при расчетах с четверной точностью, т.е. при использовании разрядной сетки 128 бит, соответствующей десятичной разрядности в 34 значащих цифры.

Кроме того, здесь же расположен выключатель, при включении которого шаг интегрирования будет корректироваться в соответствии со "стыками" интервалов интерполирования в используемых эфемеридах. Эта корректировка заключается в том, что если в пределы какого-то шага интегрирования попадает стык смежных интервалов интерполирования, то шаг на этом стыке заканчивается. Такая корректировка [3, 4] позволяет избавиться от влияния на точность численного интегрирования скачков производных второго и более высокого порядка от интерполируемых координат, которые имеют место стыках. Заметим, что корректировка шага в реализованной здесь вычислительной программе может применяться лишь при выборе переменного шага интегрирования (т.е. при положительном параметре LL) и отсутствии расчета сближений с планетами.

Далее задаются начальный шаг интегрирования (в сутках), порядок метода Эверхарта [5] численного интегрирования обыкновенных дифференциальных уравнений и параметр LL, определяющий поведение шага интегрирования: при задании отрицательного значения LL шаг будет постоянным, равным начальному шагу интегрирования; при задании положительного значения LL шаг будет переменным и его величина будет определяться из условия, чтобы локальная погрешность не превышала величины 10-LL. Но, как показывает практика использования метода Эверхарта, его реальная точность заметно превышает эту величину, поэтому оптимальное значение LL обычно составляет примерно две трети от порядка метода интегрирования. Что касается самого порядка, то эмпирически установлено, что его оптимальное значение на единицу меньше используемой разрядной сетки, т.е., например, при использовании 16-разрядной сетки (64 бита) это 15-й порядок и т.д. Однако, поскольку в реализованной программе возможен выбор только некоторых значений порядка (от 7 до 39 через 4), то для двух других разрядных сеток (19-ти и 34-разрядной) в качестве оптимальных значений используются 19-й и 27-й порядок.

После параметров интегратора задается размер разрядной сетки, точнее, размер чисел с плавающей точкой, используемых при выполнении вычислений. Возможности выбора этого размера невелики и целиком определяются возможностями используемого компилятора. Так как все программы, выполняющие на сайте расчеты, написаны на языке Фортран, современные компиляторы которого поддерживают только четыре размера чисел с плавающей точкой (32, 64, 80 и 128 бит), то и были использованы три последних размера. Размер чисел в 32 бита (так называемые числа одинарной точности, содержащие в десятичном выражении 7 значащих цифр) не использован, поскольку такие числа не обеспечивают необходимой точности астрономических вычислений.

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

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

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

Улучшение орбиты малого тела

Эта программа позволяет выполнять улучшение орбит астероидов и комет методом наименьших квадратов по данным всех видов наблюдений, выполняемых в настоящее время: позиционных (стационарных, спутниковых и роверных) и радарных. Наблюдения вводятся в основном окне интерфейса, расположенном с левой стороны и имеющем наибольший размер. Наблюдения должны иметь формат, принятый в Центре малых планет. Чтобы получить наблюдения какого-либо астероида или кометы, следует открыть на этом сайте страницу поиска и ввести в строке поиска какое-либо обозначение или название объекта, например, 2024 YR4 или Tomskuniver. Можно также вводить номер астероида или обозначение в сокращенном 7-символьном виде, принятом в Центре малых планет. На открывшейся (после нажатия кнопки Show) странице с результатами поиска следует нажать кнопку Download (прокрутив страницу немного вниз) и скопировать из открывшегося окна наблюдения объекта, для чего удобно использовать сочетание клавиш Ctrl-A (выделить все) и Ctrl-C или Сtrl-Ins (копировать). Затем следует очистить окно интерфейса от содержащихся там тестовых наблюдений, нажав в нем Ctrl-A и Del, и вставить туда из буфера наблюдения объекта (нажав Ctrl-V или Shift-Ins).

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

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

Ниже расположены переключатели, позволяющие учитывать релятивистский и негравитационный эффекты. При учете негравитационного эффекта, следует задать параметры A1 и A2 негравитационных сил, представляющие собой радиальную и трансверсальную составляющие негравитационного ускорения на расстоянии 1 а.е. от Солнца. Учитывать неграфитационные эффекты имеет смысл только для комет. Возможно также улучшение параметров A1 и A2 вместе с координатами и скоростями объекта (задается переключателем).

Крайние справа переключатели позволяют выбрать эфемериды больших планет и Луны, используемые при учете возмущений. По умолчанию используются последние выпущенные эфемериды DE440. В скобках указан охватываемый эфемеридами интервал времени в виде пары юлианских дат. Как и в программе прогнозирования движения малого тела, можно использовать сглаженный на "стыках" интервалов интерполирования вариант эфемерид, что имеет смысл лишь при расчетах с четверной точностью (на 128-битовой разрядной сетке).

Ниже расположен переключатель используемой при вычислениях разрядной сетки, точнее размера в битах чисел с плавающей точкой: 64 бита - числа двойной точности, 80 бит - расширенной точности и 128 бит - четверной точности. При переключении разрядной секи автоматически изменяется принятая по умолчангию точность сходимости итерационного процесса: 10-9 а.е., 10-12 а.е. и 10-20 а.е. соответственно. Такие значения получены эмпирически и пользователь может изменять их по своему усмотрению.

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

Ссылки

1. ftp://ssd.jpl.nasa.gov/pub/eph/planets
2. ftp://quasar.ipa.nw.ru/incoming/EPM
3. Батурин А.П. Влияние способов построения эфемерид больших планет и Луны на точность прогнозирования движения астероидов // Изв. вузов. Физика. 2014. Т. 57. № 9. С.72–79.
4. Батурин А.П. Повышение точности численного интегрирования уравнений движения астероидов при учете возмущений от больших планет и Луны из эфемерид DE // Астрон. вестн. 2018. Т. 52. № 4. С. 360–363.
5. Everhart E. An efficient integrator that uses Gauss-Radau spacings // Dynamics of comets: their origin and evolution // Proc. 83rd IAU Colloq. Rome, 11–15 June 1984 / Eds Carusi A., Valsecchi G.B. Dordrecht: D. Reidel Publ. Co., 1985. P. 185–202.

Автор, веб-мастер: А.П.Батурин (ТГУ, Томск)