Астрономический журнал, 2022, T. 99, № 3, стр. 250-264

3.1. Методика расчета движения метеороидных частиц

После того как частицы достаточно отдалились от ядра родительского тела (на ${{R}_{{esc}}}$), можно рассматривать дальнейшее их движение как движение в общем гравитационном поле Солнечной системы с учетом сил радиационного характера. Для этого используется численный подход. Модель гравитационного поля Солнечной системы включает гравитацию от Солнца и планет. Начальные условия для гравитирующих тел были взяты из базы даных NASA JPL HORIZONS²². Под силами радиационного характера подразумеваются сила давления излучения и эффект Пойтинга-Робертсона, общее действие которых можно записать следующим образом [20]:

Расчет движения осуществлялся при помощи программного комплекса REBOUND [21], предоставляющего набор инструментов для моделирования широкого класса задач динамики гравитирующих систем. Для интегрирования была выбрана гибридная численная схема MERCURIUS [22], являющаяся развитием широко известного кода MERCURY [23]. В данном модуле интегрирование ведется при помощи двух схем: основная – симплектическая WHFast с постоянным шагом [24], вспомогательная – IAS15 с переменным шагом [25], которая задействуется в случае возникновения сближений частиц с планетами. Такой подход позволяет сочетать симплектическую схему и достаточно точную обработку близких сближений. Для основной схемы был выбран шаг интегрирования по времени 12 ч. Отметим, что гибридное интегрирование часто применяется в задачах моделирования динамики метеороидных роев (см. [7, 8]).

3.2. Схема определения значений индекса массы $s$

Опишем процедуру определения значений $s$ в потоке. Поскольку мы сравниваем модель с наблюдениями метеорных потоков, логично определять значения индекса $s$ в потоке на расстоянии 1 а.е. от Солнца. На рис. 1 проиллюстрирована схема подсчета частиц в метеороидном потоке. Земля пересекает метеороидный поток. В околоземной сфере радиусом ${{R}_{{cs}}}$ (сфере подсчета) в моменты времени близ ожидаемого максимума потока подсчитывается количество частиц каждого размера, соответствующего группе $j$ (см. табл. 1). Последовательные положения сферы соответствуют изменению долготы Солнца со скважностью 1°.

Рис. 1.

Картина пересечения Землей метеороидного потока, cхематично показаны сферы подсчета частиц и орбита родительского тела (кометы) (a) и сечение метеороидного потока и проекция орбиты Земли на это сечение (b). Точками различных размеров указаны положения траекторий частиц потока. Красная точка соответствует орбите кометы.

Естественно, что Земля не обязательно пересекает самый центр потока, т.е. наиболее плотную его часть. Структура потока (здесь – распределение плотности потока в сечении) скорее всего сложнее, чем это схематично показано на рис. 1b. Что касается размеров сферы подсчета, то, конечно, было бы корректнее брать значение ${{R}_{{cs}}}$ поменьше, но это благое пожелание наталкивается на необходимость значительного (на многие порядки) увеличения числа модельных частиц. По нашим оценкам свертка по сфере радиусом порядка 0.1–0.2 толщины потока является приемлемым приближением. Толщина потока определяется по его продолжительности. Средняя оценка продолжительности сильных (ZHR > 10) метеорных дождей – от нескольких суток (например, Дракониды) до нескольких месяцев (Южные Тауриды, ассоциируемые с кометой 2P/Encke). В работе [26] характерные длительности периода наблюдений для 6 сильных метеорных потоков составляли от 7 до 30 сут, а поскольку индекс $s$ определялся каждые сутки, это позволило довольно детально построить ход изменения $s$ за время пересечения потоков.

Мы определяли значения индекса массы, учитывая весь набор размеров (масс) частиц. Минимальное число размеров (масс) для определения $s$, естественно, равно 2 (см., напр., [8], выражение (1)). В нашей модели индекс массы определяется по двум-четырем значениям плотности частиц различных групп. Понятно, что крупных частиц изначально было намного меньше, чем мелких, поэтому, чтобы как-то учесть это обстоятельство, влияющее на статистическую достоверность, мы ввели веса для разных групп частиц в процедуре определения $s$. Вес для каждой группы был выбран пропорциональным десятичному логарифму (нормированного) количества частиц.

4.1. Общая эволюция метеороидного роя

На рис. 2 проиллюстрирована эволюция метеороидного роя, выброшенного из родительского тела в модели М1 (базовой) на протяжении ~1000 лет. Показано положение частиц (метеороидов) различного размера в проекциях на плоскость эклиптики $XY$ и на момент 501.6 лет в плоскостях $XZ$ и $YZ$. Здесь и далее цветовая кодировка соответствует группам (по размеру) частиц в табл. 1.

Рис. 2.

Эволюция метеороидного потока в базовой модели М1 (комета 96P/Machholz). Показано положение частиц (метеороидов) различного размера в проекции на плоскость эклиптики $XY$ и (на момент 501.6 лет) в плоскостях $XZ$ и $YZ$. Цветовое обозначение соответствует табл. 1. Черными точками и линиями показаны планеты и их орбиты, пурпурной точкой – положение родительского тела на указанный момент времени.

В самом начале эволюции, когда первый выброс частиц уже произошел, но родительское тело прошло по орбите после выброса лишь небольшую часть периода обращения (0.2 года при периоде 5.28 лет), формируется удлиненный след из частиц, причем хорошо видно разделение траекторий движения частиц различных размеров. Естественно, что наиболее сильное влияние негравитационные факторы оказывают на частицы меньшего размера. Самые массивные “синие” частицы отходят от ядра кометы симметрично на максимальное расстояние, определяемое скоростью начального выброса. Более легкие частицы уходят от ядра гораздо дальше. Видно, например, что частицы наименьшего размера ($r = 0.005$ см) уже “забыли” про скорость начального выброса и ни о какой симметрии нет речи. Облако частиц удалилось от ядра более чем на 4 млн. км, а через 1 год на еще более далекое расстояние, порядка астрономической единицы. В целом же ансамбль метеороидов образует вытянутую структуру, которая со временем становится похожей на узкий хвост, вытянутый вдоль орбиты кометы в обоих направлениях. Отстающая часть хвоста длиннее передней. С 1983 г. такие структуры, называемые dust trails, которые образуются свежими выбросами из кометных ядер, стали объектом наблюдений. С помощью инфракрасного телескопа на спутнике IRAS были получены изображения таких длинных и узких структур в ряде активных комет (см. очень детальное описание в [16]).

Отметим, что поскольку у нас рассматриваются частицы не с непрерывным распределением по размерам, а четырех групп размеров, на рис. 2 видны отдельные нити орбит. На самом деле, распределение частиц по размерам непрерывное и общий вид потока в реальности будет более размытым (непрерывным). Количество самых массивных модельных частиц, показанных на рис. 2, относительно завышено (см. обсуждение табл. 1). В реальности таких частиц в сравнении с маломассивными намного меньше и картинка должна быть гораздо более “желтой”.

На рис. 2 на более поздних стадиях можно заметить различные структурные изменения в потоке, вызванные гравитационными факторами – сближениями с планетами и прежде всего с Юпитером. Тяжелые частицы (синий и зеленый цвет) мало подвержены действию негравитационных факторов, но сближение с Юпитером вызвало появление дополнительных ветвей в потоке этих частиц. Мы не исследуем этот вопрос детально, т.к., напомним, главной целью работы является исследование поведения спектра масс метеороидов.

Естественный и весьма важный вопрос – как оценивать вековые изменения потока. В данной работе применяется метод, позволяющий отследить эволюцию распределения орбитальных элементов всех частиц в ансамбле. Он состоит в построении диаграмм, представляющих изменение элементов орбит членов ансамбля метеороидов каждой из групп со временем. На рис. 3 показаны изменения массива значений большой полуоси $a$, эксцентриситета $e$ и наклонения $i$ орбит частиц в модели М1. Видно, что у группы массивных частиц распределение $a$ почти не меняется и в каждый момент времени концентрируется к текущему значению $a$ орбиты родительского тела. А вот у групп более легких частиц, для которых вклад негравитационных факторов более существенен, изменение значений массива $a$ происходит тем заметнее, чем меньше частицы. Самые легкие частицы вообще очень быстро “забывают” исходную орбиту родительского тела и рассеиваются в фазовом пространстве. Весьма показательно изменение наклонений орбит частиц. Более массивные частицы, естественно, в большей степени “отслеживают” эволюцию изменения наклонения орбиты ядра кометы. Самые мелкие частицы рассеиваются гораздо сильнее. Можно представить, что из самой мелкой фракции образуется дискообразная структура. Здесь с нашей точки зрения вполне уместно обсуждать возможный механизм образования дискообразного пылевого околосолнечного облака, ответственного за явление зодиакального света.

Рис. 3.

Эволюция ансамбля частиц в модели М1 (комета 96P/Machholz). Показана эволюция распределений большой полуоси $a$, эксцентриситета $e$ и наклонения $i$ для частиц разных групп. Цветовая шкала справа соответствует плотности распределения частиц, т.е. количеству частиц в интервале $\Delta a = 0.1$ а.е., либо $\Delta e = 0.005$, либо $\Delta i = 0.7^\circ $. Линия пурпурного цвета показывает эволюцию элементов орбиты родительского тела.

На вcех диаграммах рис. 3 видно, что существенные изменения в структуре потока происходят в период 500–700 лет после начала эволюции потока. По нашему мнению, это связано с тем, что при $t \sim 500$ лет наклонение орбиты ядра кометы и в целом потока становится минимальным, что существенно изменяет характер сближения с планетами (прежде всего с Юпитером), возмущения становятся гораздо более выраженными. В потоке выделяются отдельные ветви (см. обсуждение изменений структуры потока, показанных на рис. 2). Со временем (где-то после 750 лет) наклонение орбиты ядра кометы и потока в целом (за исключением самых мелких частиц) снова растет, сближения играют менее существенную роль и структура потока в фазовой плоскости $a{\kern 1pt} - {\kern 1pt} i$ выходит на некоторый “стационар”. Для более детального анализа этой сложной картины требуются дополнительные исследования.

На рис. 4 показана эволюция структуры модельного потока, образованного при испарении ядра кометы, орбита которой совпадает с современной орбитой кометы Энке. Картина эволюции сходна с показанной на рис. 2 для кометы 96P/Machholz. Конечно, эти данные требуют более углубленного анализа. Отметим лишь уже упомянутую выше деталь, которая здесь проявляется гораздо более отчетливо – образование со временем дискообразной структуры из самых мелких метеороидов. Это еще более привлекает внимание к проблеме образования и поддержания околосолнечного пылевого диска, ответственного за явление зодиакального света. Но снова отметим, что эта интересная тема не является главной в данной работе.

Рис. 4.

Эволюция метеороидного потока в базовой модели М2 (комета 2P/Encke). Показано положение частиц (метеороидов) различного размера в проекции на плоскость эклиптики $XY$ и плоскости $XZ$ и $YZ$. Цветовое обозначение соответствует табл. 1. Черными точками и линиями показаны планеты и их орбиты, пурпурной точкой – положение родительского тела на указанный момент времени.

На рис. 5 показаны изменения массива значений большой полуоси $a$, эксцентриситета $e$ и наклонения $i$ орбит частиц в модели М2 (комета 2P/Encke). Видно, что как и в модели M1, у группы массивных частиц распределение $a$ почти не меняется, и в каждый момент времени концентрируется к значению $a$ орбиты родительского тела. Тем не менее разброс по $a$ более выражен по сравнению с моделью М1. Также хорошо видны различные периодические детали на этих диаграммах. По нашему мнению, это может объясняться тем, что орбиты кометы Энке и ассоциированных с ней модельных метероидов лежат в среднем гораздо ближе к плоскости эклиптики, сближения с планетами, расположенными ближе чем Юпитер, гораздо более частые, и возмущения от таких сближений более существенны, чем в модели М1.

Рис. 5.

Эволюция ансамбля частиц в модели М2 (комета 2P/Encke). Показана эволюция распределений большой полуоси $a$, эксцентриситета $e$ и наклонения $i$ для частиц разных групп. Цветовая шкала справа отличается от той, что показана на рис. 3. Она соответствует плотности распределения частиц, т.е. числу частиц в интервалах $\Delta a = 0.008$ а.е., либо $\Delta e = 0.0025$, либо $\Delta i = 0.33^\circ $. Линия пурпурного цвета показывает эволюцию элементов орбиты родительского тела.

4.2. О зависимости результатов от параметров модели

Как отмечалось в [1], расчетные начальные скорости выброса частиц могут заметно отличаться в разных моделях. Также цитируемые в [1] авторы рассматривают разную геометрию выброса. Для исследования зависимости результатов от параметров выброса были проведены вспомогательные сравнительные расчеты, в которых условия выброса параметризовались. Основные параметры – это геометрия и скорость выброса. В базовом расчете (модель М1) моделировался выброс в сферу с терминальными скоростями, указанными в табл. 1, в модели М1.1 выброс происходил с полусферу, в модели М1.2 выброс происходил в узкий конус (угол 10°) в направлении на Солнце, а в модели М1.3 была в два раза увеличена скорость выбрасываемых частиц по сравнению с моделью M1.

На рис. 6, 7 и 8 показаны распределения частиц (всех размеров) по $a$, $e$, $i$ соответственно после 10 периодов обращения родительского тела вокруг Солнца при различных допущениях о характере выброса частиц из ядра кометы. Из рис. 6 видно, что структура распределения числа частиц по большой полуоси для выбросов в полусферу и в сферу практически не отличается для частиц всех размеров. При выбросе в узкий конус характер распределения тот же, для частиц группы 1 (т.е. наименьшего размера), но для более крупных частиц распределение по $a$ ýже. Выброс с двойной скоростью приводит к существенному изменению геометрии потока. Наиболее сильно это проявляется для самых мелких частиц, но и для более крупных частиц распределение по $a$ в разы шире, чем в случае выбросов с неувеличенной скоростью. Из анализа рис. 6 можно сделать заключение, что результаты слабо зависят от геометрии выброса, если не принимать крайние предположения (очень узкий конус выброса), а вот зависимость от скорости выброса намного существеннее. Из рис. 7 и 8 видно, что структура распределения числа частиц по эксцентриситету и наклонению для всех видов выбросов внутри данной группы отличается несущественно.

Рис. 6.

Распределение частиц (всех групп размеров) по $a$, после 10 периодов обращения вокруг Солнца при различных допущениях о характере выброса частиц из ядра кометы. Выброс в полусферу (оранжевая линия), выброс в сферу (синяя линия), выброс в узком конусе (красная линия) и выброс в сферу, но с удвоенной (по сравнению с табл. 1) скоростью (зеленая линия).

Рис. 7.

Распределение частиц (всех групп размеров) по $e$ после 10 периодов обращения вокруг Солнца при различных допущениях о характере выброса частиц из ядра кометы. Обозначения те же, что и на рис. 6.

Рис. 8.

Распределение частиц (всех групп размеров) по $i$ после 10 периодов обращения вокруг Солнца при различных допущениях о характере выброса частиц из ядра кометы. Обозначения те же, что и на рис. 6.

4.3. Изменение индекса массы $s$ в потоке

Теперь рассмотрим распределение метеороидов на расстоянии 1 а.е. в плоскости эклиптики, т.е. в той области, где движется Земля, а при входе Земли в наиболее плотную часть потока метеороидов может наблюдаться метеорный дождь. На рис. 9 на левой панели показаны распределения частиц потока в сферах подсчета (см. раздел 3) при ${{R}_{{cs}}} = 0.1$ а.е. в двадцатидневный период с центральным значением $t = 501.6$ года для модели М1 (комета 96P/Machholz). Редуцирующий фактор для учета повторных выбросов был выбран равным $f = 0.97$. Подсчеты проведены в моменты, соответствующие изменениям долготы Солнца на 1 градус (т.е. практически каждые сутки). Количество тяжелых частиц (группы 4) показано как есть, а согласно изложенному в разделе 2 приему, количество расчетных частиц из группы $j$ в сфере подсчета умножено фактор Norm$(j)$ (см. табл. 1). Напомним, что это число модельных, а не реальных частиц, и наши модели дают представления в основном об относительных характеристиках потоков. Величина $s$ как раз является одной из таких относительных величин. На правой панели рис. 9 крестиками показан соответствующий данным на левой панели профиль индекса массы $s$. Для сравнения с наблюдениями профиль нанесен на наблюдательную диаграмму (см. [26, рис. 7], отражающую результаты наблюдений потока Ариетиды в 2008 г. с радаром CMOR). Поведение $s$ при прохождении Земли сквозь поток качественно именно таково, как описано в [1]. В максимуме потока $s < 2$, на краю потока $s \sim 2$ или выше. То есть из рис. 9 видно, что есть неплохое качественное согласие между изменением во времени индекса $s$ в модели М1 и наблюдениями.

Рис. 9.

Слева – распределения числа частиц потока в сферах подсчета при ${{R}_{{sc}}} = 0.1$ а.е. в двадцатидневный период с центральным значением $t = 501.6$ года для модели М1 (комета 96P/Machholz), цвета соответствуют обозначениям в табл. 1. Справа – модельные значения индекса массы $s$, соответствующие левой панели (крестики), для сравнения приведены наблюдения (точки), построенные по данным о потоках Ариетид в 2008 г. (см. [26, рис. 7]).

На рис. 10 показано распределение индекса массы в ассоциированном с кометой 2P/Encke модельном потоке (модель М2). Сферы подсчета распределены со скважностью 1 градус по окружности радиусом 1 а.е. от Солнца в плоскости эклиптики по всему диапазону долгот. Это распределение $s$ получено на момент $t = 511.6$ лет от начала расчетов. Фактор редукции $f = 0.97$. Видно, что метеороидный поток от кометы 2P/Encke пересекает орбиту Земли дважды (на долготах ~105° и 230°). Видно, что в периоды пересечения потока с окружностью (почти совпадающей с орбитой Земли) индекс массы ведет себя так же, как в случае кометы 96P/Machholz (см. рис. 9).

Рис. 10.

Индекс массы $s$, вычисленный для модели М2 (комета 2P/Encke) в области орбиты Земли в момент $t = 511.5$ г.

4.4. Влияние времени начала интегрирования

До сих пор мы исследовали поведение модельных потоков, не делая детальной привязки к истории потока. То есть мы не учитывали первый этап в списке этапов, приведенных во Введении. Для детального исследования эволюции конкретного потока, конечно, необходимо знать не только как, но и когда образовался (начал образовываться) поток. Как уже отмечалось, вопрос о возрастах потоков весьма сложен и во многих случаях данные о возрастах противоречивы. В модели М1.4 мы прошли по классической схеме этапов (1–7), рассматривая момент времени ${{t}_{0}}$ как параметр. На рис. 11 показаны изменения индекса массы $s$ в максимуме потока, вычисленные для модели М1.4 в предположении, что истечение из кометы (комета 96P/Machholz) началось в момент ${{t}_{0}} = - 500$ или ${{t}_{0}} = - 1000$ лет. Распределение s показано в момент $t = {{t}_{0}} + 501.6$ лет. Повторные выбросы учитывались с фактором $f = 0.97$. Сравнивая эти диаграммы, а также правую диаграмму на рис. 9 (там было ${{t}_{0}} = 0$, но возраст потока тот же $t = 501.6$ лет), можно сделать вывод, что структура модельного потока одного и того же возраста значимо зависит от начального момента, когда началось истечение (выброс метеороидов) из ядра кометы. Это вполне понятно, ибо как видно из рис. 3 и 5, на временах во многие сотни и тысячи лет орбита ядра кометы меняется существенно, также существенно меняется и картина сближений с планетами, так что для различных значений ${{t}_{0}}$ эволюция потока идет по другому из‑за различной истории сближений с планетами. Однако общий характер поведения индекса масс все же сходен: в максимуме потока $s < 2$, на краю потока $s$ может превышать 2.

Рис. 11.

Слева – изменения индекса массы $s$ в максимуме потока, вычисленные для модели М1.4 в предположении, что истечение из кометы (96P/Machholz) началось в момент ${{t}_{0}} = - 500$ лет, справа – то же, но для ${{t}_{0}} = - 1000$ лет. Распределение $s$ показано в момент ${{t}_{0}} + 501.6$ лет.