Астрономический вестник, 2022, T. 56, № 2, стр. 109-131

ВВЕДЕНИЕ

Космические миссии, начиная с проектов ВЕГА-1, -2 (Sagdeev и др., 1985; 1986; 1987) и Giotto (Keller и др., 1986; 1987) по изучению кометы Галлея и заканчивая недавней миссией Rosetta по изучению кометы 67P/Чурюмова–Герасименко (El-Maarry и др., 2019; Choukroun и др., 2020; Vincent и др., 2019; Groussin и др., 2019; Weissman и др., 2020; Keller, Kührt, 2020), открыли новую эру в изучении физики ядер комет. Уже для кометы Галлея были получены доказательства того, что поверхность ядра покрыта пылевой мантией и ядро обладает значительной пористостью (Sagdeev и др., 1987; Newburn и др., 1991). В ходе выполнения миссии Rosetta различные научные инструменты на борту орбитального и посадочного модулей космического аппарата провели многосторонние исследования, нацеленные на изучение физических свойств кометного ядра, его приповерхностного слоя и кометной пыли. Огромное количество таких наблюдательных данных было получено в первую очередь с инструментами MIRO (Schloerb и др., 2015; Choukroun и др., 2015), COSIMA (Langevin и др., 2016; Hilchenbach и др., 2016), MIDAS (Bentley и др., 2016) и GIADA (Rotundi и др., 2015; DellaCorte и др., 2015). Анализ результатов наблюдений показал, что приповерхностный слой с очень высокой вероятностью состоит из сложных пористых частиц – агрегатов, размеры которых изменяются в широких пределах: от микронных (Bentley и др., 2015) до, возможно, сантиметровых (Blum и др., 2017). Ядро в целом обладает высокой пористостью ~65–85% (Pätzold и др., 2016; 2019), что накладывает ограничение сверху на пористость частиц-агрегатов. Газовая и пылевая активности наблюдаются практически везде, следуя за солнечным освещением (Keller и др., 2015; Vincent и др., 2019). Это позволяет заключить, что лед находится близко к поверхности (Fornasier и др., 2015; 2019), а связность материала приповерхностного слоя очень низкая (так как пылинки выносятся с ядра под действием потока газа). Современное качественное описание процесса сублимации льда из-под пористого нелетучего слоя выглядит следующим образом. Поскольку альбедо поверхности составляет лишь несколько процентов (Fornasier и др., 2015), почти все солнечное излучение поглощается ядром и за счет теплопроводности, включающей контактную теплопроводность твердой фазы, теплопроводность, обусловленную потоком газа через слой, и радиационную теплопроводность, энергия доставляется к областям, содержащим лед. Лед при нагреве сублимирует, и продукты сублимации, мигрирующие сквозь пористый слой наружу, создают газовую кому. При этом из-за перепада давления могут возникать условия для полного или частичного разрушения пылевой матрицы, и происходит вынос нелетучих частиц, образующих пылевую кому. Согласованная модель переноса энергии и массы в приповерхностном слое, включающая поглощение прямого солнечного излучения, тепловое излучение нагретых пылинок, диффузию пара через пористый слой и изменение его характеристик при прохождении этого слоя, играет первостепенную роль в изучении газо-пылевой активности комет на всех масштабах: от микроскопических до глобальных. В свою очередь, оценка эффективной скорости сублимации с учетом газодинамического сопротивления слоя требует аккуратного изучения транспортных свойств экранирующего поверхностного пористого пылевого слоя.

В данной работе мы исследуем такие характеристики как а) проницаемость слоя для газового потока, б) эффективный размер пор, в) среднее вертикальное смещение между столкновениями молекул с пылевой матрицей, г) распределения глубин для первого и последнего столкновения в слое. Это исследование является продолжением и дальнейшим расширением анализа, представленного в нашей предыдущей работе (Reshetnyk и др., 2021). В цитируемой работе был предложен согласованный подход для количественной оценки транспортных характеристик на основе анализа микроструктуры модельного пористого слоя. Используя этот подход, мы 1) на первом этапе создавали различные модельные слои, 2) затем, применяя метод пробных частиц, количественно анализировали поведение исследуемых характеристик, и наконец, 3) получали аппроксимирующие формулы, позволяющие включить результаты компьютерного моделирования в нестационарные тепловые модели, разработанные для анализа наблюдений (Macher и др., 2019; Skorov и др., 2017; 2020).

До сих пор мы исследовали слои, содержащие только твердые сферы – мономеры одного размера. Были проанализированы однородные и неоднородные слои с внутренними полостями и трещинами, имеющие эффективную пористость в диапазоне, согласующемся с наблюдениями (Reshetnyk и др., 2021; Skorov и др., 2021). Предположение о том, что слой состоит из одинаковых сфер, было одним из основных упрощений, которое делалось в предыдущих работах. Следует почеркнуть, что оно плохо согласуется с наблюдениями (Bentley и др., 2016; Langevin и др., 2016), явно показывающими сложную структуру кометных пылинок. В данной работе мы делаем шаг вперед и впервые представляем результаты исследований, проведенных для слоев, имеющих сложную иерархическую структуру. Говоря об иерархическом пористом слое, мы подчеркиваем именно его многоуровневую организацию: такой слой состоит не из твердых частиц, а из пористых агрегатов, которые также имеют пространственное строение. Используя новые модели слоев, мы следуем общей схеме исследования, применявшейся в наших прежних работах. Там, где это возможно и осмысленно, мы приводим сравнение с результатами, полученными для однородных и неоднородных слоев, построенных из твердых сфер. Это помогает читателю составить целостное представление о предмете исследования, а нам сформулировать общие выводы, описывающие зависимость транспортных характеристик поверхностных слоев кометного ядра от его микроструктуры.

МОДЕЛЬНЫЕ СЛОИ

Как было отмечено во “Введении”, в данной работе мы следуем общей схеме исследования, предложенной в (Reshetnyk и др., 2021). На первом этапе мы строим пористый слой, который используется в дальнейшем в качестве модельной среды. Обычно, модельные пористые слои строятся из простейших сферических непересекающихся частиц (мономеров). Более подробный обзор различных методов построения таких слоев можно найти, например, в (Skorov и др., 2011). Получаемые пористые модельные слои могут иметь как однородную, так и неоднородную структуру, а сферы могут быть либо одного размера (монодисперсные слои), либо двух размеров (бидисперсные слои). Кроме того, в слое могут быть полости (пустоты, размер которых значительно превышает размер частиц) или щели (модельный аналог микроскопических трещин, ширина которых также превышает размер модельных частиц). Все эти типы слоев нами рассматривались и систематически анализировались в (Reshetnyk и др., 2021; Skorov и др., 2021).

Основное упрощающее модельное допущение, которое используется в таких моделях, это допущение о том, что слой состоит из твердых сплошных сфер. Результаты наблюдений показывают, что реальные кометные пылинки гораздо больше похожи на сложные агрегаты, обладающие заметной пористостью и сложным неоднородным строением. Мы адресуем читателя к изображениям захваченных частиц, представленных в (Langevin и др., 2016). Именно такие несферические пористые агрегаты, скорее всего, формируют поверхностный нелетучий слой кометного ядра. Чтобы создать модельные слои, структура которых больше соответствует наблюдениям, мы предлагаем двухшаговую схему. На первом шаге мы создаем агрегаты, состоящие из твердых сферических частиц. Для построения таких комплексов используется баллистическое агрегирование – метод, который много лет успешно и широко применяется для изучения как оптических (например, Skorov и др., 2008; 2010), так и механических (например, Skorov и др., 2016; Reshetnyk и др., 2018) свойств космических пылинок. То, что баллистическое агрегирование, по-видимому, является характерным процессом образования частиц на начальном этапе их роста, было подчеркнуто в (Skorov, Blum, 2012). На втором шаге мы строим из таких агрегатов случайную пористую среду, имеющую иерархическое строение и, как следствие, флуктуации пористости на микромасштабе, сравнимом с размером используемых агрегатов. Для построения модельного слоя в данной работе мы используем два различных подхода, которые описаны подробно ниже.

Прежде, чем перейти к описанию конкретных процедур, мы сделаем одно важное общее замечание. Будем характеризовать слои степенью заполнения ψ и пористостью φ, где первая характеристика – это отношение объема, занятого частицами V_s к общему объему образца V_м, а вторая – это разность между общим объемом и полным объемом частиц, деленная на полный объем образца (V_м – V_s)/V_м или 1 – ψ. Выше отмечалось, что результаты наблюдений позволяют оценить возможный диапазон пористости ядра: средняя пористость не превосходит 85%. Очевидно, в малых масштабах возможно существование как более плотных и даже сплошных микронных и субмикронных частиц, так и очень рыхлых агломератов с пористостью, очень близкой к 100%. Однако для наших задач мы ориентируемся на среднее ожидаемое значение. В этом случае мы получаем важное ограничение на возможное строение пористых баллистических агрегатов. Эффективная пористость иерархического слоя в нулевом приближении может быть записана как (1 – ψ_s ψ_A), где в произведении стоят степени заполнения слоя и агрегата, соответственно. В этом выражении степень заполнения слоя ψ_s можно оценить, заменив агрегаты твердыми сферами подходящего размера. Если принять, что для плотно упакованного монодисперсного пористого слоя характерная степень заполнения составляет около 40%, то степень заполнения агрегатов ψ_A должна иметь близкое к этому значение. В таком приближении эффективная пористость слоя равна (1 – 0.4 × 0.4), что примерно равно наблюдательной оценке. Для достижения меньшей эффективной пористости агрегаты должны быть еще более плотными. Классические баллистические агрегаты создаются по простейшей схеме, в которой мономер, сталкивающийся с кластером-мишенью, прилипает в первой точке касания. В этом случае эффективная степень заполнения агрегата быстро уменьшается по мере его роста и асимптотически приближается к примерно 15%, что заметно ниже требуемого значения. Более плотные агрегаты получаются, если после столкновения с кластером-мишенью мономер может мигрировать, находя более “устойчивое” положение. Такие баллистические кластеры рассмотрены всесторонне в (Shen и др., 2008). Для агрегатов, содержащих от сотен до тысяч мономеров, ожидаемая степень заполнения составляет примерно 40–50% (см. рис. 3 в Shen и др., 2008), что удовлетворяет нашим требованиям. В данной работе мы опираемся на результаты цитированной работы и будем использовать именно такие агрегаты (BAM2) различного размера. Пространственные структуры агрегатов, содержащих различное количество мономеров, могут быть найдены по адресу (http://www.astro.princeton.edu/~draine/agglom.html). Важно подчеркнуть, что процедура баллистического агрегирования носит вероятностный характер: относительное начальное положение мономер-кластер выбирается случайным образом. Поэтому даже в случае фиксированного числа мономеров получаемые агрегаты имеют различную структуру, степень несферичности и даже эффективный размер (см. рис. 1 и 4 в Shen и др., 2008). Авторы цитируемой работы представили результаты для 16 различных реализаций агрегатов одной массы. Мы использовали все доступные данные в нашем исследовании, чтобы оценить влияние микро-свойств используемых агрегатов на изучаемые характеристики модельных слоев. Используя плотные баллистические агрегаты типа BAM2 в качестве модельных примитивов, мы переходим к построению иерархических слоев, применяя два различных подхода.

В первом случае модельный слой создавался с использованием компьютерного пакета YADE (https://www.yade-dem.org). YADE – это комплекс программ, использующий метод дискретных элементов для моделирования движения многих сферических частиц под действием различных физических сил, включая силы взаимодействия между частицами и внешние силы (например, гравитацию). Для расчета ускорения каждой частицы на каждом временном шаге применяются законы Ньютона. Взаимодействие между частицами описывается в рамках модели “мягких сфер” (вязкоупругие сферы), в которой частицы могут частично перекрываться. При этом имитируется деформация частиц во время контакта. Соответствующая степень перекрытия позволяет рассчитывать силы притяжения/отталкивания (в приближении силы реакции нелинейной пружины, зависящей от приведенного модуля Юнга и приведенного радиуса взаимодействующих частиц). В модели мягких сфер несколько частиц могут контактировать друг с другом, что полезно для моделирования движения плотных потоков частиц.

Для создания иерархического слоя агрегаты BAM2 считываются в YADE и группируются в объект, который называется в пакете “clump” (комок). Для моделирования комок представляет собой жесткую совокупность мономеров и динамически ведет себя в соответствии с общим поведением составляющих его мономеров. На начальном этапе комки располагаются в области моделирования так, что их координаты центров масс и пространственная ориентация задаются случайным образом. В этой начальной конфигурации системы нет контактов между комками. После этого начинается моделирование, при котором частицы падают под действием внешнего гравитационного поля. В ходе осаждения возникают контакты между комками. Моделирование заканчивается, когда изменение положения агрегатов становится малым (т.е. система комков достигает квазиустойчивого положения). На последнем временном шаге положение и тип каждого мономера и идентификатор его родительского комка сохраняются в качестве входных данных для дальнейшего моделирования. Пример слоя, построенного из агрегатов BAM2, показан на рис. 1. В данном исследовании мы использовали слои, построенные из кластеров одинаковой массы. Характеристики нескольких иерархических слоев, построенных из агрегатов BAM2 с использованием компьютерного пакета YADE, приведены в табл. 1. В первом столбце приведено название слоя, во втором – размер области моделирования в размерах сфер-мономеров, в третьем – максимальное значение средней пористости из 16 вариантов, соответствующих различным реализациям построения агрегата, в четвертом – среднее число столкновений, испытанных одной пробной частицей, прошедшей через слой.

Основным недостатком изложенного подхода с контролем контактов между элементами для нас является его длительное время, затрачиваемое на компьютерное моделирование. Поиск возможных контактов между частицами и достижение квазистационарной финальной конфигурации слоя – это затратные по времени процессы. Как видно из табл. 1, характерное число мономеров в модельном слое составляет сотни тысяч сфер. Очевидно, что прямой перебор при поиске возможных контактов крайне затруднителен и времязатратен. Для того чтобы избежать возникающих вычислительных трудностей, мы использовали в работе второй подход, который мы будем называть методом без контроля контактов. В основе этого подхода лежит простая идея, предложенная нами в (Skorov и др., 2016; Reshetnyk и др., 2018) для изучения оптических и динамических свойств иерархических пористых агрегатов.

В простейшем случае при применении метода без контроля контактов модельный слой конструируется из непересекающихся псевдомономеров, каждый из которых является баллистическим агрегатом типа BAM2 с определенным количеством сфер-мономеров. Очевидно, что при реализации этого метода построения в общем случае агрегаты не имеют контактов друг с другом. К сожалению, средняя пористость получаемого слоя обычно превышает 90%, что выходит за рамки интересующего нас диапазона значений пористости. Поэтому этот простейший вариант не использовался в нашем исследовании. Чтобы получить слой с более низкой пористостью, следует обратить внимание на два факта. Во-первых, используемые агрегаты BAM2 имеют неоднородное распределение плотности: по мере удаления от центра массы эффективная плотность агрегата падает, т.е. пористость возрастает. Анализ показал, что, выделяя центральную часть агрегата, мы можем получить объекты со средней пористостью около 40%. В таком случае пористость слоя оказывается в диапазоне наблюдаемых значений. Во-вторых, анализ, выполненный в (Shen и др., 2008), показал, что форма таких агрегатов может довольно заметно отличаться от сферической (на десятки процентов). Таким образом, помещая их внутрь контрольных сфер (псевдомономеров), мы избыточно увеличиваем пористость слоя. Если искусственно увеличить размер каждого псевдомономера, допуская их перекрытие, то пористость слоя снизится, а число перекрытий между реальными мономерами-сферами, входящими в агрегаты, скорее всего, будет небольшим. Мы увеличивали размер псевдомономера на один радиус и на три радиуса сферы-мономера. Используя этот простой прием, можно получить слои с результирующей пористостью около 70–80%. При этом доля перекрытий между мономерами обычно не превосходит нескольких процентов, что не оказывает заметного влияния на значение эффективной пористости слоя. На правой панели рис. 1 схематично показано, как строились слои без контроля за контактами между мономерами. В табл. 1 в качестве иллюстрации приведены базовые характеристики нескольких реализаций таких слоев, построенных из баллистических агрегатов BAM2.

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

МОДЕЛИРОВАНИЕ ПЕРЕНОСА ГАЗА И ИЗЛУЧЕНИЯ В СЛОЕ

Целью данной работы, также как и предыдущего исследования (Reshetnyk и др., 2021), является получение количественных оценок физических транспортных характеристик пористого слоя, имеющих важное прикладное значение. Эти характеристики должны быть использованы в моделях, описывающих тепло- и массоперенос в приповерхностной области ядра кометы. В свою очередь, эти модели должны помочь нам в анализе данных, полученных в ходе наблюдений in-situ. Прежде всего, речь идет об изучении газовой и пылевой активности ядра. Эта прикладная цель определяет выбор используемых модельных подходов и анализируемых характеристик. Последние, как и прежде, мы подразделяем на две группы: относящиеся к прохождению через слой продуктов сублимации льда (т.е. молекул разреженного газа) и относящиеся к распространению в слое излучения (как прямого солнечного, так и теплового излучения нагретых солнцем пылинок). К первой группе относятся: проницаемость слоя Ψ, среднее вертикальное смещение между последовательными столкновениями Δ, и распределение глубин, на которых имело место последнее перед вылетом из слоя столкновение пробной частицы с пылевой матрицей Λ(z). Во вторую группу входят средний размер длин свободного пробега в слое Ω и распределение глубин, на которых имело место первое столкновение пробной частицы с пылевой матрицей Θ(z). Результаты анализа этих характеристик и их влияния на перенос энергии и массы в иерархических пористых слоях представлены в следующих разделах.

Следуя нашим предшествующим работам, мы будем использовать метод трассировки лучей (или так называемый метод пробной частицы) для обеих групп. Обоснование применения данного подхода и ограничения его использования были нами представлены в (Reshetnyk и др., 2021). Здесь отметим лишь, что в кометных условиях обычно можно рассматривать течение газа в приповерхностном пористом слое как безстолкновительное, т.е. мы пренебрегаем столкновениями между молекулами, находящимися внутри слоя. На основании наших расчетов и данных о размерах пылинок (Hilchenbach и др., 2016), можно заключить, что такое упрощение не является сильной идеализацией и физически обосновано. Для анализа характеристик переноса излучения применение метода пробных частиц является значительно более сильным упрощением и возможно лишь в случаях, где применим геометрический подход.

Oписание метода пробных частиц и его применение к рассматриваемой задаче можно найти, например, в (Skorov и др., 2011). Краткое описание схемы моделирования представлено ниже. На первом шаге координаты центров всех сфер-мономеров, составляющих модельный слой, загружаются в память компьютера. После определения расчетной области, имеющей форму кубоида, мы приступаем к моделированию движения пробной частицы. Движение частицы начинается на нижней грани среды. Координаты точки входа в кубоид и значения для направляющих косинусов в соответствии с предполагаемой функцией распределения пробных частиц задаются случайным образом. При изучении диффузии газа задается так называемая полумаксвелловская функция распределения скоростей молекул. Для моделирования ослабления интенсивности прямого солнечного излучения считается, что все пробные частицы изначально движутся в одном направлении по внутренней нормали к нижней грани кубоида. Отметим, что для оценки интересующих нас характеристик (Ψ, Δ, Λ(z)) не требуется генерировать значения абсолютной величины скорости. После задания начальных условий ищется ближайшая сфера, с которой пересекается линия движения пробной частицы, и находятся координаты столкновения. Затем, в предположении, что имеет место диффузное рассеяние, генерируются новые значения направляющих косинусов. После акта рассеяния процедура повторяется до тех пор, пока пробная частица не пересечет одну из граней кубоида. Если пробная частица пересекает горизонтальные грани, то значения счетчиков испускаемых или поглощаемых частиц изменяются, и мы переходим к следующей тестовой частице. Если частица пересекает вертикальные боковые грани кубоида, то в предположении о периодических граничных условиях задаются новые координаты, и моделирование движения продолжается. Описанная процедура продолжается до тех пор, пока счетчик частиц, вылетевших через верхнюю грань кубоида, не достигнет заданного значения (обычно, 100 тысяч частиц).

СТРУКТУРНЫЕ ХАРАКТЕРИСТИКИ СЛОЯ

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

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

Для описания пористой среды существует множество подходов, начиная с простейших, в которых пористая среда описывается как система прямых цилиндрических каналов – капилляров, и заканчивая сложными геометрическими моделями, построенными, например, на основе компьютерного анализа трехмерных микроизображений реальных пористых сред (Christou и др., 2020). Мы будем описывать пористый слой как множество элементарных сфер, как это было сделано в (Skorov и др., 2011; 2021; Reshetnyk и др., 2021). В этих работах можно найти краткое описание достоинств и недостатков капиллярных и гранулированных моделей. Важной чертой моделей последнего типа является то, что движение молекул становится “внешним” по отношению к пылевой матрице, и пустоты между твердыми частицами образуют в общем случае единую глобальную сеть. В цитируемых работах мы рассматривали однородные и неоднородные слои, составленные из сфер, которые не пересекались друг с другом. Были рассмотрены случаи, когда в слое имеются полости или трещины, характерные размеры которых заметно превосходят размеры сфер-мономеров. В данном исследовании мы делаем следующий важный шаг и рассматриваем слои, построенные из сложных пористых агрегатов. Наша основная цель заключается в изучении эффектов, обусловленных таким сложным строением пористой среды, и получении аппроксимирующих функций для транспортных характеристик. Для сравнения мы будем использовать результаты, полученные в (Reshetnyk и др., 2021).

Общая схема предлагаемого анализа следующая.

– На первом шаге мы рассчитываем и анализируем функцию распределения длин отрезков между двумя последовательными столкновениями в слое Ξ. Поскольку в модели отсутствуют межмолекулярные столкновения, то это распределение можно называть распределением длин пробегов молекул.

– На основе этих распределений, мы вычисляем средние значения длин пробегов, которые используются в качестве характеристики, задающей “эффективный размер пор” Ω. Важно подчеркнуть, что мы будем говорить о “порах”, хотя, как отмечалось выше, в наших моделях пустоты образуют связанную сеть.

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

Результаты расчетов функций распределения длин пролетов и соответствующих средних значений показаны на рис. 2. Мы использовали баллистические агрегаты типа BAM2, содержащие различное количество сфер (верхний ряд N = 128; средний ряд N = 256; нижний ряд N = 512). Для каждого числа мономеров в агрегате мы проводили расчеты для 16 различных вариантов его построения. Как отмечалось выше, сам процесс баллистического агрегирования вносит статистический разброс в свойства получающегося агрегата. Этот эффект хорошо виден на правых панелях, где показаны значения средней длины пролета Ω. Для фиксированного числа мономеров в агрегате наблюдается изменение его эффективной пористости для различных вариантов построения (см. Shen и др., 2008). Поэтому мы рассчитывали функцию Ξ для всех доступных вариантов построения агрегата. Результаты расчетов представлены на левых панелях рис. 3. Видно, что поведение функций распределения длин пролетов качественно не отличаются при увеличении размера агрегата. По сравнению с распределениями длин, полученными для однородных слоев, построенных без учета контактов между мономерами (см. рис. 2 левый столбец в Reshetnyk и др., 2021), для иерархических слоев не наблюдается уменьшение значений функции распределения в области малых длин пролетов. Эта особенность обсуждалась нами, и было показано, что она связана именно со способом построения слоя. Для иерархических слоев функция распределения длин есть монотонно убывающая функция для всех значений пролетов. Во всех рассмотренных случаях изменчивость функции для различных вариантов построения агрегата увеличивается при величинах пролетов, составляющих около 0.4–0.6 размеров сферы-мономера. В этой области достигает максимума функция распределения длин пролетов в однородном монодисперсном слое (см., например, рис. 2 в Reshetnyk и др., 2021). По мере увеличения длины пролета функция распределения падает, как и в случае однородного слоя, но темп падения заметно меньше – возрастает доля длинных пролетов. Эти длинные пролеты соответствуют уже столкновениям, которые происходят при движении молекулы между агрегатами, а не внутри одного кластера. Как и следовало ожидать, относительная изменчивость функции Ξ увеличивается по мере роста длины пролета, и в области порядка десяти размеров мономера разброс значений может составлять несколько сотен процентов. Такое поведение наблюдается для всех рассмотренных значений числа сфер-мономеров N. Напомним, что эта изменчивость обусловлена только реализацией процедуры построения агрегата, т.е. она имеет стохастическую природу.

Все изученные слои были построены с контролем контактов между агрегатами (см. “Структурные характеристики…”), поэтому можно ожидать относительно небольших вариаций средней пористости слоя при использовании агрегатов фиксированной массы. Расчеты подтверждают это предположение: на рис. 1 в правом столбце по оси абсцисс показаны значения средней пористости слоев, построенных из различных реализаций кластера, содержащего заданное число сфер-мономеров. Видно, что изменения пористости составляют всего лишь 3–5%. Эти изменения связаны с различием в форме (например, степени несферичности) используемых при построении слоя агрегатов. Вместе с тем, интересно отметить, что среднее значение длины пробега изменяется в зависимости от пористости в очень хорошем согласии с теоретической зависимостью Ω ~ φ/(1 – φ) (Dullien, 1991). Эта же зависимость, как было показано в (Reshetnyk и др., 2021), выполняется для однородных слоев и слоев, содержащих неоднородности (полости или трещины). Значение коэффициента пропорциональности равно примерно 2/3 и не изменяется в пределах точности вычислений при увеличении размера агрегата.

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

На рис. 3 показаны функции распределения длин пробега (левая панель) и средние длины пробега (правая панель) для иерархических слоев, построенных из псевдомономеров без контроля за контактами между агрегатами. Мы использовали агрегаты различного размера (исходное число сфер-мономеров в агрегате N = 512, 1024) и слои различной пористости (см. “Моделирование переноса…” о повышении степени заполнения в таких слоях). Для сравнения приведены результаты, полученные для однородных монодисперсных слоев, имеющих пористость 75 и 85%, и слоя с полостями, имеющего среднюю пористость 71%.

Хорошо видно, что в иерархических слоях относительная доля дальних пролетов выше, чем в однородных слоях, имеющих более высокую пористость. Очевидно, это связано с наличием более крупных по размерам пустот между агрегатами. Напомним, что степень заполнения самих агрегатов заметно превышает значения степени заполнения для однородных слоев. Как уже отмечалось выше для слоев, построенных с контролем за контактами между агрегатами (см. рис. 2), в иерархических слоях доля коротких пролетов, которые в подавляющем большинстве имеют место внутри агрегата, во много раз превосходит соответствующие значения, рассчитанные для однородных слоев. Эта особенность связана с методом построения однородных слоев. Косвенным доказательством служит близость функций распределения длин пролетов между иерархическими слоями и слоем, содержащим полости, вплоть до значений пролетов порядка размера сферы-мономера. Поскольку размер и форма пустот между агрегатами в иерархическом слое имеют более случайный характер, чем в случае однородного слоя с полостями, имеющими одинаковый размер и сферическую форму (см. подробное описание построения такого слоя в Reshetnyk и др., 2021), вторичные максимумы в функциях распределения в иерархических слоях исчезают. Мы уже не видим характерного пика функции для пролетов с длиной, близкой к размеру полости в неоднородном слое.

Если сравнивать функции распределения в области длинных пролетов для иерархических слоев, то можно заметить, что пористость слоя играет определяющую роль: кривые, построенные для слоев, имеющих близкие значения средней пористости (см. табл. 1), образуют близкие пары. Интересно сравнить значения средних длин пролетов, приведенные на правой панели. Показаны результаты, полученные для двух исходных размеров агрегатов (512 и 1024), из которых выделена центральная часть с меньшей пористостью. Затем мы понижали пористость слоя, увеличивая размер псевдомономера на один или три радиуса сферы-мономера. Были выполнены расчеты для всех имеющихся 16 вариантов агрегатов BAM2 заданной массы. Видно, что образуются два “облака” значений пористости: более низкие значения соответствуют слоям, в которых псевдомономеры могут перекрываться на три радиуса сферы-мономера, а более высокие значения пористости – слоям, в которых псевдомономеры могут перекрываться только на один радиус сферы-мономера. Изменчивость пористости внутри “облака” связана с использованием различных реализаций исходных агрегатов BAM2. Важным и нетривиальным является тот факт, что для всех рассмотренных модельных слоев по-прежнему очень хорошо работает аппроксимация вида Ω ~ φ/(1 – φ). Значение коэффициента пропорциональности равно примерно 2/3 и не изменяется в пределах точности вычислений.

Подводя итог анализу представленных результатов, можно заключить, что

– функция распределения длин пролетов для иерархических слоев имеет монотонно убывающий характер;

– форма агрегата оказывает влияние на эту функцию, и для пролетов, сравнимых с размером агрегата, вариации могут составлять сотни процентов;

– по сравнению со слоями, построенными из твердых сфер, наблюдается относительное увеличение доли длинных пролетов (т.е. пролетов, превышающих размер агрегата);

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

Перейдем теперь к анализу результатов моделирования проницаемости иерархических слоев. Рассмотрим сначала слои, построенные методом осаждения баллистических агрегатов во внешнем гравитационном поле (рис. 4, левая панель). Как отмечалось выше, в моделях этого класса агрегаты могут вращаться после образования первого контакта. Кроме того, форма агрегатов отличается от сферической (см. рис. 1 в Shen и др., 2008). Это приводит к вариациям средней пористости слоев и функции распределения длин пролетов (см. рис. 2). В результате изменяется эффективная проницаемость слоя: для наиболее крупных из рассмотренных частиц вариации могут составлять 20–30%. Разница в проницаемости между слоями, построенными из агрегатов, содержащих 128 и 512 сфер, отличается примерно вдвое для толщин, превосходящих размер агрегата. Отметим, что средняя пористость этих агрегатов отличается не больше чем на 10% (большие агрегаты имеют большую пористость). Это означает, что в увеличении проницаемости слоев, содержащих крупные кластеры, играет роль увеличение среднего размера пустот между кластерами. Различия в средних значениях длин пробегов (примерно 1.8 для BAM2.128 и 2.6 для BAM2.512) подтверждают это предположение. В то же время общее поведение проницаемости имеет такой же характер, как в случаях однородных и неоднородных слоев, состоящих из сфер-мономеров (см. рис. 7 и 8 в Skorov и др., 2021). Это позволяет надеяться, что простая аппроксимация, использованная нами ранее для оценки проницаемости, будет применима и для этих слоев.

На правой панели рис. 4 показаны результаты для слоев, сконструированных из псевдомономеров. В этом случае несферичность используемых агрегатов играет еще большую роль. Псевдомономер – это сфера, содержащая наш баллистический агрегат, поэтому размер этой сферы может заметно превышать так называемый эффективный размер агрегата, пропорциональный размеру сферы-мономера и корню кубическому из числа мономеров (Shen и др., 2008). Для таких слоев, в которых в общем случае нет контактов между агрегатами, проницаемость изменяется в разы при использовании различных вариантов агрегатов, состоящих из заданного числа мономеров. На рис. 4 такие области изменения проницаемости показаны различными видами штриховки. Для слоев, толщина которых превышает размер агрегата, проницаемость может меняться в 5–8 раз в рассмотренных случаях. Важно отметить, что средняя пористость изменяется от примерно 70% (512s-3r) до примерно 90% (1024s-r), т.е. степень заполнения этих слоев отличается также значительно: от 30 до 10%. Более чем вдвое отличаются значения средней длины пролета: они равны примерно 1.6 и 4.1 размера мономера, соответственно. Интересно, что в значительном интервале значений толщин проницаемость слоев, содержащих более мелкие агрегаты, может превышать проницаемость слоев, построенных из более крупных, но более плотных агрегатов. Расчеты показали, что, несмотря на значительный разброс значений проницаемости при заданной толщине слоя, для всех рассмотренных случаев общее поведение зависимости проницаемости от толщины слоя сохраняет прежний вид: проницаемость монотонно убывает с увеличением толщины слоя. Опираясь на этот результат, для этого класса слоев мы будем искать аппроксимирующую функцию в том же виде, что и для слоев, построенных с контролем контактов между агрегатами.

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

где h = L/D_m – это безразмерная толщина слоя, выраженная в размерах сфер-мономеров.

Выбирая такой вид функции, мы хотим оценить, насколько хорошо простейшая формула описывает рассматриваемые в данной работе иерархические слои, и проанализировать, как сильно изменяются значения найденного коэффициента b для различных слоев. Таким образом, мы хотим проверить, насколько универсальной является формула, полученная нами в (Gundlach и др., 2011) на основе анализа конкретных экспериментальных данных. На рис. 5 в качестве иллюстрации приведена проницаемость как функция толщины слоя для слоев, составленных из агрегатов BAM2 (N = 512), в классе моделей с контролем контактов (левая панель) и для слоев, составленных из псевдомономеров двух размеров (правая панель). Пористость слоев на левой панели примерно равна 80%, пористость слоев на правой панели равна примерно 76% (группа слоев 512s-r) и 85% (группа слоев 1024s-3r). Показаны аппроксимирующие функции, со средними значениями коэффициента b, рассчитанного для всей выборки реализаций слоев (т.е. 16 вариантов для каждого числа мономеров в агрегате).

Можно отметить, что тестируемая аппроксимирующая формула начинает хорошо работать, начиная с толщин слоя, превосходящих размер агрегата (т.е. с толщин, больших примерно 5 размеров сфер-мономеров). Это особенно хорошо видно на левой панели, где разброс значений проницаемости для слоев заданной толщины меньше. Такой вывод не является неожиданным, так как наши иерархические слои имеют высокую степень неоднородности на масштабе, сравнимом с размером агрегата. Кроме того, на проницаемость тонких слоев, оказывает влияние граничное условие, используемое в модели осаждения агрегатов. Для результатов, представленных на правой панели, разброс значений проницаемости намного больше. Распределение агрегатов вблизи нижней грани образца более однородное для класса моделей слоев, построенных из псевдомономеров. Поэтому область, в которой аппроксимирующая формула работает удовлетворительно, несколько шире. Для всех рассмотренных случаев предлагаемая аппроксимация хорошо работает для слоев толщиной в несколько размеров агрегатов. Средние значения коэффициента b и ошибки вычисления для рассмотренных слоев и соответствующие значения средних длин пролетов приведены в табл. 2. Напомним, что усреднение выполнялось для выборки из всех слоев, содержащих агрегаты, построенные из одинакового числа сфер-мономеров (т.е. 16 вариантов для каждого сорта баллистического агрегата).

Приведенные результаты доказывают, что при использовании предложенной в (Gundlach и др., 2011) простейшей формулы для аппроксимации проницаемости необходимо учитывать среднюю длину пробега в слое, которая, в свою очередь, нелинейно зависит от средней пористости слоя и от размера агрегата (см. правые панели на рис. 2 и 3). Таким образом, мы заключаем, что значение коэффициента, приведенное в цитируемой работе (где лучшее совпадение аппроксимирующей функции с результатами эксперимента для образцов лед–пыль было получено при значении коэффициента b, равном 7.31 ± 0.20, для слоя со средней пористостью чуть выше 40%), можно использовать только для слоев, имеющих структурные микрохарактеристики, близкие к характеристикам экспериментальных образцов. Для слоев, имеющих неоднородности (полости или трещины), а так же для слоев, имеющих иерархическую структуру, полученная зависимость не применима. Несмотря на значительный разброс в значениях коэффициента в формуле для аппроксимации проницаемости, предлагаемый подход последователен и легко применим на практике. Мы показали, что среднее значение длины пролетов есть однозначная функция средней пористости слоя. Также мы показали, что для относительно толстых слоев (толщиной более нескольких размеров агрегата) при известном значении средней длины пролета хорошо работает простейшая аппроксимационная формула, имеющая лишь один модельный коэффициент. Эта схема позволяет использовать данный подход для описания газодинамического сопротивления пористого слоя, имеющего иерархическую структуру, в моделях тепломассопереноса в приповерхностной области кометного ядра, что и является одной из конечных целей данного исследования.

В конце этого раздела мы коротко остановимся на результатах анализа модельного параметра извилистости τ. Эта характеристика пористой среды обычно используется при описании пористости как связки цилиндрических капилляров. Такой подход традиционно применялся в теплофизических моделях кометного ядра (см., например, Mekler и др., 1990; Enzian и др., 1997; Skorov и др., 1999; 2001). Во всех этих моделях используется формула для проницаемости, предложенная в (Fanale, Salvail, 1984), в которой проницаемость среды обратно пропорциональна квадрату извилистости. Более подробно мы рассматривали использование данного подхода в (Skorov и др., 2011; Reshetnyk и др., 2021). Было показано, что в моделях, где извилистость определяется как отношение средней длины пролета к среднему вертикальному смещению в направлении диффузии между последовательными столкновениями, ее значение близко к 2. Такое значение извилистости получается теоретически в модели цилиндрических пор с изотропным распределением их ориентации в пространстве (Markin, 1965). Систематические расчеты, выполненные для иерархических слоев различных классов, пористость которых изменялась от примерно 70% до примерно 90%, показали, что значения извилистости по-прежнему изменяются в окрестности значения 2. В качестве иллюстрации на рис. 6 мы приводим результаты расчетов для слоев двух классов: 1) построенных с контролем контактов между агрегатами (левая панель), 2) построенных из псевдомономеров (правая панель). Хорошо видно, что как в случае моделей однородных слоев и слоев со щелями извилистость слабо растет с ростом пористости (Reshetnyk и др., 2021). Лишь для самых пористых случаев (>90%) ее значение превышает теоретическое значение более чем на 10%. Мы заключаем, что можно не учитывать незначительные изменения извилистости с изменением пористости и/или размеров агрегата. Это означает, что в моделях, описывающих пористую среду в виде связки параллельных цилиндрических капилляров, учет извилистости сводится к умножению на константу, как это предполагалось в пионерской работе (Fanale, Salvail, 1984).

КИНЕТИЧЕСКИЕ ХАРАКТЕРИСТИКИ ВЫЛЕТЕВШИХ МОЛЕКУЛ

Следуя структуре исследования, предложенной в (Reshetnyk и др., 2021), мы переходим к анализу кинетических характеристик молекул, вылетевших из иерархического слоя. Как было показано, например, в (Macher и др., 2019; Skorov и др., 2017; 2020), из-за крайне низкой тепловой инерции, разница температур на поверхности ядра (т.е. на верхней границе пористого слоя) и на фронте сублимации льда (т.е. на нижней границе слоя) может достигать примерно двухсот кельвинов. Это означает, что большое значение имеет определение глубины, на которой молекула испытала последнее столкновение перед вылетом из слоя. Именно эта характеристика в первом приближении определяет тепловую скорость вылетевшей молекулы (в модели с диффузным рассеянием с полной тепловой аккомодацией). Оценка этой величины требуется, например, для расчета ускорения пылинок потоком газа вблизи ядра (Skorov и др., 2016; Reshetnyk и др., 2018) и для интерпретации наблюдений распределения температуры во внутренней коме, сделанных инструментом MIRO (Marschall и др., 2019).

Как было отмечено в (Reshetnyk и др., 2021), в общем случае для получения аккуратной оценки значений температуры вылетевших молекул необходимо детально описать акт элементарного столкновения молекулы и частицы слоя. Это предполагает, например, введение в рассмотрение коэффициента тепловой аккомодации при столкновении. Стараясь сохранить максимально возможную простоту модели, мы предполагаем, как и прежде, что тепловая скорость молекулы после столкновения полностью определяется температурой пылинки в точке столкновения. Тогда нахождение распределения абсолютного значения скорости вылетевших молекул сводится к определению распределения глубины последнего столкновения в слое перед вылетом тестовой частицы. В этом случае мы получаем максимальный нагрев газа более горячим пылевым слоем, т.е. полученный результат является оценкой сверху. Так как в модели зеркального рассеяния сохраняется абсолютное значение скорости после столкновения с частицами слоя, то эта модель рассеяния дает нам оценку снизу на величину скорости вылетевших молекул. Более сложные модели рассеяния дадут значения, лежащие между этими предельными значениями.

Результаты расчетов относительной доли вылетевших из слоя частиц, испытавших последнее столкновение на данной глубине, показаны на рис. 7. Моделирование выполнялось для иерархических слоев, построенных из баллистических агрегатов BAM2, содержащих различное количество сфер-мономеров (N = 128, 256, 512). При построении слоев контролировались контакты между агрегатами. Для всех рассмотренных случаев относительная доля вылетевших частиц сначала возрастает по мере приближения к верхней границе слоя, а затем падает непосредственно вблизи верхней границы слоя. Как можно ожидать, общее поведение кривых не отличается от поведения, полученного для монодисперсных слоев (см. рис. 7 в Reshetnyk и др., 2021). Как и ранее, с увеличением пористости слоя возрастает относительная доля частиц, испытавших последнее столкновение на заданной глубине вне самой верхней зоны слоя. Вблизи верхней границы слоя всегда наблюдается область падения относительного числа вылетевших молекул. Мы уже объясняли этот эффект асимметричным распределением точек столкновений по поверхности рассеивающей сферы вблизи верхней грани нашего слоя (Reshetnyk и др., 2021). Для иерархических слоев размер этой зоны заметно меньше, чем для монодисперсных слоев. Почти во всей области глубин для функции распределения удовлетворительной является простейшая экспоненциальная аппроксимация (показана линиями на левой панели). При использовании такой аппроксимации мы вводим параметр эффективной глубины вылета, как глубины, на которой число вылетевших частиц уменьшается в e раз. Значения соответствующих эффективных глубин показаны на правой панели. Для каждого числа мономеров в агрегате (N = 128, 256, 512) мы проводили расчеты для всех имеющихся 16 реализаций построения баллистического агрегата. Полученные значения эффективной глубины (показаны маркерами различного типа) образуют три выраженных кластера, каждому из которых соответствует свой интервал изменений пористости. Это говорит о том, что, как и в случае анализа проницаемости, пористость позволяет оценить эффективную глубину вылета молекулы, а значит (с учетом сделанных модельных упрощений) и абсолютное значение ее скорости. Прямое сравнение результатов, полученных для иерархических и однородных монодисперсных слоев, имеющих близкие значения средней пористости, показывают, что в первом случае эффективная глубина вылета примерно вдвое больше, что обусловлено присутствием в иерархических слоях больших по размеру полостей между агрегатами. Такие пустоты повышают вероятность вылета из более глубоколежащих слоев.

Используя параметр эффективной глубины вылета, мы можем легко оценить толщину подслоя, из которого вылетает две трети всех молекул. Учитывая, что в дневное время, когда поверхность освещена Солнцем, температура приповерхностных областей может значительно превышать температуру сублимирующего льда, весовой вклад этих молекул становится определяющим в температуру газа, прошедшего через неизотермический пористый слой. Следуя нашей схеме исследований (Reshetnyk и др., 2021), мы предполагаем, что тепловой поток, переносимый газом, очень мал, а время диффузии много меньше времени, при котором распределение температуры изменяется существенно. Тогда можно считать, что температура в пылевом слое есть линейная функция глубины и, задавая значения граничных температур пылевого слоя, можно получить весовое распределение температуры вылетевших молекул и среднее значение температуры газа. Детали схемы расчетов можно найти в цитируемой работе. В качестве иллюстрации рассматривались слои, в которых температура нижней границы равна температуре сублимирующего водяного льда и равна 200 К, а температура поверхности – 300 К. Выбранная температура льда является репрезентативной для гелиоцентрических расстояний, где наблюдается заметное газопроизводство водяного льда (т.е. <2–3 а. е.; см., например, Weissman и Kieffer (1981), Fanale и Salvail (1984)). Примеры модельных расчетов приведены в табл. 3. В первом столбце приведено название слоя, во втором столбце – соответствующее значение пористости. В следующих столбцах приведены средние температуры газа для трех значений толщины слоя. Как видно из табл. 3, во всех рассмотренных случаях температура газа заметно превышает температуру сублимирующего льда. Этот эффект может иметь значение для анализа распределений скоростей и температур газа вблизи поверхности ядра. Кроме того, увеличение средней температуры газа, очевидно, приводит к росту силы тяги, действующей со стороны газового потока на пылевые частицы, ускоряемые в зоне связанного двухфазного течения (Skorov и др., 2016; Reshetnyk и др., 2018). Таким образом, нами показано, что средняя температура вылетевших молекул изменяется в зависимости от разницы температур на границах слоя, его толщины и пористости.

Кроме температуры для описания вылетевших молекул имеет значение их угловое распределение (Skorov и Rickman, 1995). Для однородных и неоднородных слоев, сконструированных из твердых сфер, результаты расчетов были представлены в (Reshetnyk и др., 2021). Для иерархических слоев угловое распределение очень мало отличается от распределений, рассчитанных ранее. Как и прежде, средние значения угла между направлением вектора скорости молекулы и нормалью к поверхности составляли 40°–41°. Небольшие различия наблюдались между моделями с зеркальным и с полностью диффузным рассеянием. Эти различия имели такой же порядок величин, что и различия, полученные для разных реализаций слоев. Таким образом, можно заключить, что во всех рассмотренных моделях наблюдается незначительная “фокусировка” вылетевших молекул: их скорости в большей степени направлены по нормали к поверхности. Моделирование показало, что изменение абсолютной скорости при прохождении неизотермического слоя имеет намного большее значение, чем изменение углового распределения молекул (как при диффузном, так и при зеркальном рассеянии).

Как было отмечено в начале этого раздела, модель рассеяния молекул на пылевых частицах может заметно влиять на результаты. Для оценки этого влияния мы провели сравнение результатов моделирования для слоев с зеркальным и диффузным рассеянием. В качестве тестовой характеристики была выбрана проницаемость иерархического слоя. Полученные результаты отличаются заметно для всех рассмотренных размеров агрегатов, типов слоев и их толщин. Иллюстрирующие примеры приведены на рис. 8. В верхнем ряду показаны результаты, полученные для иерархических слоев, построенных из баллистических агрегатов типа BAM2, содержащих различное количество сфер-мономеров: слева – модель слоя с контролем контактов между агрегатами, справа – модель слоя без контроля контактов между агрегатами. Рассматривались полностью диффузная и зеркальная модели рассеяния. Относительная проницаемость слоя показана как функция его безразмерной толщины, выраженной в размерах сфер-мономеров. В нижнем ряду показаны отношения рассчитанных проницаемостей к проницаемости однородного слоя с пористостью 85% в случае зеркальной модели рассеяния.

Видно, что для иерархических слоев, построенных с контролем контактов между агрегатами (левый столбец), заметные различия в проницаемости наблюдаются для всех рассмотренных размеров агрегатов и для всех толщин слоя, примерно превышающих размер агрегата. Схожие результаты получены для слоев, построенных из псевдомономеров (правый столбец). Эти слои более неоднородные в малом масштабе, что проявляется в легком волнообразном поведении проницаемости. Тем не менее отличия между двумя моделями рассеяния остаются хорошо заметными. Различие в проницаемостях для иерархических слоев с различными моделями рассеяния составляет 20–30% почти во всем диапазоне рассмотренных толщин. Для их качественной характеристики мы построили отношения рассчитанных проницаемостей к проницаемости однородного монодисперсного слоя с пористостью 85% (нижний ряд). Также мы привели отношение, полученное для идентичных однородных слоев с различными моделями рассеяния. Поскольку средние пористости иерархических слоев (см. табл. 1) отличаются от пористости контрольного однородного слоя с пористостью 85%, к сравнению нужно подходить осторожно. В то же время эти результаты хорошо иллюстрируют ожидаемую степень неопределенности в оценке проницаемости слоя, которая может быть вызвана как локальными изменениями пористости, так и неустранимой сейчас неопределенностью в моделировании, обусловленной выбором модели рассеяния молекул на пылевых частицах.

ОСЛАБЛЕНИЕ ПРЯМОГО ИЗЛУЧЕНИЯ В ПОГЛОЩАЮЩЕМ СЛОЕ

В заключительной части нашего исследования мы представим результаты моделирования, относящиеся к оценке ослабления прямого солнечного излучения в пористом иерархическом слое. Гипотеза о том, что солнечное излучение поглощается не на поверхности кометного ядра, а в некотором эффективном объеме – приповерхностной области пористого слоя, была высказана в (Davidsson, Skorov, 2002; 2004) и позже получила подтверждение в лабораторных экспериментах (Kaufmann и др., 2006; Kaufmann, Hagermann, 2017). В (Reshetnyk и др., 2021) мы обсуждали этот процесс и его важность для моделирования переноса энергии и массы в приповерхностной области кометного ядра. Мы не будем повторять здесь проделанный анализ, а перейдем к краткому описанию модели и обсуждению полученных результатов. Как и прежде, предполагается, что пылевая среда только поглощает солнечное излучение, т.е. процессы рассеяния не рассматриваются. Мы также не рассматриваем собственное тепловое излучение пылевых частиц. Эти модельные упрощения обусловлены конечной целью данного исследования, которая состоит в получении профиля поглощения падающего солнечного света, который можно использовать в моделях тепломассопереноса для изучения макроскопических процессов кометной активности. Для моделирования ослабления прямого излучения используется метод пробных частиц, примененный выше для изучения свободно-молекулярной диффузии. В качестве аппроксимирующей функции мы будем использовать экспоненциальную зависимость, предполагая выполнение закона Бугера–Ламберта–Бера.

Нами были проведены систематические расчеты для модельных иерархических слоев, построенных из баллистических агрегатов разного размера. Исследовались слои, сконструированные с учетом и без учета контактов между агрегатами. Пористость слоев изменялась от примерно 70 до 90%, т.е. был рассмотрен весь ожидаемый диапазон значений пористости ядра кометы Чурюмова–Герасименко. Моделирование выполнялось следующим образом. На горизонтальной контрольной поверхности (например, на нижней границе слоя) случайным образом задавались начальные координаты тестовой частицы. Все частицы двигались вертикально вверх, поскольку солнечное излучение рассматривалось как плоскопараллельный пучок. Определялись координаты первого столкновения с пылевыми частицами. После этого процедура повторялась для следующей тестовой частицы. В результате многократного повторения этой процедуры была получена функция распределения длин первого пролета.

Результаты расчетов для иерархических слоев различной пористости приведены на рис. 9 (левая панель – слои, построенные с контролем контактов между агрегатами, правая панель – слои, построенные без контроля контактов). На этом же рисунке приведены аппроксимационные кривые, коэффициенты которых приведены в табл. 4. Хорошо видно, что по мере увеличения пористости (т.е. уменьшения степени заполнения) излучение эффективно проникает в более глубокие области.

Как и следовало ожидать, полости между агрегатами позволяют прямому излучению проникать в более глубокие области иерархического слоя, чем в случае монодисперсных и бидисперсных слоев, рассмотренных ранее. Если для однородных слоев около 30% прямого солнечного излучения проникало на глубину в 3–4 размера пылевой частицы, а полная толщина поглощающего слоя составляет около 10 размеров сфер-мономеров, то для иерархических слоев соответствующие значения в 3–4 раза больше при сравнимых значениях средней пористости. Теперь уже около 10% прямого излучения сохраняется на глубинах 15–30 размеров сфер-мономеров в зависимости от строения слоя. Интересно отметить, что для иерархических слоев не наблюдается четко выраженный вторичный максимум в функции распределения длин первого пролета, который был обнаружен для неоднородных слоев с полостями фиксированного размера (рис. 10 в Reshetnyk и др., 2021). Как и при расчете функции распределения длин пролетов (см. рис. 2–3), это объясняется тем, что размеры и форма полостей в модельных иерархических слоях изменчивы. Этим же объясняется заметно больший разброс значений функций для глубин, сравнимых или превосходящих размеры агрегатов, для всех рассмотренных модельных слоев. Как уже отмечалось выше, слои, построенные без контроля за контактами между агрегатами, в целом более неоднородные. Этим объясняется более сложный характер поведения функции распределения для этих слоев (правая панель).

Анализируя результаты, мы обнаружили, что выбор положения контрольной стартовой поверхности имеет значение для расчета длин первых пролетов. Этот эффект подробно обсуждался в (Reshetnyk и др., 2021). Для иерархических слоев в случае выбора нижней границы слоя в качестве контрольной поверхности наблюдается немонотонное поведение функции распределения длин первых пролетов в приграничной области. Если же контрольная поверхность размещается внутри слоя на некоторой глубине (тестировались случаи, когда эта глубина была не меньше 5 размеров мономеров), то эта немонотонность исчезает, и экспоненциальная аппроксимация удовлетворительно работает, как в случаях однородных слоев (см. рис. 9 в Reshetnyk и др., 2021). Средние значения коэффициента ослабления в экспоненциальной функции для различных слоев приведены в табл. 4. Видно, что можно говорить об увеличении глубины проникновения с увеличением пористости. Кроме того, можно говорить о большей “прозрачности” слоев, построенных из псевдомономеров. Похожий результат мы получили при нахождении аппроксимации для проницаемости слоев. Результаты представлены в табл. 2. В первом столбце таблицы указан сорт слоя, во втором – средняя пористость, в третьем – коэффициент ослабления C в размерах сфер-мономеров, D_m. Хотя средние пористости слоев (512s-r) и (BAM2.512) очень близки, однако их функции проницаемости и поглощения излучения отличаются заметно, что связано с различиями в микростроении используемых агрегатов и их организации в слое.

ЗАКЛЮЧЕНИЕ

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

В предыдущих работах нами было показано, что структурные характеристики пылевого слоя (такие как пористость, толщина, размер и строение слагающих его частиц) влияют на перенос энергии и продуктов сублимации, определяя эффективное газопроизводство и свойства газового потока над поверхностью ядра. Оценка и анализ транспортных характеристик пористого приповерхностного слоя требуются для создания новых теплофизических моделей и является актуальной задачей кометной физики. Ранее нами уже были рассмотрены различные типы модельных пористых слоев, сконструированных из непересекающихся сфер, упакованные случайным образом, с использованием различных методов построения. Были проанализированы как однородные, так и неоднородные по структуре слои, содержащие полости и микротрещины. Такие типы неоднородностей, по-видимому, характерны для приповерхностного слоя кометного ядра (Groussin и др., 2019; El-Maary и др., 2019). Рассматривая свободно-молекулярное течение и используя метод тестовых частиц, мы выполнили систематическое компьютерное моделирование переноса продуктов сублимации в слое. Полученные результаты были опубликованы в (Reshetnyk и др., 2021; Skorov и др., 2021). Главное модельное упрощение, принятое во всех предыдущих моделях, связано с использованием твердых сфер-мономеров в качестве строительных элементов, из которых создавались слои. Это допущение не согласуется с результатами анализа кометных частиц (Bentley и др., 2016; Langevin и др., 2016) и носит вынужденный характер, так как построение пористых случайных слоев из твердых сфер – это хорошо изученный и относительно несложный процесс (см., например, обзор в Skorov и др., 2011). В данной работе мы делаем важный шаг вперед в развитии моделей поверхностного кометного пылевого слоя и впервые рассматриваем иерархические слои, созданные не из единичных твердых сфер, а из пористых агрегатов. В качестве строительных единиц мы использовали баллистические агрегаты различной массы (следовательно, размера). Слои строились на основе двух различных подходов. В первом случае пористые агрегаты осаждались во внешнем гравитационном поле, и в ходе моделирования контролировались контакты между агрегатами. Во втором случае слои создавались из псевдомономеров, которые сами являлись агрегатами. Использование различных методов построения и различных размеров агрегатов позволило нам получить модельные слои со значениями средней пористости, полностью покрывающими диапазон ожидаемых значений для пористости кометного ядра, полученный из наблюдений. Таким образом, мы получили широкий набор модельных слоев, необходимых для изучения транспортных свойств поверхностного слоя и свойств, вылетающих газовых молекул.

Структура проведенного исследования совпадает со структурой ранее проведенных исследований (Reshetnyk и др., 2021). Это позволяет сравнивать полученные результаты и выделять те или иные особенности, связанные с отличиями в строении модельных слоев.

В данной работе мы:

• Разработали модели различных иерархических слоев, обладающих пористостью в диапазоне значений, полученных в наблюдениях кометы 67P/Чурюмова–Герасименко.

• Получили количественные оценки основных транспортных характеристик, таких как: распределение длин пролетов между столкновениями в слое и их среднее значение, извилистость слоя и его проницаемость для свободно-молекулярного течения.

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

• Получили оценки для объемного поглощения солнечного света в поглощающем слое.

• Используя результаты компьютерного моделирования, представили простые аппроксимационные формулы, позволяющие включить рассматриваемые эффекты в согласованные модели тепломассопереноса в кометном ядре.

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

Показано, что для всех рассмотренных слоев средняя длина пролета с достаточной точностью связана с пористостью слоя известным теоретическим соотношением (Dullien, 1991) для различных значений пористости и различных методов построения иерархических слоев. Подобный вывод был ранее сделан для слоев, сконструированных их твердых сфер и имеющих как однородную, так и неоднородную структуру. Мы делаем вывод, что эта характеристика может успешно использоваться в качестве универсальной характеристики “среднего размера пор”. Это позволяет строить общие зависимости для проницаемости слоя, включающие только два параметра: толщину слоя и средний размер пор. Анализ показал, что учет пористости в оценке проницаемости важен, и ее влияние носит нелинейный характер. Результаты, полученные в лабораторных экспериментах, не могут прямо обобщаться на слои, свойства которых отличаются от изученных в опытах образцов. В то же время параметрические зависимости могут быть построены на основе систематических компьютерных расчетов, покрывающих широкий диапазон изменения модельных параметров. Нами показано, что для ожидаемых значений средней пористости приповерхностного слоя проницаемость слоя заданной толщины может отличаться в разы. Для иерархических слоев имеет значение не только средняя пористость, но и размер агрегатов, из которых построен слой. Суммируя результаты данной работы и предшествующих, мы заключаем, что нами получены простые аппроксимационные формулы для самых различных типов пористых слоев, позволяющие количественно оценить проницаемость экранирующего слоя. Важно подчеркнуть следующее. Тот факт, что мы не знаем, как именно происходит процесс столкновения молекул с нелетучей фракцией кометного ядра, вносит в теоретические модели неустранимую неопределенность. Расчеты показали, что, например, проницаемости, рассчитанные для простейших случаев зеркального и полного диффузного рассеяния, отличаются на 20–30%. Этот разброс значений задает общее представление о допустимой точности оценок и аппроксимаций.

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

Кроме изучения структурных характеристик иерархических пористых слоев и кинетических свойств вылетевших молекул, мы продолжили изучение ослабления прямого солнечного света в приповерхностном слое. Были получены распределения интенсивности как функции глубины (расстояния от верхней границы слоя) и оценки характерной глубины ослабления излучения для полностью поглощающей среды. Как и в случае однородных слоев, последняя характеристика связана со средней длиной пролетов (а тем самым с пористостью слоя). Общий ход зависимости удовлетворительно описывается законом Бугера–Ламберта–Бера: интенсивность излучения экспоненциально убывает с глубиной. В иерархических слоях наблюдаются заметные локальные вариации пористости на масштабах, сравнимых с размерами агрегатов. Поэтому наблюдается заметный разброс значений интенсивности на близких глубинах. Для значений пористости, характерных для поверхности комет, заметная доля излучения проникает на глубину в несколько средних размеров пор. Однако для иерархических слоев средний размер “пор” существенно больше, т.е. излучение проникает в более глубокие слои ядра, а его поглощение имеет более выраженный объемный характер. Полученные результаты усиливают заключение, которое мы ранее сформулировали для слоев, сделанных из сфер-мономеров: анализ наблюдений с инструментами MIRO и VIRTIS следует проводить с учетом объемного поглощения излучения в слое.

Следует сделать еще одно важное замечание. Одной из намеченных целей нашего исследования является изучение процесса выноса пылинок с поверхности кометного ядра. До сих пор обычно предполагается, что основной механизм связан с тягой газового потока, проходящего через пористый слой. В (Skorov и др., 2017) нами были проанализированы возникающие подъемные силы, и было показано, что эти силы в абсолютном большинстве модельных случаев слишком малы, чтобы преодолеть связанность пылевой матрицы. В данной работе, как и в цитированной, мы рассматривали нейтральное незаряженное вещество. В то же время есть веские основания считать, что частицы на поверхности кометного ядра имеют электрический заряд (Popel, Gisko, 2006). Возникающие электростатические силы могут оказывать важное влияние на физику пылинок на поверхности кометы. Учет этих процессов крайне желателен при будущем моделировании пылевой активности комет.

Результаты, представленные в данной работе для различных иерархических слоев, собранных из пористых баллистических агрегатов, вместе с результатами, полученными ранее для однородных и неоднородных слоев, собранных их твердых сфер (Reshetnyk и др., 2021; Skorov и др., 2021), позволяют на новом уровне конструировать согласованные модели тепломассопереноса в приповерхностном слое кометного ядра. В этих моделях впервые непротиворечиво будут учитываться микроструктура моделируемой области, проницаемость слоя, роль радиационной теплопроводности, объемный характер поглощения солнечной энергии и т.д. Такие модели необходимы для анализа наблюдений, выполненных с инструментами MIRO и VIRTIS (см., например, Macher и др., 2019), для изучения пылевой активности (см., например, Skorov и др., 2017) и физических свойств поверхности, для лучшего понимания глобального поведения газопроизводительности с изменением гелиоцентрического расстояния (см., например, Skorov и др., 2020). Эти исследования будут нами продолжены в следующих работах.