Космические исследования, 2020, T. 58, № 1, стр. 27-39

Исследование характеристик вероятной области соударения астероида Апофис с Землей в 2036 г.

В. В. Ивашкин 1 2*, П. Гуо 1, К. А. Стихно 2

1 Институт прикладной математики им. М.В. Келдыша РАН
г. Москва, Россия

2 Московский государственный технический университет им. Н.Э. Баумана
г. Москва, Россия

* E-mail: ivashkin@keldysh.ru

Поступила в редакцию 26.12.2018
После доработки 26.03.2019
Принята к публикации 25.04.2019

Полный текст (PDF)

Аннотация

В настоящей работе определено и исследовано возможное место соударения опасного астероида Апофис с Землей в 2036 г. Разработаны алгоритмы для поиска всего вероятного множества попадающих в Землю траекторий астероида Апофис с помощью методов Монте-Карло и градиентного спуска. Для данного множества произведен поиск нескольких вероятных семейств столкновительных с Землей в 2036 г. траекторий астероида, соответствующих подмножествам значений перигейного расстояния орбиты Апофиса от минимально возможного расстояния, ~2069 км, до максимального, соответствующего касанию поверхности Земли, ~6376 км. Получено отображение всего множества попадающих траекторий на поверхность Земли, а также построена на карте мира столкновительная полоса. Изучены особенности этой зоны соударения и свойства столкновительных траекторий Апофиса. Дано сравнение результатов данного исследования с другими работами.

1. ВВЕДЕНИЕ

Астероид (99942) Апофис (лат. Apophis) потенциально угрожает Земле в XXI в. С момента его открытия в 2004 г. теоретическая вероятность его столкновения с Землей в 2029 г. и связанные с резонансными возвратами к Земле возможные его соударения с нашей планетой многократно изучались отечественными и зарубежными учеными [116]. По результатам обработки наблюдений в 2005 г. было выявлено, что астероид Апофис тесно сближается с Землей в 2029 г. на расстояние ~38 ± 2.3 тыс. км от центра Земли и имеет некоторую малую вероятность его столкновения с Землей при следующем сближении в 2036 г. Однако эта вероятность практически исключена NASA после существенного уточнения его орбиты из новых наблюдений в 2012–2013 гг., когда Апофис опять сближался с Землей [8, 13]. Несмотря на это, изучение и выявление характеристик возможной зоны падения астероида Апофис на Землю в 2036 г. важны как сами по себе, для анализа структуры и характеристик попадающих траекторий астероида, так и для разработки методик и алгоритмов, которые могут понадобиться в любой момент при тесном сближении с Землей других опасных небесных тел. Поэтому исследование характеристик зоны падения Апофиса представляется важной и актуальной задачей в комплексной проблеме астероидно-кометной опасности.

В настоящей работе определена и исследована вероятная область соударения опасного астероида Апофис с Землей в 2036 г., а также проанализированы характеристики этой области падения астероида. Задача выполнялась в два этапа. На первом этапе работы разработаны алгоритмы для поиска траекторий опасного астероида Апофис, попадающих в Землю в 2036 г. с учетом влияния реального гравитационного поля и давления Солнечного света. Произведен поиск нескольких семейств попадающих в Землю траекторий Апофиса, для каждого из которых перигейное расстояние траектории близко к некоторому фиксированному перигейному расстоянию из диапазона rπ = {2069–2100; 2500 ± 1; 3000 ± 1; …; 6360 ± 1; 6375 ± 1} км – от минимального, ~2069 км, до максимального, ~6376 км. Далее, на втором этапе, определено множество точек пересечения этих траекторий с поверхностью Земли и на основе этого построена столкновительная полоса на карте мира. Исследованы геометрические, временные и энергетические характеристики столкновительных траекторий Апофиса и полосы точек соударения астероида с Землей, в частности, ее ширина и длина, ориентация в пространстве, углы входа астероида в атмосферу, время и скорость столкновения, направление подлета астероида к Земле. Выполнено сравнение с результатами других авторов по данной проблеме. Данная работа развивает [6].

2. МЕТОДИКА ОПРЕДЕЛЕНИЯ МНОЖЕСТВА ПОПАДАЮЩИХ ТРАЕКТОРИЙ АПОФИСА

Для поиска попадающих в Землю траекторий астероида Апофис использованы методы Монте-Карло и градиентного спуска. При этом траектории Апофиса определялись численным интегрированием уравнений движения астероида с учетом возмущений от реального гравитационного поля и давления солнечного света при варьировании начальных данных в допустимом множестве [16]. Разработаны также алгоритмы для поиска множества попадающих в Землю траекторий Апофиса с заданными значениями перигейного расстояния.

2.1. Модель движения астероида. Чтобы повысить точность определения траекторий астероида Апофис, была использована достаточно полная модель его движения. Применены уравнения движения астероида в гелиоцентрической геоэкваториальной прямоугольной невращающейся системе координат (СК) OXYZ на эпоху J2000 с учетом возмущений от притяжения всех больших планет и Луны, как точечных тел, сжатия Земли (с точностью до второй зональной гармоники J2), а также давления солнечного света:

(1)
$\begin{gathered} \frac{{{{d}^{2}}{{{\mathbf{r}}}_{{\text{A}}}}}}{{d{{t}^{2}}}} = - \frac{{{{\mu }_{{\text{S}}}} + {{\mu }_{{\text{A}}}}}}{{r_{{\text{A}}}^{3}}}{{{\mathbf{r}}}_{{\text{A}}}} + \sum\limits_i {{{\mu }_{i}}\left( {\frac{{{{{\mathbf{r}}}_{i}} - {{{\mathbf{r}}}_{{\text{A}}}}}}{{{{{\left| {{{{\mathbf{r}}}_{i}} - {{{\mathbf{r}}}_{{\text{A}}}}} \right|}}^{3}}}} - \frac{{{{{\mathbf{r}}}_{i}}}}{{{{{\left| {{{{\mathbf{r}}}_{i}}} \right|}}^{3}}}}} \right)} + \\ + \,\,{{\Delta }_{{\text{E}}}} + {{{\mathbf{a}}}_{{SP}}}, \\ \end{gathered} $
где t – время, rA – гелиоцентрический радиус-вектор астероида, μS и μA – гравитационные параметры Солнца и астероида, ΔE – возмущающее ускорение из-за полярного сжатия Земли, aSP – возмущающее ускорение из-за давления солнечного света, ri, μi – гелиоцентрические радиус-векторы и гравитационные параметры больших планет и Луны, их координаты определяются по JPL-эфемеридам DE405.

В силу большого объема вычислений интегрирование выполнено методом предиктор-корректора 8-ого порядка, эффективным по быстродействию, разработанным и используемым в ИПМ им. М.В. Келдыша РАН [17]. В качестве номинальных начальных данных r0n, v0n для уравнений движения на начальный момент t0 в эпоху JD2453400.5 (30.I.2005 г.) и их возможных отклонений были использованы результаты ИПА РАН по обработке оптических и радиолокационных наблюдений астероида [1], причем:

(2)
$3{{\sigma }_{0}}\left( {x,y,z} \right) = 3\,{\text{км}},\,\,\,\,3{{\sigma }_{0}}\left( {{{{v}}_{x}},{{{v}}_{y}},{{{v}}_{z}}} \right) = 2{\text{ мм}}/{\text{с}},$
где σ0(x, y, z), σ0(vx, vy, vz) – среднеквадратичные отклонения (СКО) по компонентам начальных радиус-вектора r0 и вектора скорости v0 астероида.

2.2. Поиск попадающих в Землю траекторий астероида Апофис. Чтобы как можно более полно и за приемлемое время оценить множество опасных траекторий, не перебирая при этом с очень мелким малом все множество начальных состояний, был применен комплексный метод их поиска.

Распределение начальных данных для траекторий из эллипсоида рассеивания (2) было нами принято нормальным с диагональной ковариационной матрицей. В соответствии с этим, случайным образом смоделировано ~104 начальных векторов состояния астероида z = (r0, v0). При этом нормальное распределение моделировалось на основе равномерно распределенных случайных величин на отрезке [0, 1], согласно [18].

Каждая траектория была проинтегрирована до сближения с Землей в 2036 г. Затем из общего множества траекторий было выделено подмножество из около 50 траекторий с перигейным расстоянием в 2036 г. вблизи сферы действия Земли.

Для этих траекторий их перигейное расстояние в 2036 г. рассматривалось как минимизируемый параметр, зависящий от начального состояния астероида z = (r0, v0):

(3)
${{r}_{\pi }}({\mathbf{z}},{{t}_{k}},\Delta t) = \mathop {\min }\limits_{t \in ({{t}_{k}} - \Delta t,{{t}_{k}} + \Delta t)} \left| {{{{\mathbf{r}}}_{{AE}}}(t)} \right|,$
здесь rAE – геоцентрический радиус-вектор астероида в окрестности (±Δt) момента tk его тесного сближения с Землей. Далее методом градиентного спуска было получено соответствующее количество траекторий с перигейным расстоянием меньше радиуса Земли. При этом градиент функции rπ(z) численно определялся методом конечных разностей. Далее делается шаг в направлении антиградиента:
(4)
$\tilde {z} = {\mathbf{z}} - {{\vartheta }_{k}} \cdot {\text{grad}}({{r}_{\pi }}({\mathbf{z}})),$
где ${{\vartheta }_{k}}$– коэффициент коррекции шага. На основании интегрирования новой траектории и определения ее перигейного расстояния принималось решение о продолжении или остановке процесса минимизации. Блок-схема алгоритма итерационного процесса приведена на рис. 1.

Рис. 1.

Блок-схема алгоритма поиска одной попадающей траектории.

Таким образом было получено соответствующее количество (~50) траекторий с перигейным расстоянием меньше радиуса Земли [25]. Для начальных состояний найденных траекторий был посчитан диапазон интеграла энергии h0 задачи двух тел гелиоцентрического движения астероида:

(5)
${{h}_{0}} = {v}_{0}^{2} - 2{{{{\mu }_{S}}} \mathord{\left/ {\vphantom {{{{\mu }_{S}}} {{{r}_{0}}}}} \right. \kern-0em} {{{r}_{0}}}},$
где r0 и v0 – начальные гелиоцентрическое расстояние и скорость астероида.

Затем снова методом Монте-Карло проводился перебор достаточно большого множества траекторий (порядка 107) из начального эллипсоида (2) и выбирались траектории с величиной интеграла энергии в начальный момент времени, близкой к его значениям в указанном диапазоне для определенных выше столкновительных траекторий. Последующее поверочное интегрирование до 2036 г. проводилось только в том случае, когда интеграл энергии исследуемой траектории не выходил далеко за рамки диапазона, задаваемого ранее найденными попадающими траекториями. Так было найдено уже ~5000 траекторий с пролетом внутри сферы действия Земли. Для получившегося множества траекторий снова варьированием начального состояния и минимизацией перигейного расстояния в 2036 г. было сгенерировано соответствующее множество попадающих траекторий.

Найденные попадающие в Землю в 2036 г. траектории Апофиса были значительно разбросаны вдоль всего эллипсоида рассеивания начальных данных. На рис. 2 приведены “сферически-симметричное” множество исследовавшихся точек, соответствующих исходному множеству (2), и множество точек, образующих “отрезок прямой” и соответствующих попадающим в Землю орбитам астероида. Здесь dr = |r0r0n|, dv = |v0v0n|, где r0n, v0n соответствуют номинальной траектории [1], r0, v0 соответствуют траектории в эллипсоиде (2).

Рис. 2.

Множество начальных данных и множество точек, соответствующих начальным данным для столкновительных траекторий (~5000) в проекции на плоскость (dr, dv).

2.3. Детализация множества попадающих траекторий. Главной характеристикой каждой траектории построенного семейства является ее перигейное расстояние в 2036 г. Для всех траекторий этого множества оно варьируется в диапазоне от минимального значения, равного ~2069 км [26] до максимального, соответствующего предельной величине для “касания” траекторией Земли, ~6376 км. Для дальнейшего изучения всего множества попадающих траекторий были построены его подмножества, семейства, объединяющие траектории с одинаковым перигейным расстоянием rπ в 2036 г. для нескольких значений rπ из указанного диапазона. При этом было выявлено, что каждое данное семейство с перигейным расстоянием, большим минимального, состоит из двух несвязных частей, для которых место падения на Землю лежит по разные стороны от места падения для траектории с минимальным перигейным расстоянием.

Эти траектории с заданным значением перигейного расстояния получались из полученных ранее попадающих траекторий как порождающих. Перигейное расстояние rπ в 2036 г. рассматривалось как функция начального состояния астероида rπ(z). Затем, в отличие от [25], выбиралось два направления – вдоль градиента интеграла энергии задачи двух тел и противоположное ему:

(6)
${{{\mathbf{g}}}_{1}} = {\text{gra}}{{{\text{d}}}_{{({{{\mathbf{r}}}_{0}},{{{\mathbf{v}}}_{0}})}}}({{h}_{0}}),\,\,\,\,{{{\mathbf{g}}}_{2}} = - {\text{gra}}{{{\text{d}}}_{{({{{\mathbf{r}}}_{0}},{{{\mathbf{v}}}_{0}})}}}({{h}_{0}}).$

Вдоль этих направлений производился поиск нуля функции fi = rπ(r0, v0) – ci, где ci – требуемое значение rπ на изолинии. Подобный однопараметрический поиск быстро приводил к траектории, удовлетворяющей требуемым свойствам. Затем каждая траектория семейства отображалась на поверхность Земли в момент столкновения в 2036 г., т.е. для нее определялась точка падения (ее координаты) на поверхности Земли.

Минимальное перигейное расстояние rπ множества попадающих траекторий астероида Апофис в 2036 г. оценено, как ~2069 км. Для нашей задачи было получено 11 семейств попадающих траекторий, для каждого из которых значение перигейного расстояния меняется в некотором диапазоне из множества rπ = {2069–2100; 2500 ± 1; 3000 ± 1; …; 6360 ± 1; 6375 ± 1} км. Для каждой попадающей траектории получены векторы состояния на момент прохождения условного перигея в 2036 г.:

(7)
${{t}_{\pi }},{{{\mathbf{r}}}_{\pi }}({{x}_{\pi }},{{y}_{\pi }},{{z}_{\pi }}),{{{\mathbf{v}}}_{\pi }}({{{v}}_{x}}_{\pi },{{{v}}_{y}}_{\pi },{{{v}}_{z}}_{\pi }).$

3. АНАЛИЗ ХАРАКТЕРИСТИК ЗОНЫ СОУДАРЕНИЯ И СТОЛКНОВИТЕЛЬНЫХ ТРАЕКТОРИЙ АСТЕРОИДА АПОФИС

Рассматривается задача отображения множества сталкивающихся траекторий астероида Апофис на поверхность Земли в момент столкновения. Разработаны алгоритмы для определения географических координат точек соударения с Землей. Анализируются характеристики множества точек соударения астероида Апофис, в том числе, зоны падения, угол входа в атмосферу, время и скорость столкновения, а также характеристики столкновительной полосы на Земле. Выполнен анализ влияния прецессии и нутации на определение точки падения. На данном этапе анализа нами не учитывалось влияние атмосферы Земли на точки падения астероида (некоторые вопросы анализа характеристик движения астероида типа Апофис после его входа в атмосферу Земли ранее изучены в работе [16]).

3.1. Определение точек соударения Апофиса с Землей. Интегрируя уравнения движения астероида (1) с начальными данными (7), полученными в п. 2, можем получить на любой момент времени t радиус-вектор rA и вектор скорость vA астероида. Отметим, что в данном разделе применен метод интегрирования Рунге–Кутты 7(8) порядка, с использованием эфемерид DE421. При этом было предварительно выполнено сравнение с методом интегрирования в п. 2. Оно показало хорошую близость результатов этих обоих методов. Отсюда определяем радиус-вектор и вектор скорость астероида относительно центра масс Земли в геоэкваториальной системе координат r(x, y, z), v(vx, vy, vz):

(8)
${\mathbf{r}} = {{{\mathbf{r}}}_{{\text{A}}}} - {{{\mathbf{r}}}_{{\text{E}}}},\,\,\,\,{\mathbf{v}} = {{{\mathbf{v}}}_{{\text{A}}}} - {{{\mathbf{v}}}_{{\text{E}}}},$
где rЕ, vЕ − радиус-вектор и вектор скорость Земли в гелиоцентрической геоэкваториальной системе координат на момент времени t из эфемериды DE421. Далее переводим вектор состояния астероида в гринвичскую вращающуюся систему координат, и находим, в частности, его радиус-вектор rG. На первом этапе этого перехода мы не учитывали влияние прецессии и нутации на определение координат точки падения, в силу малости их влияния. Их влияние учтено далее, в п. 3.5.

Для перевода во вращающуюся систему координат сначала находим среднее Гринвичское звездное время на момент времени t [1923]. Используя его, переходим к Гринвичской вращающейся системе координат, находим радиус-вектор rG, геоцентрические широту φ' и долготу λ астероида. Далее с учетом сжатия Земли вычисляем географическую широту φ и высоту H астероида над поверхностью земного эллипсоида [1925]. В качестве формы Земли принят земной эллипсоид, характеризуемый средним экваториальным радиусом Rэ = 6378.137 км и полярным сжатием αе = 1/298.25.

Считаем, что при H = 0 астероид достигает поверхность Земли. Обозначим соответственно время попадания, географические широту и долготу точки падения через tc, φс и λс. При H = 100 км астероид проходит границу атмосферы. Обозначим соответственно время входа в атмосферу через tатм.

В данной работе для определения времен tатм и tc, соответствующих H = 100 км и H = 0, построен итерационный процесс, использующий интегрирование, при заданной точности по высоте 1 мм. В качестве начального времени интегрирования выбран момент времени прохождения перигея tπ для каждой попадающейся траектории. Интегрируя уравнения движения астероида до времен tатм и tc, получим радиус-векторы и векторы скорости астероида в геоцентрической СК: rатм, vатм и rc, vc, переводим их в географическую СК. Если Ω − угловая скорость вращения Земли, то компоненты скорости астероида vr(vxr, vyr, vzr) относительно поверхности Земли:

(9)
${{{v}}_{{xr}}} = {{{v}}_{x}} + \Omega y;\,\,\,\,{{{v}}_{{yr}}} = {{{v}}_{y}}--\Omega x;\,\,\,\,{{{v}}_{{zr}}} = {{{v}}_{z}}.$

При t = tc, получаем скорость соударения vcr астероида с Землей с учетом вращения Земли. Вводя орты (ξ, η, ζ) по нормали к поверхности от центра Земли, по параллели на Восток, по меридиану на Север, определим соответствующие составляющие относительной скорости vξ, vη, vζ, ее проекцию на поверхность vηζ, угол наклона проекции к параллели ψr. При t = tатм получим относительную скорость астероида при входе в атмосферу Земли vатмr.

(10)
$\begin{gathered} \sin {{\psi }_{r}} = \frac{{{{{v}}_{{r{\zeta }}}}}}{{\sqrt {{v}_{{r{\eta }}}^{2} + {v}_{{r{\zeta }}}^{2}} }}, \\ cos{{\psi }_{r}} = \frac{{{{{v}}_{{r{\eta }}}}}}{{\sqrt {{v}_{{r{\eta }}}^{2} + {v}_{{r{\zeta }}}^{2}} }},\,\,\,\,0 \leqslant {{\psi }_{r}} < 2\pi . \\ \end{gathered} $

Для всех полученных траекторий соударения определяем указанные характеристики соударения, в частности, время и скорость tc, vc = |vc|, vcr, географические широту и долготу точки падения φс и λс, углы входа в атмосферу θa, θr:

(11)
$\begin{gathered} {{\theta }_{а}} = \frac{\pi }{2} - \arccos \left( {\frac{{{{{\mathbf{r}}}_{{{\text{атм}}}}} \cdot {{{\mathbf{v}}}_{{{\text{атм}}}}}}}{{{{r}_{{{\text{атм}}}}} \cdot {{{v}}_{{{\text{атм}}}}}}}} \right), \\ {{\theta }_{r}} = \frac{\pi }{2} - \arccos \left( {\frac{{{{{\mathbf{r}}}_{{{\text{атм}}}}} \cdot {{{\mathbf{v}}}_{{{\text{атм}}\;r}}}}}{{{{r}_{{{\text{атм}}}}} \cdot {{{v}}_{{{\text{атм}}\;r}}}}}} \right), \\ \end{gathered} $
здесь углы θа и θr отрицательны, когда астероид проходит через границу атмосферы.

3.2. Построение столкновительной полосы на поверхности Земли. Были построены 11 семейств попадающих траекторий An (n = 0, 1, …, 9), Bn (n = = 1, 2, …, 10), для каждого из которых значение перигейного расстояния меняется в некотором узком диапазоне: A0 = {rπ: 2069–2100 км}, (A1, B1) = = {rπ: 2500 ± 1 км}, (A2, B2) = {rπ: 3000 ± 1 км}, (A3, B3) = {rπ: 3500 ± 1 км}, (A4, B4) = {rπ: 4000 ± 1 км}, (A5, B5) = {rπ: 4500 ± 1 км}, (A6, B6) = {rπ: 5000 ± ± 1 км}, (A7, B7) = {rπ: 5500 ± 1 км}, (A8, B8) = {rπ: 6000 ± 1 км}, (A9, B9) = {rπ: 6360 ± 1 км}, B10 = {rπ: 6375 ± 1 км}, и получены точки соударения астероида с Землей для каждого из всех указанных семейств An, Bn.

На рис. 3 в плоскости географических углов (φ, λ) приведена область точек падения астероида для попадающих траекторий семейства A0 с перигейным расстоянием rπ из диапазона [2069; 2100] км. Центр этой области (ее также обозначим через A0) отмечен крестиком. При этом связна только часть семейства A0 = {rπ: 2069 км ≤ rπ ≤ 2100 км}, примыкающая к орбите T0 c минимальным rπ, для этой части 2069 км ≤ rπ ≤ 2085 км. При rπ > 2085 км множества {rπ = const} несвязны, они состоят из двух подмножеств. Для одного: точки падения лежат правее центра области A0, оно соответствует уменьшению гелиоцентрической энергии h до соударения по сравнению с траекторией T0. Для другого точки падения лежат левее области А0, оно соответствует увеличению энергии по сравнению с траекторией T0.

Рис. 3.

Точки соударения с Землей астероида Апофис в области A0.

В качестве примера на рис. 4 приведены области точек падения астероида для семейств А9 и B9 с rπ из диапазона [6360 ± 1] км на плоскости (δλ cosφ, δφ), где δλ, δφ – долгота и широта точки падения относительно координат центра этой области падения. Здесь же указаны изолинии rπ = = const. Каждая область падения А9 и B9 разделена этими изолиниями на четыре подмножества по перигейному расстоянию rπ, так что на каждом из этих подмножеств rπ меняется на 0.5 км. Например, для крайних темных областей rπ меняется в диапазоне [6360.5; 6361] км.

Рис. 4.

Точки падения астероида для областей A9 (слева) и B9 (справа) (rπ = 6360 ± 1 км).

Все множество столкновительных орбит астероида порождает столкновительную полосу на Северном полушарии Земли. Она приведена на рис. 5. Для лучшего анализа полосы падения астероида полученные результаты сохраняются в файле формата “.kml”, который позволяет отображать трехмерные геопространственные данные и другую полезную информацию в программе “Google Earth”, а также на двухмерных картах Google.

Рис. 5.

Столкновительная полоса астероида Апофис на северном полушарии Земли.

На рис. 6 приведена эта полоса, построенная нами на карте мира с помощью системы “Google my maps” (здесь учтены прецессия и нутация в соответствии с анализом в п. 3.5), ее можно посмотреть на сайте по ссылке (https://drive.google.com/open?id=1VXX1dxJO2e6hn51FFTVf_wbk9gs&usp=sharing). Увеличивая масштаб этого изображения, можно просмотреть все точки данной полосы. Ее начало (Start, семейство А9) – вблизи Омска, Томска, Новосибирска. Далее она последовательно проходит через Красноярский край, Якутию, Дальний восток и Камчатку России, вблизи Алеутских островов, через Тихий океан, Центральную Америку, Южную Америку (Колумбия, Венесуэла), и уходит дальше в Атлантику, до конца (The End, семейство В10) у западного побережья Африки, вблизи Дакара (см. рис. 6).

Рис. 6.

Вероятная столкновительная полоса астероида Апофис на поверхности Земли.

Ниже, в табл. 2 указаны угловые координаты, широты и долготы (с учетом прецессии и нутации), а также местонахождение на Земле центров областей падения An, Bn. Область падения для траекторий семейства A0 с наименьшими значениями rπ находится примерно в центре столкновительной полосы – в Тихом океане, в 1000 км от западного побережья Мексики. Для каждого семейства (An, Bn) при 2085 км < rπ < 6364 км получим две области падения, находящиеся примерно симметрично относительно центра области А0, к западу (An) и востоку (Bn) от него, см. рис. 5 и 6. Отметим, что ранее эту столкновительную полосу привели R. Schweickart, et al. [9], D.B. Gennery [10]. Сравнение данных [9, 10] с нашими результатами показывает, что они практически полностью совпадают, что говорит, в частности, о совпадении отечественных и зарубежных результатов в определении орбиты Апофиса и в его модели движения.

На рис. 7 приведены некоторые энергетические, временные и геометрические характеристики столкновительных траекторий для центров областей падения An, Bn. Области Аn порождаются орбитами астероида с большей энергией перед соударением, чем в области А0, h ∈ (−799.7689, −799.7554] км22, здесь большая полуось a ∈ (а0, а0 + 2784] км, а0 ≈ 165 938 501 км, см. рис. 7, кривая h(λ). Области Bn порождаются орбитами с меньшей энергией, h ∈ [−799.7824, −799.7689) км22, a ∈ [а0−2817, а0) км. Из-за несферичности Земли и наклонности траекторий и полосы к экватору, расстояния от центра Земли до центров областей А9 и В9 различны и составляют ~6364 и ~6376 км, т.е. траектории астероида из семейства А9 почти касаются поверхности Земли, а траектории из семейства В9 еще не касаются Земли, так как в этих областях rπ = 6360 ± 1 км. Поэтому построено еще семейство В10 с rπ из диапазона 6375 ± 1 км, создающее восточную область В10.

Рис. 7.

Характеристики столкновения астероида с Землей в 2036 г. в зависимости от долготы точки соударения.

Для кривой времени tctc0 (λ) на рис. 7 опорное время tc0 = 9 ч, 13.IV.2036 г., это, примерно, среднее время падения. Для всего множества траекторий время столкновения tc находится в диапазоне от tc0 – 17 мин до tc0 + 13.8 мин, см. рис. 7, кривая tctc0 (λ). Абсолютная скорость столкновения астероида с Землей составляет Vca ≈ 12.6 км/c, аналогично [2, 47, 10]. Учет вращения Земли приводит к увеличению скорости соударения в западной части полосы до ~12.9 км/c и к ее уменьшению на восточной части до ~12.2 км/c, рис. 7, кривая vcr (λ). На рис. 3, 6 для областей A9, A0 и B10 показаны скорости vηζ, дающие направления подлета астероида к Земле.

На рис. 8 приведены угол входа астероида в атмосферу θr, Δθ = θa θr и перигейное расстояние rπ в зависимости от долготы центров областей падения An, Bn; где Δθ разность углов входа в атмосферу для абсолютной и относительной скоростей. Наибольший угол входа, около 60° – в области А0. С увеличением rπ угол входа монотонно уменьшается, до ~8° для областей A9 и B10.

Рис. 8.

Угол входа астероида в атмосферу θr, Δθ и перигейное расстояние rπ орбиты астероида в зависимости от долготы точки соударения.

3.3. Геометрия столкновительной полосы. Рассмотрим геометрию столкновительной полосы в гринвичской вращающейся СК. Радиус-векторы центров трех областей падения Аi, А0 и Вi (i = 1, …, 9) задают в этой СК плоскости Пi с векторами нормали к ним ni. Углы между этими векторами nj (j = 1, …, 8) и n9 = [0.1909, –0.3842, 0.9033]Т малы, до 2.7°. Поэтому можно считать, что полоса столкновений лежит примерно в одной плоскости П9. Угол наклона вектора n9 к плоскости экватора Земли составляет 64.6°, и ~87.6° – к плоскости эклиптики (в момент tc0), т.е. плоскость полосы примерно параллельна плоскости эклиптики в этот момент.

Выполнен анализ по оценке длины и ширины полосы на поверхности Земли. Считаем приближенно, что пересечение плоскости полосы и Земного шара является окружностью, ее радиус RП = = 5370 км. Тогда центральный угол полосы ϕ ≈ ≈ 4.534 рад = 259.8°, а ее длина LП ≈ ϕ RП = 24 347 км. Для оценки ширины столкновительной полосы построены системы координат oixiyizi, начало каждой системы находится в центре ее области падения. При этом ось oixi направлена вдоль продольной оси полосы на восток, и ось oizi направлена по нормали к местному горизонту. Ось oiyi, направленная “поперек” полосы, дополняет СК до правой. Плоскость oixiyi соответствует местному горизонту. Преобразуя координаты точек падения в СК oixiyi, и определяя диапазон координаты yi для точек области, получим оценку ширины полосы в этой области. Ширина области А0 максимальна, она оценена в ~27 км. Для областей Ai и Bi (i ≥ 1) ее оценка дала 10–22 км. В среднем, ширина полосы составляет ~20 км. Отметим, что ранее ширина полосы падения Апофиса оценивалась в 30–50 км [14, 15], т.е. близко. Выявлено, что ширина полосы зависит от принятых разбросов начального вектора состояния (1). Для анализа этого фактора рассмотрены СКО начального состояния, в два раза большие и в два раза меньшие прежних значений (2), для множества траектории с rπ = 6000 км, A8 и B8 (см. табл. 1). Анализ показал, что ширина полосы примерно прямо пропорциональна СКО начального состояния.

Таблица 1.  

Ширина зоны падения для семейств траекторий А8, B8 с rπ = 6000 ± 1 км в зависимости от СКО начального состояния

СКО начального состояния Ширина зоны падения, км
А8 В8
σ = σ0/2 7.21 7.51
σ = σ0 14.32 14.96
σ = 2σ0 28.22 29.52
Таблица 2.  

Географические координаты точки падения и соответствующее название места падения

Область падения ${\tilde {\varphi }}$с, град $\tilde {\lambda }$с, град Название места падения
А9 54.3 75.2 Омская обл, около г. Черлак
А9–А8     Около Томска, р. Ангара
А8 60.0 114.6 Около г. Хамра (Якутия, Сибирь)
А7 58.5 140.7 г. Курун-Урях (Хабаровский край, Дальний Восток)
А6 55.3 158.3 Камчатский край РФ, ~300 км от Петропавловска-Камчатского
А5 51.6 171.7 Северная часть Тихого океана, ~ к югу от Алеутских островов
А4 47.6 –177.2
А3 43.3 –167.1 Северная часть Тихого океана, к югу от Аляски,
у Западного побережья США
А2 38.7 –157.3
А1 33.2 –146.1
A0 23.7 –126.0 Тихий океан, 1000 км от западного побережья Мексики
B1 16.3 –106.8 Тихий океан между Гавайскими островами и Центральной Америкой
B2 13.5 –96.8
B3 11.8 –88.4
B3–B4     Пересекает Центральную Америку по Никарагуа, вблизи Манагуа
B4 10.7 –80.6 Карибское море
B5 10.1 –72.7 г. Мачике (Венесуэла)
B6 10.0 –64.4 г. Барселона (Венесуэла)
B7 10.4 –54.8 Атлантический океан
B8 11.7 –42.3
B9 14.5 –23.4 Около г. Прая (Кабо-Верде), острова Зеленого Мыса
B10 15.1 –19.4 Атлантический океан, вблизи от Дакара (Сенегал, Западная Африка)

Определена также ориентация продольной оси полосы вдоль всей зоны падения. Угол χ наклона продольной оси полосы относительно восточной параллели меняется в диапазоне от +26° для семейства А9 (+ соответствует повороту вокруг оси + oizi против часовой стрелки) до –32° для семейств А2, А3, т.е. продольная ось полосы близка к параллели, см. рис. 9.

Рис. 9.

Угол ориентации продольной оси полосы относительно параллели в зависимости от географической долготы λ.

3.4. Качественный анализ подлета к Земле столкновительных траекторий Апофиса. Чтобы лучше понять качественные особенности множества орбит соударения астероида с Землей, рассмотрим векторы прицельной дальности b и векторы скорости “на бесконечности” v гиперболических орбит Апофиса перед столкновением с Землей [24, 25]. Численный анализ показал, что для всех столкновительных орбит Апофиса концы векторов b, отложенных из одной точки, лежат практически на одной прямой L, единичный вектор направления которой L0 = [–0.2443, 0.8776, 0.4124]T. Векторы v почти постоянны и образуют практически плоскость, проходящую через L, она переходит затем в плоскость полосы. Среднее значение векторов v составляет вектор v = = [5.5849, 0.7923, 1.6077]T км/c, |v| = 5.87 км/c. На рис. 10 показано геоцентрическое движение астероида Апофис для всего множества траекторий столкновения с Землей в 2036 г. Здесь приведены векторы прицельной дальности и векторы скорости “на бесконечности”, а также столкновительная полоса на поверхности Земли. Изолинии rπ = = const будут примерно перпендикулярны соответствующему вектору b. Так, для области А0, при rπrπmin ≈ 2069 км изолиниями rπ = const будут дуги, идущие вдоль оси полосы (см. на рис. 3), примерно перпендикулярно среднему для области A0 вектору прицельной дальности b0. С увеличением перигейного расстояния rπi векторы bi поворачиваются, поэтому наклон изолиний rπi = const к продольной оси столкновительной полосы монотонно возрастает. На рис. 3 виден небольшой начальный поворот этих изолиний при смещении от центра А0. На краях полосы, в областях А9, В9, В10, изолинии, будучи перпендикулярными к вектору bi, становятся почти перпендикулярными к оси полосы, что видно на рис. 4. Для других областей эти повороты изолиний можно проследить на рис. 6 при увеличении его масштаба с помощью системы “Google my maps” по ссылке (https:// drive.google.com/open?id=1VXX1dxJO2e6hn51FFTVf_ wbk9gs&usp=sharing).

Рис. 10.

Вероятное множество геоцентрических траекторий астероида Апофис при падении на Землю в 2036 г.

3.5. Определение точки падения с учетом прецессии и нутации. Для уточненного перехода от геоэкваториальной инерциальной системы координат на эпоху J2000 (вектор r) к гринвичской системе координат текущего момента (вектор rG) применяется модель преобразования координат точки с учетом прецессии, нутации и движения полюса Земли [1923]:

(12)
${{{\mathbf{r}}}_{{\text{G}}}} = {\text{WRNP}} \cdot {\mathbf{r}},$
где W − матрица движения полюса; R − матрица суточного вращения Земли; N − матрица нутации; P − матрица прецессии. При прецессии земная ось описывает конус вокруг Северного полюса эклиптики со скоростью ~50.37″ в год. Нутация представляет собой колебания с амплитудой ~9.21″, накладываемые на прецессионное движение.

Матрицы прецессии и нутации определяются в соответствии с моделью прецессии-нутации МАС2000 [19, 20, 23]. Отметим, что по теории нутации МАС2000 разложения нутации в долготе и в наклоне экватора включают 1365 членов (678 лунно-солнечных и 687 планетных), здесь мы учитывали лишь 10 главных членов для вычисления лунно-солнечной нутации [23]. Так как изменения координат полюса достаточно малы, то влияние движения земных полюсов на преобразование координат точки падения в нашей задаче не учитывалось. При этом преобразование (12) учитывает вращение Земли, прецессию и нутацию. Для проверки матрица перехода (12) определена также по сайту (http://hpiers.obspm.fr/eop-pc/ index.php?index=matrice_php&lang=en). Сопоставление результатов показывает, что они хорошо совпадает, точность по элементам матрицы составляет ~10–5.

На основании этой методики выполнен численный анализ определения координат точек падения Апофиса в рамках данной модели. В табл. 2 приведены географические широта и долгота ${{{\tilde {\varphi }}}_{с}},$ ${{\tilde {\lambda }}_{с}}$ точки падения Апофиса в центрах областей падения астероида Апофис на Землю с учетом прецессии и нутации. Изменение величины географической широты и долготы (по сравнению со случаем без учета прецессии и нутации) составило Δφ ≈ –0.22°–0.15°, Δλ ≈ 0.4°–0.8°, изменение координат точек составляет примерно 50 км. Полоса соударения астероида Апофис с Землей в 2036 г., помещенная на сайте [26], также сделана с учетом прецессии и нутации.

ВЫВОДЫ

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

Авторы признательны д.ф.-м.н. В.В. Чазову, к. ф.-м. н. В.А. Степаньянцу, сотруднику НПО “Комета” Ю.П. Кулешову за обсуждение настоящей работы и ряд полезных советов.

Список литературы

  1. Ягудина Э.И., Шор В.А. Орбита АСЗ (99942) Apophis = 2004 MN4 из анализа оптических и радарных наблюдений // Всероссийская конференция “Астероидно-кометная опасность-2005 (АКО-2005)”. Санкт-Петербург. 3–7 октября 2005 г. Материалы конференции. СПб: ИПА РАН, 2005. С. 355–358.

  2. Ivashkin V.V., Stikhno C.A. An Analysis of the Correction Problem for the Near-Earth Asteroid (99942) Apophis = 2004 MN4. // 2007 Planetary Defense Conference. G. Washington University, Washington, D.C., USA. March 5–8, 2007.

  3. Ивашкин В.В., Стихно К.А. Анализ проблемы коррекции орбиты астероида Apophis. В кн. “Околоземная астрономия 2007” // Труды Международной конференции. 3–7 сент. 2007. Терскол. КБНЦ РАН, Нальчик. 2008. С. 44–48.

  4. Ивашкин В.В., Стихно К.А. О проблеме коррекции орбиты сближающегося с Землей астероида (99942) Apophis // Доклады Академии наук. 2008. Т. 419. № 5. С. 624–627.

  5. Ивашкин В.В., Стихно К.А. О предотвращении возможного столкновения астероида Apophis с Землей // Астрономический вестник. 2009. Т. 43. № 6. С. 502–516.

  6. Ивашкин В.В., Стихно К.А., Гуо П. О структуре множества вероятных траекторий соударения астероида Апофис с Землей в 2036 г. // ДАН. 2017. Т. 475. № 4. С. 389–394.

  7. Соколов Л.Л., Башаков А.А., Борисова Т.П. и др. Траектории соударения астероида Апофис с Землей в XXI в. // Астрон. вестн. 2012. Т. 46. № 4. С. 311–320.

  8. Соколов Л.Л., Кутеева Г.А. Возможные соударения астероида Апофис после уточнения его орбиты // Вестник СПБГУ. 2015. Сер. 1. Т. 2(60). Вып. 1. С. 148–156.

  9. Russell Schweickart et al. Threat Characterization: Trajectory Dynamics (White Paper 39). B612 Foundation. 2006 // URL: https://www.researchgate.net/ publication/2176734_Threat_Characterization_Trajectory_Dynamics

  10. Gennery D.B. Scenarios for dealing with Apophis // 2007 Planetary Defense Conference. Washington, DC. 2007.

  11. Giorgini J.D., Benner L.A.M., Ostro S.J. et al. Predicting the Earth encounters of (99942) Apophis // Icarus. 2008. V. 193. № 1. P. 1–19.

  12. Chesley S.R. Potential Impact Detection for Near-Earth Asteroids: The Case of 99942 Apophis (2004 MN4) // Proc. 229th IAU Symp., 2005 / Eds. Ferraz Mello, S., Lazaro, D., and Fernandes, J. New York, Melbourne, Madrid: Cambridge Univ. Press, 2006. P. 215–228.

  13. Wlodarczyk I. The potentially dangerous asteroid (99 942) Apophis // Monthly Notices of the Royal Astronomical Society. 2013. V. 434. P. 3055–3060.

  14. GALSPACE Малые тела Солнечной системы. Столкновение Земли с астероидом – проблема Апофиса http://galspace.spb.ru/index129.html

  15. Wikipedia / B612 Foundation https://en.wikipedia.org/ wiki/B612_Foundation

  16. Ferrier L., Vérant J.-L., Moschetta J.-M. AeroThermodynamical Study for the Entry of an Apophis-Like Asteroid // 49th AIAA Aerospace Sciences Meeting including the New Horizons Forum and Aerospace Exposition, Aerospace Sciences Meetings. AIAA 2011-1037.

  17. Степаньянц В.А., Львов Д.В. Эффективный алгоритм решения системы дифференциальных уравнений движения. // Математ. модел. 2000. Т. 12. № 6. С. 9–14.

  18. Математика в техническом университете / Под ред. Зарубина В.С. Выпуск XVI. Теория вероятностей. М.: МГТУ им. Н.Э. Баумана, 1999.

  19. Kaplan G.H. The IAU Resolutions on Astonomical Reference Systems, Time Scales, and Earth Rotation Models. Explanation and Implementation. U.S. Naval Observatory. Circular № 179. Washington, D.C. 20392. 2005.

  20. Аксенов Е.П., Чазов В.В. Модель движения ИСЗ. Главная проблема. Основные алгоритмы. М.: ГАИШ МГУ, 2011.

  21. Аджян А.П., Аким Э.Л., Алифанов О.М. и др. Ракетно-космическая техника. Машиностроение. Энциклопедия. Т. IV-22 / Под ред. Легостаева В.П. М.: Машиностроение, 2012. С. 15–26.

  22. Жаров В.Е. Сферическая астрономия. М.: Фрязино, 2006.

  23. Крылов В.И. Координатно-временные преобразования в геодезии: учебное пособие. М.: МИИГАиК, 2014.

  24. Аким Э.Л., Энеев Т.М. Определение параметров движения космического летательного аппарата по данным траекторных измерений // Космич. исслед. 1963. Т. 1. Вып. 1. С. 5–50.

  25. Абалакин В.К., Аксенов Е.П., Гребенников Е.А. и др. Справочное руководство по небесной механике и астродинамике. М.: Наука, 1976.

  26. Сайт рисунка столкновительной полосы Апофиса на поверхности Земли в 2036 г.: https://drive.google. com/open?id=1VXX1dxJO2e6hn51FFTVf_wbk9gs& usp=sharing

Дополнительные материалы отсутствуют.