CC BY
42
DOI: 10.18721/JPM.11109 УДК 541.11
О МОДЕЛИРОВАНИИ ПРОцЕССА ОБРАЗОВАНИЯ КЛАСТЕРОВ ПРИ СФЕРИЧЕСКОМ РАСШИРЕНИИ ВОДЯНОГО ПАРА В ВАКууМ
Н.Ю. Быков
Санкт-Петербургский политехнический университет Петра Великого, Санкт-Петербург, Российская Федерация
Методом прямого статистического моделирования выполнен расчет истечения водяного пара в вакуум с учетом процессов конденсации в объеме. Рассмотрены две модели образования кластеров, основанные на кинетическом подходе и на выводах модифицированной классической теории нуклеации. Проведен анализ физической адекватности моделей и особенностей их алгоритмической реализации. Выполнены расчеты сферического расширения пара от испаряющейся поверхности в диапазоне чисел Кнудсена, соответствующем околоконтинуальному и переходному режимам течения. Проанализировано влияние процесса конденсации на газодинамические параметры течения. Продемонстрирован эффект «замораживания» мольных концентраций кластеров по мере удаления от испаряющейся поверхности.
Ключевые слова: прямое статистическое моделирование; конденсация; кластер молекул воды; переходный режим течения
Ссылка при цитировании: Быков Н.Ю. О моделирование процесса образования кластеров при сферическом расширении водяного пара в вакуум // Научно-технические ведомости СПбГПУ. Физико-математические науки. 2018. Т. 11. № 1. С. 86 - 101. DOI: 10.18721/JPM.11109
on simulation of cluster formation process under water
VAPOR SPHERICAL EXPANSION INTO VACuuM
N.Yu. Bykov
Peter the Great St. Petersburg Polytechnic University, St. Petersburg, Russian Federation
The calculation of a water vapor outflow into vacuum has been performed using the direct simulation Monte Carlo with taking into account a condensation process. Two cluster-formation models were employed for simulation. The former is based on a kinetic approach, the latter is based on conclusions drawn from the modified classical nucleation theory. The analysis of physical adequacy of models and features of their program implementation was carried out. The simulation of spherical vapor expansion from an evaporating surface was performed over a range of Knudsen numbers corresponding to transient and near-continual flow regimes. The influence of a condensation process on gasdynamic flow parameters was analyzed. The effect of freezing of cluster mole-fractions when receding from the evaporating source was demonstrated.
Key words: direct simulation Monte Carlo; condensation; water cluster; transient flow regime
Citation: N.Yu. Bykov, On simulation of cluster formation process under water vapor spherical expansion into vacuum, St. Petersburg Polytechnical State University Journal. Physics and Mathematics. 11 (1) (2018) 86 - 101. DOI: 10.18721/JPM.11109
Введение
Процессы конденсации водяного пара в расширяющихся течениях широко распространены в природе и технике. Среди природных феноменов можно выделить атмосферные явления, в том числе образование кластеров воды в потоках газа, формирующих околоядерную атмосферу комет [1]. К техническим приложениям относятся процессы конденсации водяного пара в соплах и струях [2], в том числе при работе двигателей космических летательных аппаратов на больших высотах [3].
В зависимости от характерного числа Кнудсена (Кп = 1/Л, где Х — средняя длина свободного пробега молекул, Л — характерный размер течения) истечение пара в вакуум является разреженным на периферии струи или во всей области течения. Для расчета разреженных течений традиционно используется метод прямого статистического моделирования (ПСМ) [4].
Интерес к моделированию процессов формирования и роста кластеров в разреженных течениях методом ПСМ возник в начале 2000-х годов в связи с бурным развитием вакуумных технологий нанесения различного рода покрытий, в том числе технологий лазерной абляции материалов. Для расчета таких течений использовались кинетические модели [5 — 8], в которых вероятности ассоциации частиц при взаимных столкновениях являлись функциями размера кластера (в простейшем случае были равны единице), а мономолекулярный распад кластера моделировался формулой теории Райса — Рамспергера — Касселя (РРК). Работы отличались детальностью описания первого этапа кластерообразования — процесса димеризации и алгоритмической реализацией РРК-формулы для мономолекулярного распада. Дальнейшее развитие методики прямого статистического моделирования течений с процессами конденсации, основанной на применении кинетической модели, имело место в работе [9], где вероятности процессов ассоциации получали, применяя подходы молекулярной динамики (МД), а также в работе [10], где был предложен точный алгоритм реализации
формулы РРК мономолекулярного распада для метода ПСМ. Скорости прямого процесса ассоциации частиц и обратного процесса распада во всех указанных работах оказываются не связанными между собой. В работе [11] предложена модель конденсации, в рамках которой вероятности роста/ распада кластеров получены на основе данных о скоростях прямых и обратных реакций соответствующих процессов. Константы скоростей прямых и обратных реакций связаны через константы равновесия, которые считаются известными. Отличительными особенностями работы [11] являются оригинальные выражения для вероятности мономолекулярного распада кластеров (в частном случае они совпадают с теорией РРК) и вероятности формирования димера при трехчастичном столкновении мономеров. Указанная кинетическая модель [11] используется при получении основных результатов данной работы.
Отдельное место занимают работы [12, 13], в которых выполнен расчет ПСМ конденсации паров воды в дальнем разреженном поле струи ракетного двигателя. В данных статьях использована классическая теория нуклеации (КТН), выводы которой были адаптированы для алгоритма ПСМ.
Стремление использовать КТН для практических расчетов имеет основания. Во-первых, выводы теории довольно давно используются для описания течений с конденсацией. Расчет необходимых термодинамических параметров (давление и температура пара, степень пересыщения, скорость зародышеобразования и критический размер кластера) можно реализовать в рамках ПСМ. Во-вторых, моделирование процесса роста кластеров, начиная с момента образования критического зародыша, может дать вычислительные преимущества, так как нет необходимости в затратах вычислительного времени на моделирование процесса формирования докритических кластеров, в том числе образования димера при тройной рекомбинации мономеров.
Однако КТН имеет известные недостатки: использует понятие поверхностного натяжения по отношению к кластерам малого размера, предсказывает размер критическо-
го кластера меньший единицы при больших степенях пересыщения пара. В работах [12, 13] игнорируется эффект выделения энергии связи при ассоциации частиц, а следовательно, не учитывается эффект повышения внутренней температуры кластера, сопровождающий данный процесс; не ясна также реализация алгоритма для случая больших степеней пересыщения.
С целью устранения указанных недостатков КТН, в работе [14] предложено использовать модифицированную классическую теорию нуклеации (МКТН). В рамках модели МКТН не требуется использования понятия поверхностного натяжения [15], размер критического кластера при больших степенях пересыщения равен единице. Последнее обстоятельство требует дополнительных предположений при алгоритмической реализации модели. Модель МКТН была реализована для пространственно-однородного случая конденсации пара меди [16]. Однако прямого сравнения результатов расчетов одних и тех же течений с использованием обеих моделей (на базе кинетического подхода и подхода МКТН) выполнено не было. Кроме того, не был проведен полный анализ алгоритмических возможностей реализации модели МКТН для метода ПСМ.
Одной из целей настоящей работы является сравнение результатов расчетов (рассчитывались размерные распределения получаемых кластеров, распределения газодинамических параметров по координате) с использованием двух моделей (кинетической и МКТН) на примере одного объекта исследования.
В статье также обсуждаются особенности алгоритмической реализации подхода МКТН и параметры модели МКТН, при которых указанные расчетные результаты оказываются близкими к результатам, полученным при использовании кинетического подхода.
В качестве объекта исследования в работе выбрано одномерное стационарное течение водяного пара от испаряющейся сферической поверхности. Моделирование истечения нереагирующего газа с поверхности сферы в вакуум выполнено в рабо-
тах [17, 18]. В указанных статьях рассмотрено влияние степени разреженности на газодинамику течения, проанализированы особенности формирования кнудсеновско-го слоя у поверхности сферы и процессы релаксации поступательных степеней свободы.
В работе [19] проведен расчет параметров кластеров в области X < 5 Я (X — радиальная координата, Я — радиус сферы) и пара на звуковой линии при истечении водяного пара от сферического источника. Для расчетов использована упрощенная модель конденсации, в которой, в частности, вероятности ассоциации частиц в ходе реакций кластеризации предполагались равными единице, а процесс испарения моделировался согласно теории РРК. Для используемых параметров модели [19] наблюдалась значительная доля кластеров в течении (до 19 %), показано существенное влияние процесса конденсации на параметры звуковой линии.
Второй целью настоящей работы является проведение серии расчетов течения от сферического источника с использованием уточненной модели конденсации, предложенной в статье [11], изучение влияния процесса кластеризации на параметры течения пара и исследование параметров формирующихся кластеров.
Модель формирования кластеров
Согласно данным работы [11], кластер воды Ш характеризуется числом мономеров g, числом поверхностных мономеров g0, массой т,, диаметром й,, поступательной скоростью внутренней энергией Е ш, числом вращательных (£ ) и колебательных степеней свободы, энергией испарения одного мономера с поверхности е^
Для кластеров размера g > 2 рассматриваются следующие реакции образования/ распада:
Ш + М__1 > Шg,
(1)
ш + Шg_l,
(2)
где р , р— вероятности ассоциации и распада соответственно.
Реакция (1) описывает присоединение мономера к кластеру при парном столкновении, а реакция (2) — мономолекулярный распад кластера. Время мономолекулярного распада кластера имеет порядок 1Д0, где v0 — характерная частота колебаний мономеров в кластере. Прямой процесс (1) характеризуется константой скорости К—1, обратный процесс (2) — константой скорости К — . £
Упругое столкновение мономер-кластер происходит с вероятностью (1—р ):
Ml+ Mg
>M, + M„.
В рамках рассматриваемой модели д^и^ меры воды образуются в результате трехча-стичных столкновений:
Д(3)
~^M2 + M1,
М1 + М1 + М1
2 1 (4)
где р1(3) — соответствующая вероятность.
Обратная реакция происходит по схеме:
Plc
> M + M + M,
м2 + м ■
- 2 1 1 1 1 (5)
где р2с — вероятность столкновительного
распада.
Скорости прямой (4) и обратной (5) реакций обозначим как К1(?) и К—с соответ-
ственно.
Модель на базе модифицированной классической теории нуклеации
Отличительными чертами модели на базе МКТН являются необходимость расчета в каждой вычислительной ячейке размера критического кластера и скорости процесса зародышеобразования (нуклеации) /, которые являются функциями те-плофизических параметров.
В рамках модели МКТН скорость ну-клеации определяется следующим образом [2, 15]:
J = у I Х-^ R- = K-/ -=2
Лп1п-
(6)
где пе — концентрация --меров; Kg— — кон
станты скоростей процессов (1) и (4) для £ > 2 (в случае формирования димеров воды в ходе реакции (4) К1 = К^д); верхний индекс е означает равновесную концентра-
цию, n1e = п1.
Заданной степени пересыщения S соответствует определенный критический размер кластера -* [2, 15]. В диапазоне S < Smin критический размер оказывается равным критическому размеру, определенному согласно классической теории нуклеации (КТН) [2]; в диапазоне S . < S < S этот
\ / l j ? ^ mm max
размер равен координационному числу Nc, а для случая S > S — = 1. Величины S
^ J max max
и Smin определены в работах [14, 15]. Знание критического размера кластера необходимо для определения величины J, так как кластеры околокритического размера вносят основной вклад в сумму в выражении (6).
После определения значений J и -* вычислительный алгоритм на базе подхода МКТН сводится к вбросу кластеров в ячейку со скоростью J. Определение размера вбрасываемых кластеров является первой отличительной особенностью алгоритма, базирующегося на МКТН. В работах [12, 13] в ячейку вбрасывались кластеры критического размера -*. Вброс кластеров размера -* > 2 в ячейку приводит к «усе-ченности» реального размерного распределения кластеров. В этом случае предполагается, что концентрации докритических кластеров в ячейке соответствуют равновесному распределению. Такая ситуация требует дополнительного моделирования распределения докритических кластеров в рамках метода ПСМ (в работах [12, 13] эта проблема не рассматривалась) и изменения вычислительного алгоритма. С другой стороны, при S >> 1, в рамках подхода МКТН, -* = 1 (классическая теория ну-клеации предсказывает в этом случае, что -* < 1). Вброс мономеров в ячейку лишен смысла. Поскольку в рамках КТН скорости образования кластеров любого размера равны скорости нуклеации J [2], возможен вброс в ячейку кластеров околокритического размера или докритического размера -a [14]. Оптимальным является выбор размера кластера ga = 2, так как, с одной стороны, в любой ячейке должно воспроизводиться распределение кластеров по размерам (начиная с димера), а с другой — решается проблема вброса кластеров при S >> 1. Именно такой подход используется в на-
стоящей работе.
Второй отличительной особенностью алгоритма на базе МКТН является исключение реакций образования кластеров размера ga (при ассоциации мономера и кластера размера (^а — 1)) и их распада из вычислительного алгоритма, так как оба процесса входят в скорость зародышеобра-зования / [2].
Третьей отличительной особенностью вычислительного алгоритма на базе МКТН является необходимость определения температуры пара для розыгрыша процесса испарения. В стандартном случае для розыгрыша столкновений и передвижения молекул, в рамках метода ПСМ, определения макроскопических параметров на каждом временном шаге не требуется. Модели, базирующиеся на МКТН или КТН, являются однотемпературными, т. е. предполагают, что поступательные и внутренние температуры мономеров и кластеров равны. Такой подход требует использования одной температуры для определения констант скоростей роста и распада кластеров. Это может быть либо поступательная температура мономеров Т1, либо температура Тш — усредненная по поступательным и внутренним энергиям частиц в ячейке [11]. Но с одной стороны, при п1 >> пусредненная температура Ttot оказывается близкой к поступательной температуре мономеров Т1, и так же, как и Т1, не отражает реальной внутренней энергии индивидуального кластера g. С другой стороны, реализация алгоритма, связанная с использованием величины Тш, требует внесения некорректных (с физической точки зрения) изменений в основной столкновительный блок программы, так как число столкновений в ячейке зависит от температуры. С учетом изложенных обстоятельств, в настоящей работе для определения скоростей реакций используется температура мономеров
Тг
Следует отметить, что возможна реализация смешанного алгоритма (она также рассматривается в настоящей работе), при котором димеры вбрасываются в поле течения согласно алгоритму МКТН, а даль-
нейший процесс испарения кластеров рассматривается аналогично кинетическому подходу (вероятность испарения мономера из кластера зависит от внутренней температуры индивидуального кластера). Однако такую вычислительную модель следует признать некорректной с точки зрения основ КТН.
Кинетический подход к построению модели
В рамках кинетического подхода не требуется дополнительных предположений о размере вбрасываемых в вычислительную ячейку кластеров и вычисления скорости нуклеации. Процесс кластеризации стартует с формирования димеров в соответствии с реакцией (4).
Вероятности протекания реакций в случае известных констант скоростей определяются с использованием модели полной энергии столкновения (Total Collision Energy, TCE) [4, 20, 21]. Основой данной модели является предположение о том, что соответствующие вероятности могут быть заданы в виде
(Etot s)Y/2-1
p = L
Z ш /2-1
Jtot
(7)
где Еш — полная энергия; £ш — полное число степеней свободы; е — энергия активации; Ь, У — константы, зависящие от параметров модели упругих столкновений частиц и параметров констант скоростей реакций.
Для определения параметров Ь и У в выражении (7) необходимо задать константы скоростей реакций в аррениусовском виде:
К = АТ%хр(—е/кТ
(к — постоянная Больцмана); константы А, В и энергия активации е полагаются известными.
Следует обратить внимание на то обстоятельство, что кинетическая модель лишена недостатков модели МКТН. Она не требует определения макроскопической температуры в ячейке и ввода дополнительных предположений о размере кластеров, «вбрасываемых» в ячейку, и их параметров.
Постановка задачи
В настоящей работе ставится задача расчета сферического стационарного расширения водяного пара в вакуум с поверхности сферы радиуса R. Расчет выполнен с помощью метода прямого статистического моделирования и учитывает процесс конденсации пара в объеме.
Рассматриваемая задача является одномерной, с единственной радиальной координатой X. Однако реализация метода ПСМ требует розыгрыша всех трех компонент скорости частиц: радиальной vx и двух перпендикулярных ей компонент — v0 и уф. Температура источника T0 полагается постоянной. С поверхности испаряются только мономеры (молекулы) воды. Для описания испарения используется закон Герца — Кнудсена, в соответствии с которым поток F + испаряющихся молекул равен
F + = ПЛ/4, где n0 — концентрация насыщенного пара (соответствует равновесному давлению, определенному по T0); v0 — средняя тепловая скорость испаряющихся частиц;
vo = (8k7o/nm1)1/2
(m1 — масса мономера).
Равновесное давление пара p0(T) опре-
деляется аналогично тому, как это сделано в работе [22]. Остальные параметры водяного пара соответствуют приведенным в статье [23]. Предполагается, что функция распределения для молекул воды, улетающих с поверхности, является полумаксвел-ловской [19].
Для рассматриваемого диапазона температур источника колебательные степени свободы испаряющихся молекул воды считаются замороженными [24]. Внутренняя энергия улетающих с поверхности мономеров равна вращательной энергии
Емд = (С,д/ 2)кГ0.
Для розыгрыша энергообмена между поступательными и внутренними степенями свободы используется модель Ларсена — Боргнакке [4]. Вероятность обмена энергией между вращательными и поступательными степенями свободы молекул воды при взаимных столкновениях принята равной единице. Вероятности обмена энергией между поступательными и внутренними степенями свободы при столкновении класса мономер-кластер приняты равными 0,1 [10]. Формирование кластеров в поле течения происходит согласно вышеописанной модели конденсации. Константы скоростей реакций, а также параметры моделей столкновений взяты из работы [25].
Таблица 1
Варианты расчета сферического стационарного расширения водяного пара в вакуум
Номер База Учет Число
варианта модели конденсации Кнудсена
1 - - 10-4
2 КП + 10-4
3 - - 10-2
4 КП + 10-2
5 МКТН + 10-4
6 МКТН + КП + 10-4
(для испарения)
Примечания. 1. Температура источника одинакова для всех вариантов: 70 = 35 0 К. 2. Число Кнудсена Кп = Х0/(2Я), где \ — средняя длина свободного пробега, Я — радиус сферы.
Сокращения: МКТН — модифицированная классическая теория нуклеации, КП — кинетический подход.
Молекулы и кластеры, достигающие внешней границы области, исключались из расчета (гиперзвуковое граничное условие). Также исключались из расчета частицы, пересекающие границу источника в обратном направлении (условие полной конденсации на поверхности).
В рассмотренной постановке газодинамическая картина течения определяется числом Кнудсена
Кп = Х0 / Б = Х0 / (2Я),
где — средняя длина свободного пробега молекул воды, соответствующая параметрам п0 и Т0; Б — диаметр сферы, а также характерными скоростями реакций (1), (2), (4), (5).
Варианты расчета и значения определяющих параметров представлены в табл. 1.
Расчеты выполнены методом ПСМ со схемой столкновений без временного счетчика [4], с помощью разработанного программного комплекса на кластере «Политехник — РСК Торнадо» в Суперкомпьютерном центре «Политехнический».
Результаты расчетов и их обсуждение
На рис. 1 представлены результаты моделирования — данные об изменении нормализованных величин концентрации, скорости и температуры пара в поле течения. Плотность пара падает по мере удаления от поверхности, асимптотически стремясь к нулю; при этом скорость газа растет (рис. 1, а). В случае истечения в режиме сплошной среды (Kn = 0) скорость стремится к своему термодинамическому пределу стационарного расширения
u= (2/(у -1))1/2«0,
где у — показатель адиабаты; a0 — скорость звука, определенная по температуре T0 [26].
В случае свободномолекулярного разлета (Kn ^ ж) скорость пара стремится к тепловой скорости, определенной по равновесной температуре пара [27]:
u = к.
max 0
Теплосодержание пара падает до нуля для сплошно-среднего течения [26].
Рис. 1. Расчетные распределения плотности (кривые 1 — 4) и скорости (5 — 8) (а), а также температуры (Ь) пара мономеров. Представлена температура поступательных (9, 11, 13, 15) и внутренних (10, 12, 14, 16) степеней свободы. Расчеты проведены с учетом (1, 3, 5, 7, 9, 10, 13, 14) и без учета (2, 4, 6, 8, 11, 12, 15, 16) процессов кластеризации. Кп = 10-4 (1, 2, 5, 6, 9 — 12)
и 10-2 (3, 4, 7, 8, 13 - 16)
В случае свободномолекулярного течения температура газа падает до величины [27]:
Т / То = 1-8/(3 п)
Уменьшение степени разреженности течения (уменьшение числа Кнудсена) приводит к большему числу столкновений между частицами и более эффективному переходу тепловой энергии в энергию направленного движения частиц. В результате скорость газа для более плотного режима течения (Кп = 10-4) возрастает по отношению к более разреженному варианту (Кп = 10—2), а плотность и температура уменьшаются по мере удаления от источника.
Следует обратить внимание на результат, который показывает, что увеличение числа Кнудсена ведет к нарушению равновесия между поступательными и внутренними степенями свободы [17, 18]. Для Кп = 10-4 различие между соответствующими компонентами температуры в рассматриваемой области незначительно. Однако уже для чисел Кнудсена порядка 0,01 различие между компонентами становится существенным.
Учет конденсации в расчетах приводит к двум эффектам:
выделение энергии скрытой теплоты конденсации в поток. В рассматриваемой модели энергия связи при ассоциации частиц выделяется во внутренние степени свободы кластеров размера g > 2 (см. реакцию (1)), а также в поступательные и внутренние степени свободы мономера и димера при тройной рекомбинации мономеров в ходе реакции (4);
естественная убыль мономеров в процессе формирования и роста кластеров в течении.
Процесс конденсации, ввиду указанных выше факторов, способствует росту скорости, температуры и более быстрому падению плотности пара в течении. На рис. 1 данный эффект хорошо заметен для варианта расчета (Кп = 10-4), для которого объемная доля кластеров составляет 5,4 %. Для числа Кнудсена 0,01 объемная доля кластеров в течении оказывается менее 0,2 % и влияние процесса конденсации на газодинамику течения отсутствует.
У испаряющейся поверхности влияние процесса кластерообразования на распределения плотности и скорости даже для значения Кп = 10-4 является незначительным (см. рис. 1). Так, параметры пара на звуковой линии для вариантов расчета с учетом и без учета процесса класте-рообразования отличаются на величину менее 3 % (табл. 2). Влияние процесса конденсации на параметры звуковой линии оказалось слабее, по сравнению с результатами, полученными в работе [19]. Положение звуковой линии в вариантах с учетом и без учета конденсации составляет 48,8Х0 и 44,4Х0 от испаряющейся поверхности соответственно (для Кп = 10-4). С увеличением расстояния от поверхности влияние процесса конденсации на параметры течения (особенно на температуру) возрастает.
Особенностью рассматриваемого класса разреженных течений является его поступательная неравновесность, возрастающая с увеличением числа Кнудсена. При любом
Таблица 2
Результаты расчета параметров пара на звуковой линии
Номер варианта Учет конденсации п/п0 т/т; Доля обратного потока, %
1 — 0,316 0,793 19,3
2 + 0,307 0,803 19,5
Примечание. Представлены результаты при значении числа Кнудсена Кп = 10—4. Номера вариантов соответствуют номерам табл. 1
Обозначения: п/п0 — нормализованная концентрация пара, Т/Т0 — нормализованная температура пара.
значении числа Кнудсена области поступательной неравновесности (область, в которой компоненты Тх, Те, Тф поступательной температуры Тг различны) находятся как непосредственно у испаряющейся поверхности (кнудсеновский слой с толщиной порядка нескольких длин свободного пробега молекул), так и на некотором удалении от поверхности, когда локальный средний размер свободного пробега частиц становится большим по отношению к характерному размеру течения [18]. Данные рис. 2, а демонстрируют размеры слоя Кнудсена у поверхности. Размер слоя можно определить, исходя из необходимости одновременного выполнения следующих двух условий [18, 19]:
Т - т I
-^ • 100% < 3%;
Тг Тм\ • 100% < 3%,
где Тм — температура внутренних степеней свободы.
Полученная толщина слоя составляет 40А,0 и практически не зависит от учета процесса конденсации.
В случае свободномолекулярного разлета областью поступательной неравновесности является вся зона течения [27]. Для варианта расчета при Кп = 10-2 различие между компонентами температуры Тх и Т0 = Тф наблюдается также в большей области течения (рис. 2, Ь), а для случая Кп = 10-4 — только в дальнем поле струи. В итоге учет процесса конденсации заметного влияния на степень неравновесности течения не оказывает.
Наличие зон неравновесности (по поступательным и внутренним степеням свободы) служит одним из аргументов в пользу кинетического подхода при выборе метода расчета, так как в нем не требуется вычисления температуры пара для определения вероятностей рассматриваемых реакций.
Отдельного рассмотрения заслуживает случай расчета течения, в котором учитывается конденсация с использованием подхода МКТН. На рис. 3 представлены расчетные данные об изменении степени пересыщения, критического размера кластера и скоростей нуклеации в зависимости от координаты, полученные для расчетного варианта 5 (расчет, учитывающий процесс кластерообразования, для кото-
Рис. 2. Поступательная неравновесность у поверхности сферы (а) и в поле течения (Ь). Приведены данные для компонент температуры Т (1, 3, 5, 7); Те = Тф (2, 4, 6, 8) при значениях Кп = 10-4 (1 — 4) и 10-2 (5 — 8). Расчеты проведены с учетом (1, 2, 5, 6) и без учета (3, 4, 7, 8)
процессов кластеризации
Рис. 3. Зависимости степени пересыщения 5 (1) и критического размера кластера g, (2, 3) (а), а также скорости зародышеобразования (4, 5, 7) и Я2 (6) (Ь) от координаты для чисел Кнудсена Кп = 104 (модель МКТН). Представлены величины g, (оценка по КТН и на базе МКТН — кривые 2 и 3 соответственно), (4), 1К(5), расчетное значение / (7)
рого использована модель на базе МКТН, см. табл. 1). У поверхности источника 1,00 < Х/Я < 1,05 степень пересыщения пара мала (см. рис. 3, а), размер критического кластера g, >> 1 и совпадает для случая совместного применения КТН и МКТН. По мере удаления от источника происходит быстрый рост степени пересыщения, и, как следствие, — уменьшение критического размера кластеров. КТН предсказывает для рассматриваемого случая значения степени пересыщения менее единицы уже на расстояниях порядка одного радиуса от источника.
Подход МКТН предсказывает величину g, = Лс на расстояниях 1,05 < Х/Я < 1,90 (кривая 3 на рис. 3, а) от поверхности источника. Далее с увеличением координаты Х/Я величина g, становится равной единице.
Соответствующие рассматриваемым вариантам скорости зародышеобразова-ния представлены на рис. 3, Ь. Довольно большой размер критического кластера непосредственно у поверхности источника при малых значениях равновесной концентрации п* предопределяет малые значения скорости зародышеобразования. Далее ве-
личина скорости зародышеобразования, определенная по формуле (1) для МКТН, совпадает со скоростью процесса димери-зации в ходе тройных столкновений мономеров. Скорость зародышеобразования /. , определенная по формуле КТН [2], и
ё
Рис. 4. Распределения кластеров по размерам в поле течения, полученные по расчетным вариантам 2 (кривая 1), 5 (2) и 6 (3)
скорость с поправкой Куртни 1К = /Б существенно отличаются от значений /, полученных в рамках модели МКТН.
Интегральные по рассматриваемому объему течения распределения кластеров по размерам для варианта расчета Кп = 10-4 представлены на рис. 4. Для вариантов расчета, использующих модели, базирующиеся на МКТН и на кинетическом подходе, получены принципиально различные результаты. При кинетическом подходе частота испарения мономера из кластера определяется индивидуально для каждого кластера и зависит от значения внутренней энергии этого кластера. В ходе реакции (4) часть, а для реакции (1) вся энергия связи перераспределяются во внутренние степени свободы кластера. Для поведения кластеров, полученного расчетом при кинетическом подходе, характерно довольно интенсивное испарение, а их максимально достижимый размер в представленных расчетах ограничен четырнадцатью мономерами. Согласно расчетам на базе МКТН, кластеры испаряются с температурой, равной поступательной температуре пара мономеров. При этом скорость роста кластеров оказывается
существенно выше скорости их испарения. В итоге число кластеров с размерами g >> 1 слабо меняется с ростом величины ^
На рис. 4 также приведены данные расчетного варианта 6 (см. табл. 1), для которого вброс критических кластеров в ячейку был реализован на базе модели МКТН, а процесс испарения каждого кластера зависел от его индивидуальной внутренней температуры (как для кинетического подхода). Видно, что использование такого подхода существенно сближает результаты с таковыми, полученными на основе кинетической модели.
Распределение плотности и мольной доли кластеров по координате приведено на рис. 5 для числа Кнудсена 10-4. Диме-ры (^ = 2) формируются в непосредственной близости от поверхности источника. Максимумы на распределениях пятимеров (^ = 5) смещены от поверхности. Расстояние от максимума на распределении до источника увеличивается с ростом размера кластера. На рис. 5, Ь приведены данные о распределениях по координате мольной доли кластеров размером g = 2, 3 и 5. Представленные результаты свидетельствуют о
Рис. 5. Распределение концентрации (а) и мольной доли (Ь) кластеров по радиальной координате для кластеров разного размера g. Расчеты проведены на базе кинетического подхода (1, 2, 7), на базе МКТН (3, 4) и согласно варианту 6 (5, 6); g = 2 (1, 3, 5), 3 (7) и 5 (2, 4, 6).
Число Кнудсена Кп = 10-4
замораживании мольных долей, начиная с некоторой координаты, зависящей от размера кластера. В силу малости скорости нуклеации / вблизи источника, концентрации димеров у поверхности, получаемые в рамках прямого статистического моделирования на базе МКТН, малы по отношению к получаемым концентрациям на базе кинетического подхода, а максимум распределения димеров оказывается смещенным от поверхности на величину 0,25Л.
Концентрации пятимеров в большей части течения оказываются существенно меньшими аналогичных концентраций, полученных при использовании модели, основанной на кинетическом подходе. В случае варианта 6 (модель МКТН, а испарение рассматривается в соответствии с кинетическим подходом) распределения концентраций димеров и пятимеров оказываются ближе к результатам расчета на основе кинетического подхода (вариант 5), однако совпадения данных не наблюдается.
Следует отметить, что результаты, полученные с использованием модели [11], отличаются от таковых работы [19]. Размерные
Рис. 6. Распределения внутренней температуры частиц по радиальной координате для кластеров двух размеров: димеров (кривые 1, 3) и тримеров (2, 4). Расчеты проведены на базе кинетического подхода (1, 2) и на базе МКТН (3, 4)
распределения кластеров имеют максимум, соответствующий размеру димера.
Поступательные скорости формирующихся кластеров для варианта расчета при Кп = 10-4 совпадают со скоростью мономеров. Внутренние температуры кластеров представлены на рис. 6. Эти температуры димеров и тримеров оказываются выше поступательной температуры за счет выделения скрытой теплоты конденсации в поток в ходе процесса роста кластеров. Процесс роста сопровождается процессом испарения кластеров, в ходе которого энергия связи поглощается. В рамках модели на основе кинетического подхода скорость испарения каждого кластера зависит от его собственной внутренней температуры и указанная скорость оказывается сопоставимой со скоростью роста. В результате процесс испарения приводит к умеренному различию поступательных и внутренних температур кластеров и ограниченному максимальному размеру кластера (см. рис. 4). Для модели МКТН скорость испарения определяется поступательной температурой мономеров. Определяющим процессом при использовании модели МКТН становится процесс роста. Скорость процесса испарения кластеров, уменьшающая их внутреннюю энергию, мала по отношению к скорости их роста. В итоге наблюдается рост кластеров больших размеров ^ >> 1), а энергия, которая выделяется в ходе этого процесса, накапливается во внутренних степенях свободы. Внутренняя температура кластеров, полученных в расчетах с использованием модели МКТН, оказывается большей, чем внутренние температуры кластеров, которые предсказывает кинетическая модель конденсации (см. рис. 6).
Заключение
В работе методом прямого статистического моделирования выполнен расчет истечения водяного пара в вакуум с поверхности сферического источника с учетом образования кластеров в поле течения. Для расчетов использованы две модели процесса кластерообразования:
основанная на выводах модифицированной классической теории нуклеации;
базирующаяся на кинетическом подходе.
Одним из основных недостатков модели МКТН является ее однотемпературность. В рамках данного подхода не удается правильно описать процессы испарения мономеров из кластеров. Реализация модели в вычислительном алгоритме является сложной и требует введения дополнительных предположений. Результаты расчетов основных характеристик процесса кластерообразования (пространственное распределение концентрации кластеров, распределение кластеров по размерам в поле течения), полученных с использованием модели МКТН, существенно отличаются от результатов, полученных с помощью модели на основе кинетического подхода.
Кинетическая модель свободна от всех недостатков модели МКТН. В рамках указанного подхода возможно получение вероятностей роста/испарения кластеров как функций индивидуальных параметров сталкивающихся/распадающихся частиц, что полностью соответствует основам метода ПСМ. Для моделирования процесса конденсации не требуется определять температуру пара в ячейке. С методической точки зрения, использование модели, основанной на кинетическом подходе, является оправданным.
Расчеты процесса конденсации водяного пара при истечении со сферической поверхности в вакуум на основе кинетической модели в околоконтинуальных и переходных по числу Кнудсена режимах показывают отсутствие влияния процесса конденсации на газодинамические параметры течения при степенях кластеризации
менее 1 %. Для околоконтинуальных режимов (Кп = 10-4) степень кластеризации становится ощутимой (5,4 %). В этом случае процесс выделения скрытой теплоты конденсации в поток приводит к более быстрому росту скорости и температуры газа, падению его плотности. Особенно сильное влияние процесс кластерообразования оказывает на распределения поступательных и внутренних температур пара. Однако даже для варианта расчета, соответствующего Кп = 10-4, процесс конденсации не оказывает существенного влияния на параметры звуковой линии, размер кнудсеновского слоя и степень поступательной неравновесности по отношению к соответствующим результатам расчета, в котором не учитывается процесс конденсации. Максимальная концентрация кластеров в поле течения соответствует размеру димера.
В работе показано наличие эффекта «замораживания» мольных долей кластеров по мере удаления от поверхности источника. Для кластеров большего размера эффект наблюдается на больших расстояниях от испаряемой поверхности. Скорости и поступательные температуры кластеров, формирующихся в течении, оказываются близкими к соответствующим параметрам мономеров. Выделение энергии связи в процессе ассоциации частиц приводит к большим значениям внутренних температур кластеров по отношению к внутренней температуре мономеров.
Работа поддержана Проектом
№ 16.8549.2017/8.9 Министерства образования и науки России.
СПИСОК ЛИТЕРАТУРЫ
1. Crifo J.F. Water clusters in the coma of comet Halley and their effect on the gas density, temperature, and velocity // Icarus. 1990. Vol. 84. No. 2. Pp. 414-446.
2. Стернин Л.Е. Основы газодинамики двухфазных течений в соплах. М.: Машиностроение, 1974. 212 с.
3. Платов Ю.В., Семенов А.И., Филиппов Б.П.
Конденсация продуктов сгорания в выхлопной струе ракетных двигателей в верхней атмосфере // Геомагнитизм и аэрономия. 2011. Т. 51. № 4. C. 556-562.
4. Bird G.A. Molecular gas dynamics and the
direct simulation of gas flows. Oxford: Clarenton Press, 1994. 458 р.
5. Briehe В., urbassek H.M. Monte Carlo simulation of growth and decay processes in a cluster aggregation source // J. Vac. Sci. Technol. A. 1999. Vol. 17. No. 1. Pp. 256-265.
6. Zeifman M.I., Garrison B.J., Zhigilei L.v. Combined molecular dynamics direct simulation Monte Carlo computation study of laser ablation plume evolution // J. Appl. Phys. 2002. Vol. 92. No. 4. Pp. 2181-2193.
7. Быков Н.Ю., Лукьянов Г.А. Прямое статистическое моделирование импульсной лазер-
ной абляции металлов с процессами кластеризации в испаренном облаке // Теплофизика и аэромеханика. 2006. Т. 13. № 4. С. 569-582.
8. Itina T.E., Gouriet K., Zhigilei L.v., Noel S., Hermann J., Sentis M. Mechanisms of small clusters production by short and ultra-short pulse laser ablation // Appl. Surf. Sci. 2007. Vol. 253. No. 19. Pp. 7656-7661.
9. Li Z., Zhong J., Levin D.A., Garrison B.J. Development of homogeneous water condensation models using molecular dynamics // AIAA Journal. 2009. Vol. 47. No. 5. Pp. 1241-1251.
10. Jansen R., Wysong I., Gimelshein S., Zeifman M., buck u. Nonequilibrium numerical model of homogeneous condensation in argon and water vapor expansions // J. Chem. Phys. 2010. Vol. 132. No. 24. P. 244105.
11. bykov N.Yu., Gorbachev Yu.E. Mathematical models of water nucleation process for the Direct Simulation Monte Carlo method //Applied Mathematics and Computation. 2017. Vol. 296. No. (1 March). Pp. 215-232.
12. Zhong J., Zeifman M.I., Gimelshein S.F., Levin D.A. Modeling of homogeneous condensation in supersonic plumes with the DSMC method // AIAA Paper. 2004. 42nd Aerospace Sciences Meeting. Jan. 5 - 8, Reno, Nevada. 2004-0166.
13. Zhong J., Zeifman M.I., Gimelshein S.F., Levin D.A. Direct simulation Monte Carlo modeling of homogeneous condensation in supersonic plumes //AIAA Journal. 2005. Vol. 43. No. 8. Pp. 1784-1796.
14. Быков Н.Ю., Горбачев Ю.Е. Прямое статистическое моделирование процессов формирования кластеров в газовой фазе: классический подход с поправкой на размер кластера // Теплофизика высоких температур. 2015. Т. 53. № 2. C. 291-300.
15. Zhukhovitskii D.I. Size-corrected theory of homogeneous nucleation // J. Chem. Phys. 1994. Vol. 101. No. 6. Pp. 5076-5080.
16. bykov N.Y., Gorbachev Y.E. Comparative analysis of condensation models within DSMC // AIP Conference Proceedings. 2014. Vol. 1628. No. 1. Pp. 139-147.
17. Булгакова Н.М., Плотников М.Ю., Ре-бров А.К. Моделирование стационарного расширения газа с поверхности сферы в вакуум // Механика жидкости и газа. 1997. № 6. C. 137-143.
18. Лукьянов Г.А., Ханларов Гр.о. Стационарное расширение паров воды с поверхности сферы в вакуум // Теплофизика и аэромеханика. 2000. Т. 7. № 4. C. 511-521.
19. Быков Н.Ю. Моделирование процесса конденсации при сферическом расширении водяного пара в вакуум // Теплофизика и аэромеханика. 2009. Т. 16. № 2. С. 189-199.
20. bird G.A. Simulation of multi-dimensional and chemically reacting flows // Rarefied Gas Dynamics. Ed. By R. Campargue. Commissariat a Lenergie Atomique. Paris. Vol. 1. 1979. Pp. 365-388.
21. Gimelshein S.F., Gimelshein N.E., Levin D.A., Ivanov M.S., Wysong I.J. On the use of chemical reaction rates with discrete internal energies in the direct simulation Monte Carlo method // Physics of Fluids. 2004. Vol. 16. No. 7. Pp. 2442-2451.
22. Goldman N., Fellers R.S., Leforestier C., Saykally R.J. Water dimers in the atmosphere: equilibrium constant for water dimerization from the VRT (ASP-W) potential surface // J. Phys. Chem. A. 2001. Vol. 105. No. 3. Pp. 515-519.
23. Wolk J., Strey R. Homogeneous nucleation of H2O and D2O in comparison: the isotope effect // J. Phys. Chem. B. 2001. Vol. 105. No. 47. Pp. 11683-11701.
24. Xu X., Goddard W. Bonding properties of water dimer: a comparative study of density functional theories // J. Phys. Chem. A. 2004. Vol. 108. No. 12. Pp. 2305-2313.
25. bykov N.Yu., Gorbachev Yu.E., Zakharov v.v. Simulation of small cluster formation in water vapor plumes // AIP Conference Proceedings. 2016. Vol. 1786. No. 1. Pp. 050001-1-050001-8.
26. Черный Г.Г. Газовая динамика. М.: Наука, 1988. 424 с.
27. Кошмаров Ю.А., Рыжов Ю.А. Прикладная динамика разреженного газа. М.: Машиностроение, 1977. 184 с.
Статья поступила в редакцию 30.11.2017, принята к публикации 08.02.2018.
СВЕДЕНИЯ ОБ АВТОРЕ
БыКоВ Николай Юрьевич — кандидат физико-математических наук, начальник научно-исследовательского отдела Санкт-Петербургского политехнического университета Петра Великого, Санкт-Петербург, Российская Федерация.
195251, Российская Федерация, г. Санкт-Петербург, Политехническая ул., 29 nbykov2006@yandex.ru
REFERENCES
[1] J.F. Crifo, Water clusters in the coma of comet Halley and their effect on the gas density, temperature, and velocity, Icarus. 84 (2) (1990) 414-446.
[2] L.E. Sternin, Osnovy gazodinamiki dvukhfaznykh techeniy v soplakh [Fundamentals of gas dynamics of two-phase flows in the nozzles]. Mashinostroyeniye, Moscow, 1974.
[3] Yu.v. Platov, B.P. Filippov, A.I. Semenov, Condensation of combustion products in the exhaust plumes of rocket engines in the upper atmosphere, Geomagnitism and Aeronomy. 51 (4) (2011) 550-556.
[4] G.A. Bird, Molecular gas dynamics and the direct simulation of gas flows, Clarenton Press, Oxford, 1994.
[5] B. briehe, H.M. urbassek, Monte Carlo simulation of growth and decay processes in a cluster aggregation source, J. Vac. Sci. Technol. A. 17 (1) (1999) 256-265.
[6] M.I. Zeifman, B.J. Garrison, L.v. Zhigilei, Combined molecular dynamics direct simulation Monte Carlo computation study of laser ablation plume evolution, J. Appl. Phys. 92 (4) (2002) 2181-2193.
[7] N.Yu. bykov, G.A. Lukyanov, Direct simulation Monte Carlo of pulsed laser ablation of metals with clusterization processes in vapor, Thermophysics and Aeromechanics. 13 (4) (2006) 569-582.
[8] T.E. Itina, K. Gouriet, L.v. Zhigilei, et al.,
Mechanisms of small clusters production by short and ultra-short pulse laser ablation, Appl. Surf. Sci. 253 (19) (2007) 7656-7661.
[9] Z. Li, J. Zhong, D.A. Levin, B.J. Garrison, Development of homogeneous water condensation models using molecular dynamics, AIAA Journal. 47 (5) (2009) 1241-1251.
[10] R. Jansen, I. Wysong, S. Gimelshein, et al., Nonequilibrium numerical model of homogeneous condensation in argon and water vapor expansions, J. Chem. Phys. 132 (24) (2010) 244105.
[11] N.Yu. bykov, Yu.E. Gorbachev, Mathematical models of water nucleation process for the Direct Simulation Monte Carlo method, Applied Mathematics and Computation. 296 (1 March) (2017) 215-232.
[12] J. Zhong, M.I. Zeifman, S.F. Gimelshein, D.A. Levin, Modeling of homogeneous condensation in supersonic plumes with the DSMC method, AIAA Paper (2004), 42nd Aerospace Sciences Meeting, Jan. 5 - 8, Reno, Nevada, 2004-0166.
[13] J. Zhong, M.I. Zeifman, S.F. Gimelshein, D.A. Levin, Direct simulation Monte Carlo modeling of homogeneous condensation in supersonic plumes,
Received 30.11.2017, accepted 08.02.2018.
AIAA Journal. 43 (8) (2005) 1784-1796.
[14] N.Y. Bykov, Y.E. Gorbachev, Direct statistical simulation of the processes of clusters formation in the gas phase: classical approach with cluster size correction, High Temperature. 53 (2) (2015) 291-300.
[15] D.I. Zhukhovitskii, Size-corrected theory of homogeneous nucleation, J. Chem. Phys. 101 (6) (1994) 5076-5080.
[16] N.Y. bykov, Y.E. Gorbachev, Comparative analysis of condensation models within DSMC, AIP Conference Proceedings. 1628 (1) (2014) 139-147.
[17] N.M. bulgakova, M.Yu. Plotnikov, A.K. Rebrov, Modeling steady-state gas expansion from a spherical surface into a vacuum, Fluid Dynamics. 32 (6) (1997) 137-143.
[18] G.A. Lukyanov, Gr.o. Khanlarov, Stationary water vapor expansion from a spherical surface into vacuum, Thermophysics and Aeromechanics. 7 (4) (2000) 511-521.
[19] N.Yu. bykov, Modelling of condensation at spherical expansion of water vapor into vacuum, Thermophysics and Aeromechanics. 2009. 16 (2) (2009) 189 -199.
[20] G.A. bird, Simulation of multi-dimensional and chemically reacting flows, In: Rarefied Gas Dynamics, ed. by R. Campargue, Commissariat a Lenergie Atomique, Paris, 1 (1979) 365-388.
[21] S.F. Gimelshein, N.E. Gimelshein, D.A. Levin, et al., On the use of chemical reaction rates with discrete internal energies in the direct simulation Monte Carlo method, Physics of Fluids. 16 (7) (2004) 2442-2451.
[22] N. Goldman, R.S. Fellers, C. Leforestier, R.J. Saykally, Water dimers in the atmosphere: Equilibrium constant for water dimerization from the VRT (ASP-W) potential surface, J. Phys. Chem. A. 105 (3) (2001) 515-519.
[23] J. Wolk, R. Strey, Homogeneous nucleation of H2O and D2O in comparison: the isotope effect, J. Phys. Chem. B. 105 (47) (2005) 11683-11701.
[24] X. Xu, W. Goddard, Bonding properties of water dimer: a comparative study of density functional theories, J. Phys. Chem. A. 108 (12) (2004) 2305-2313.
[25] N.Y. bykov, Y.E. Gorbachev, v.v. Zakharov, Simulation of small cluster formation in water vapor plumes, AIP Conference Proceedings, 1786 (1) (2016) 050001-1-050001-8.
[26] G.G. Cherniy, Gazovaya dinamika [Gas dynamics], Nauka, Moscow, 1988.
[27] Yu.A. Koshmarov, Yu.A. Ryzhov, Prikladnaya dinamika razrezhennogo gaza [Applied rarefied gas dynamics], Mashinostroyeniye, Moscow, 1977.
THE AUTHOR
BYKov Nikolay Yu.
Peter the Great St. Petersburg Polytechnic University
29 Politechnicheskaya St., St. Petersburg, 195251, Russian Federation
nbykov2006@yandex.ru
© Санкт-Петербургский политехнический университет Петра Великого, 2018
