УДК 004.021
С.С. Равдин1'2, А.В. Пругло1'2, В.Н. Лагуткин1'2, Ф.А. Козин1'2
1 Московский физико-технический институт (Национальный исследовательский университет) 2 Межгосударственная акционерная корпорация «Вымпел»
Моделирование пространственных распределений параметров сверхзвукового
потока газа с частицами
Разработана компьютерная модель для расчёта пространственных распределений параметров сверхзвуковой струи газа с частицами, истекающей в спутный воздушный поток. Модель основана на численном решении методом сеток в сочетании с методом расщепления системы уравнений в частных производных, описывающих течение гетерогенного потока. В модели учтены все основные физические факторы, влияющие на формирование потока. Проведено моделирование течения сверхзвукового потока газа с частицами при умеренно нерасчётных и сильно нерасчётных режимах. Полученные результаты расчётов соответствуют общим физическим закономерностям изменения структуры и параметров течения вдоль струи в зависимости от характеристик спутного потока.
Ключевые слова: моделирование, сверхзвуковая струя газа с частицами, гетерогенный поток, спутный поток, турбулентность.
I. Введение
В настоящее время благодаря исследованиям, опубликованным в ряде монографий и статей, общие закономерности изменения структуры и параметров сверхзвуковых струй в зависимости от давления, плотности и скорости набегающего потока достаточно хорошо известны [1-5]. Эти закономерности определяются главным образом следующими физическими факторами:
— сверхзвуковая скорость истечения продуктов сгорания из сопла,
— неизобаричность струи, вызванная тем, что давление на срезе сопла отличается от давления окружающей атмосферы,
— многофазность струи, содержащей газ и частицы,
— многокомпонентность газа струи,
— распределение частиц по фракциям (размерам),
— неравновесность газа по степеням свободы и частиц по фракциям,
— смешение продуктов сгорания со спутным воздушным потоком, параметры которого зависят от высоты и скорости полета,
— сильная вытянутость вдоль оси.
В работах [1-4, 6, 7] приведены описания различных методов и алгоритмов численного моделирования сверхзвуковых струй, включая реагирующие и гетерогенные струи в атмосфере. В зависимости от назначения моделей эти методы и алгоритмы по степени детализации газодинамической структуры струи условно разделяются на методы сквозного счета и методы с выделением разрывов, а по свойствам используемых моделей физико-химических процессов — на равновесные и неравновесные. Алгоритмы различаются по способу построения вычислительной схемы, реализу-
ющей переход от дифференциальных уравнений к соответствующим разностным аппроксимациям. Для получения вычислительных алгоритмов используют метод интегральных соотношений, явные и / или неявные схемы.
Очевидно, что увеличивающая сложность математических моделей отражает рост возможностей вычислительной техники. Однако применительно к проблеме моделирования течения реагирующей неравновесной гетерогенной струи в спут-ном воздушном потоке с учётом всей совокупности физико-химических и кинетических процессов ограничения вычислительной техники в настоящее время остаются весьма существенными, что приводит к необходимости при разработке компьютерной модели с высокими функциональными возможностями выбирать соответствующий уровень сложности и детализации составляющих частных моделей.
Кроме того, для проведения корректных расчётов пространственных распределений характеристик выхлопной струи необходимо задать исходные данные по большому списку параметров. В прикладных исследованиях точное определение исходных данных составляет серьезную проблему, при этом типичным является случай, когда имеется значительная неопределённость по ряду входных параметров и констант. В этом случае нет необходимости в использовании сложных детальных моделей перечисленных выше явлений и можно воспользоваться более простыми методами расчёта течения выхлопной струи, которые тем не менее учитывают все вышеперечисленные физические факторы и при надлежащем выборе параметров моделей физико-химических процессов обеспечивают хорошую методическую точность.
В данной работе описана компьютерная программа (модель) расчёта параметров сверхзвуко-
вых струй, разработанная на основе упрощённых моделей физико-химических процессов, и приведены результаты расчётов для различных входных данных.
II. Входные и выходные данные, основные блоки и банк данных модели
Входными данными модели являются следующие данные.
1. Параметры спутного воздушного потока:
давление, плотность, температура и скорость
потока.
2. Параметры струи в выходном сечении сопла: радиус выходного сечения сопла, давление, температура газа и частиц, скорость газа и частиц, состав продуктов сгорания, размеры частиц.
Для параметров струи в выходном сечении сопла используются известные из опубликованных материалов значения и учитываются следующие основные компоненты продуктов сгорания:
— газовые компоненты: O, O2, H, H2, HO, H2O, CO, CO2, N2, NO, HCl, Cl, CI2;
— частицы: AI2O3 (окись алюминия), C (сажа).
Выходными данными модели являются распределения в поперечных сечениях струи следующих параметров:
1) давление, плотность, температура, относительные массовые концентрации и колебательные температуры компонент газовой фазы,
2) плотность и температура частиц каждой фракции.
Основные блоки предлагаемой компьютерной модели показаны на рис. 1. Расчёт основан на численном решении методом сеток системы уравнений в частных производных, описывающих течение гетерогенной струи (газ + частицы). Задачи, решаемые основными блоками модели, определены в табл. 1.
Рис. 1. Схема модели
Таблица 1
Задачи, решаемые основными блоками модели
№ Наименование блока Решаемая задача
1 Управляющий блок Управление работой модели
2 Расчёт невязкого течения Рассчитываются параметры струи и спутного потока в текущем поперечном слое без учёта эффектов вязкости, релаксации, влияния частиц.
3 Расчёт поперечного переноса в зоне смешения Рассчитываются изменения параметров потока в зоне смешения, вызванные эффектами турбулентного и молекулярного переноса.
4 Расчёт физико-химических превращений Рассчитываются изменения параметров потока в зоне смешения, вызванные химическими реакциями.
5 Расчёт колебательной релаксации Рассчитываются изменения параметров потока, обусловленные процессами колебательной релаксации.
6 Расчёт параметров потока частиц Рассчитываются параметры потока частиц с учётом кристаллизации.
7 Расчёт взаимодействия потоков газа и частиц Рассчитываются изменения параметров потоков газа и частиц, обусловленные процессами динамического и теплового взаимодействия газа и частиц.
8 Расчёт границы струи Рассчитывается граница контактного разрыва между струей и спутным потоком.
9 Расчёт параметров турбулентности Рассчитываются параметры, определяющие зависимость турбулентной вязкости от продольной и поперечной координат.
10 Расчёт структуры течения с учётом эффектов разреженного газа Рассчитываются поперечные профили ударных волн и контактного разрыва в случае разреженного спутного потока.
И Блок формирования и отображения результатов расчётов Формируются массивы распределений параметров струи в поперечных сечениях, отображаются поперечные и продольные профили основных параметров струи.
Банк данных модели включает:
— параметры моделей турбулентной вязкости;
— параметры моделей молекулярной вязкости, теплопроводности, диффузии в газах;
— параметры моделей колебательной релаксации молекул;
— параметры моделей динамического и теплового взаимодействия частиц и газа;
— параметры температурных зависимостей констант скоростей химических реакций;
— значения теплоты образования и теплоёмкости молекул.
III. Система уравнений и соотношений для турбулентного течения сверхзвуковой струи газа с частицами в спутном воздушном потоке
Сверхзвуковая струя в спутном воздушном потоке представляет собой гетерогенный поток с взаимопроникающим движением составляющих потока и обменом массой, импульсом и энергией. Составляющими потока являются газовые компоненты, соответствующие различным молекулам с учётом квантовых состояний молекул, и частицы различных размеров.
Решение проблемы определения параметров струи, истекающей в спутный поток воздуха, основано на совместном моделировании газодинамических и физико-химических процессов в гетерогенном потоке. Эти процессы описываются законами и уравнениями термодинамики, газодинамики, динамики (для частиц), химической и колебательной кинетики.
Пусть поток состоит из K газовых компонент и M групп частиц различных размеров (K + M = N). Рассмотрим стационарное осе-симметричное течение: производные по времени равны нулю, x, y — продольная и поперечная координата, u, v — продольная и поперечная скорости. Система дифференциальных уравнений для осреднённых параметров турбулентного потока (знак усреднения опущен) при неравновесном протекании химических реакций и неравновесном возбуждении колебательных степеней свободы в параболическом приближении имеет следующий вид [1-4].
Уравнение сохранения массы газа:
дри dpv pv dx dy y
(1)
Уравнения сохранения для компонент вектора импульса газа:
д (ри2 + р) дрпу рпу дх ду у
1д_
У ду
du dy
м
-53 psjf
sx],
(2)
dpuv d(pv2 + p) pv2 dx dy y
4 д
dv
Зуду \У^Нду
м
— 53 psj fsyj ■ j=1
(3)
Уравнение сохранения энергии газа: dpHu dpHv pHv
dx
dy
y
1д
y dy \ Pr
У л
Vt
dH . du
— + Pr-l}u — dy t dy
м
-pQ - 53 psj qsj ■
j=i
(4)
Уравнения сохранения для массовых концентраций газовых составляющих С.
дрС.и дрС\у рС.V 1 д ( рь дСЛ —-------- —' ' >■'' ■
dy y y dy Smt dy
K
i = 1, ■■■, K, Ci = 1.
k=i
(5)
Уравнение для среднего числа колебательных квантов молекул а.:
дрсци драцу pa.iV
dx
1 d f pt da.
у ду у Smt ду
dy y
+ pEvi - pQvi, i = 1, ■■■, K■
(6)
Уравнение сохранения массы частиц:
3psj usj dpsj vsj psj vsj
dx
■ +
+
ду ' у у ду \У Smt3j ду J ' j = 1, ■■■, M. (7)
Уравнения сохранения для компонент вектора им-
пульса частиц:
2
dPsj'Ksj dpsjUsjVsj PsjU.sjV.sj
dx
dy
1 d dusj
-TT Wt-zr
y dy dy
+ psj fsxj ,
(8)
3psj usj vsj . dpsj vnj psj
dx
-+
dy
»+LiL*=Pijf. j I......M.
(9)
Уравнение сохранения энергии частиц:
dpsj^sj'Usj | д p sjt sjV sj ^ psj€ sjVsj
dx
1 d f pt de
уду V PrtSj dy
dy y
^ +psj qsj -psj Qsj, j = 1,■■■,M■
(10)
Здесь использованы следующие обозначения: i
p = i=i pi — плотность газа,
p = i=i pi — давление газа,
y
2
РаяТ^М=1 =2^1 РаЗ !вз — суммарная сила воздействия частиц на газ,
раЗ = Раз 1раз — сила воздействия газа на частицы 3-й фракции, отнесенная к плотности частиц 3-й фракции раз,
Н = 70, £ + е., + /?.п + Щг + Щ--полная энталь-
7о — 1 р и 0 2 2
пия,
7о — показатель адиабаты при условии, что колебательные степени свободы молекул заморожены,
ev
= ¿=i evi — колебательная энергия газа,
ho = J2¿=1 hoi — теплота образования газа, Cj = — массовая концентрация i-й компоненты газа,
pt — эффективный динамический коэффициент вязкости (сумма молекулярного и турбулентного коэффициентов вязкости),
Prt = cppt/Xt — число Прандтля, Smt = pt /pDt — число Шмидта (предполагаем, что число Льюиса Let = Prt/Smt = 1 [1]),
Xt,Dt — коэффициенты теплопроводности и диффузии,
Prtpj, Smtpj — эффективные числа Прандтля и Шмидта для частиц (определяются аналогично), esj = esj (Tsj) — внутренняя энергия частиц j-й фракции,
qsj — поток теплоты между газом и частицами j-й фракции,
Ri — скорости химических превращений, Evi — скорости колебательной релаксации, Qvi — радиационные потери колебательной энергии г-го компонента,
Q,Qsi — радиационные потери энергии газа и частиц.
Система уравнений (1)— (10) дополняется термическим уравнением состояния, связывающим давление, плотность и поступательную температуру газа:
Р = P(Ro/P)T, к
P = J2 Pi, i=1
1 P
к к
i=1
£i
Pi'
к
Р = Pi i=1
калорическими уравнениями состояния, определяющими внутреннюю энергию газовых компонент и газа в целом:
к
е = ^2 Сгвг, вг = е0 (Т) + еуг (Туг ),
0/т\ 0 гт I и о 1г Во
е,: (Т) = СУ{Г + /?,>;„ с-ы = - —
2 Цг
(с°г — теплоёмкость поступательных и вращательных степеней свободы при постоянном объёме,
— число поступательных и вращательных степеней свободы молекул, ког — энтальпия образования молекул), и уравнениями состояния для частиц:
еаз = еаЗ (ТаЗ ), 3 = 1,...,М.
Систему (1)— (10) можно записать в векторном виде:
С2А + 2Ё. + £ -(У(1 Ё (и)
дх ду у у ду у г ду !
где
А = [ри,ри2 + р,риу,рНи,рСги,раи,раЗ иаз ,раз у?аз,
]т
раЗ иаЗ УаЗ ,раЗ еаЗ иаЗ] ?
В = [ру,риу,(ру2 + р),рНу,рСгУ,рагУ,раЗ УаЗ, 2 "Г
раЗ иаЗ УаЗ ,раЗ У аЗ ,раЗ еаЗ УаЗ \ ,
C = [pv,puv,pv2 ,pHv ,pCiV,paiV,psj Vsj ,Psj u
sj sj
2 ] -psj vsj ,psj esj vsj J
Psj
n ГП 4 H , n Ci ai
D= --1-(1 — Pi") ———
3 Prt t 2 Smt Smt Smtsj
sj Prtsj
M
M
Е = [0, раЗ 1ахЗ, раЗ 1ауЗ ,
З=1 З=1
М
РЯ РаЗ ЧаЗ ,рЯъ,рЕьг - рЯуг,0,
З=1
раЗ ^ахЗраЗ !'ауЗраЗ (ЧаЗ — ЯаЗ ^ , г = 1,...,К, з = 1,...,М.
Для системы (11) можно записать эквивалентную систему интегральных соотношений, учитывающих к тому же соотношения на разрывах. Пусть решение системы (11) необходимо найти в криволинейной полосе:
х > 0, Ушш (х) < У < Утах(х).
Разобьем полосу на N элементарных равномерных полосок переменной (по х) толщины 1гу(х) = _ Верхней и нижней грани-
цей п-й полоски являются соответственно кривые
Уп(х) = Ут1п(х) +пку (х), Уп—1(х) = Ут1п(х) + (п - !)Ну (х), П =1, ..., N. Интегрируя уравнения (11) по у в пределах полосок, получим N векторных интегральных соотношений
д_
дх
Ady = y'xA - B + pt
dD ду
+
Уп-1
у дУ у
= 1, ..., N.
У
У
Уп- 1
У
Уп— 1
Система интегральных соотношений (12) используется для построения разностной схемы сквозного расчёта параметров струи.
Отметим, что система (12) имеет вид
дР
Ро
(13)
где оператор Ь представляет собой композицию операторов, описывающих воздействия соответствующих физических факторов (невязкой газодинамики, диффузии, физико-химических превращений, колебательной релаксации, взаимодействия газа и частиц):
Ь = £
8=1
(14)
Для численного решения такой задачи применяется метод расщепления по факторам, суть которого заключается в следующем [6]. Дискретная аппроксимация системы (13) первого порядка точности по Дх имеет вид
р¿+1 _ р. Ах
Ь(Р1), Р0 = Р0.
С учётом (14)
Б
р.+1 = (1 + Дх^Ь3)(р ¿), 8=1
где I — единичный (тождественный) оператор.
Эту разностную систему с точностью до первого порядка по Дх можно преобразовать к виду, содержащему композицию новых операторов, зависящих от Ь8 порознь:
Лш = Ц(! + ДхЬ8)().
8=1
В результате алгоритм вычисления р¿+1 расщепляется на последовательность этапов:
<38 = (1 + ДхЬ3)С3-1, 8 = 1,...,Б,
30 = р\ р.+1 = бз.
Таким образом, расчёт параметров газовой фазы струи на одном шаге по х разделяется на следующие этапы:
1) расчёт невязкой газодинамики;
2) расчёт процессов вязкого смешения выхлопной струи и спутного потока;
3) расчёт химической кинетики в зоне смешения;
4) расчёт динамического и теплового воздействия потока частиц на газ;
5) расчёт колебательной релаксации и потерь на излучение.
Расчёт параметров потока частиц на одном шаге по х разделяется на следующие этапы:
1) расчёт переноса массы и энергии потока частиц;
2) расчёт динамического и теплового воздействия потока газа на частицы;
3) расчёт процессов кристаллизации частиц и потерь на излучение.
Ниже дано краткое описание моделей перечисленных процессов.
IV. Краткое описание используемых моделей физико-химических процессов
Расчёт невязкой газодинамики осуществляется на основе решения автомодельной задачи о взаимодействии двух равномерных сверхзвуковых потоков (стационарный аналог схемы С.К. Годунова [9]). Для расчёта волн сжатия-разряжения вне оси струи используется обобщённое разрывное решение (соотношения на разрывах) [9-11]. Для определения параметров газового потока на оси струи при va > 0 используются соотношения автомодельной задачи о волне разряжения [12]. Расчёт параметров частиц осуществляется на основе разностной аппроксимации дифференциальных уравнений.
При расчёте струйных течений одним из определяющих факторов является модель турбулентной вязкости. В программе реализована возможность расчётов параметров турбулентного потока для различных моделей турбулентной вязкости, а именно, локальных моделей на основе первой и второй формул Прандтля, а также однопара-метрической и двухпараметрической дифференциальных моделей.
В зоне смешения со спутным воздушным потоком происходит процесс догорания компонентов выхлопа в атмосферном кислороде, в ходе которого изменяются концентрации веществ в выхлопе и температура. Символически догорание можно описать следующими эффективными (составными) реакциями:
2Н2 + 02 ^ 2Н2 0,
2С0 + О2 ^ 2С02,
4НС1 + О2 ^ 2Н2О + 2С12,
N2 + О2 ^ 2Ш.
Реально каждая из этих составных реакций представляет собой набор элементарных реакций между компонентами из вышеприведённого списка и другими промежуточными компонентами. Пример полного набора элементарных реакций можно увидеть в [13]. В данной работе используется минимально достаточная система реакций химической кинетики, определённая на основе проведённого анализа. Банк данных по константам скоростей элементарных бинарных реакций формируется на основе анализа данных о зависимости скоростей реакций от температуры из достоверных источников (например, [13, 14]).
Б
Таблица 2
Прямые реакции Обратные реакции
А В Еа А В Еа
н + о2<->он + о 1,2 Е17 -0,91 69,1 1,8 Е13 0 0
о + н2<->он + н 1,5 Е7 2 31,6
он + н2 < - > н2о + н 1,5 Е8 1,6 13,8 4,6 Е8 1,6 77,7
0 + Н20 < - > он + он 1,5 ЕЮ 1,14 72,2
н + н + м<->н2 + м 6,4 Е17 -1 0 2,2 Е14 0 402
Н + ОН + М<-> Н20 + м 1,4 Е23 -2,0 0 1,6 Е17 0 478
о + о + м<->о2 + м 1,0 Е17 -1 0 1,2 Е14 0 451
о + со + м<-> со2 + м 6,02 Е14 0 3
со + он < - > со2 + н 4,76 Е7 1,23 70
С0 + 02 < - > С02 + 0 2,5 Е12 0 47,8
В табл. 2 приведены значения параметров формулы Аррениуса К(Т) = А ■ (Т)в ■ ехр(-Еа/ВТ) для констант скоростей моделируемой системы химических реакций (размерности: А [(см3/моль)п—1 с— 1], Е а [кДж/моль]; п — порядок реакции).
Изменения концентраций компонент в элементарном объёме движущегося газа описываются системой обыкновенных дифференциальных уравнений вида
¿С,,
(З)
1, К,
(15)
где з — индекс химической реакции, г — индекс вещества,
ЕЗЗ (С) = ±СЙ(З)СКЗ), (знак «+» — для реакций образования, знак «—» — для реакций исчезновения),
КЗ — константа скорости реакции.
Сг — концентрация вещества,
М — концентрация инертных молекул (в тройных реакциях),
г(3) =0 — для бинарных реакций, г(з) = 1 — для тройных реакций.
Концентрации веществ подчиняются, кроме того, условиям сохранения химических элементов.
Для численного решения системы нелинейных дифференциальных уравнений (15) используется неявный метод простых итераций в виде
Сг^т+1 — Сг,т = [Вг((рт)] +
+Н(1 - Б)[Вг(Скт+1)],
где Н к —
шаг по времени, индекс итерации, т — номер шага во времени, Б — параметр, лежащий в интервале [0,1]. При Б = 1 метод становится явным.
Сравнение неявного метода простых итераций с параметром Б = 0,5 с явным методом Рунге-Кутта четвёртого порядка с автоматическим выбором шага показало, что неявный метод
для рассматриваемой задачи работает надежнее и устойчивее.
В ходе создания модели было проведено сравнение с результатами расчётов химической кинетики, проведёнными в других работах. В частности, было проведено сравнение с расчётами температуры горения кислородно-водородной смеси, описанными в [17]. На рис. 2 представлены результаты расчётов, приведённые в реферируемой работе (а), и результаты расчётов по представленному алгоритму (б) с использованием схемы реакций и набора параметров из той же статьи. Сопоставление позволяет говорить о хорошем согласии результатов.
0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 11ше(тз)
Рис. 2. Результаты расчёта химической кинетики: а) взято из [15]; б) расчёт предложенным алгоритмом
Присутствие в выхлопных продуктах мелкодисперсных частиц (конденсированной фазы) влияет на динамические и термодинамические характеристики струи. Для двухфазных течений характерны неравновесные процессы, обусловленные отставанием частиц по скорости и температуре от соответствующих параметров газа, а также кинетикой возможной кристаллизации жидкой конденсированной фазы при охлаждении частиц.
Взаимное влияние дисперсной фазы частиц и несущей газовой фазы определяется динамическим и тепловым взаимодействием между газом и частицами [5, 15].
Удельная сила динамического взаимодействия газа и частиц 2-й фракции определяется соотношением
-> _ IV - IV^
Ьз ~ >
'рз
где 1р,1м8з — векторы скорости газа и частиц 2-й фракции соответственно,
трз — время динамической релаксации частиц 2-й фракции.
Удельный теплопоток между газом и частицами 2-й фракции определяется соотношением
Ъз =
Т — Т-
г9 -1 «3
т т , (ги—гйд^ )2
где 1гд = 1 + -—ос--температура торможения
относительно потока,
Т8з — температура частиц 2-й фракции, тьз — время тепловой релаксации частиц 2-й фракции.
Времена динамической и тепловой релаксации частиц зависят от плотности и теплоёмкости вещества частиц, диаметра частиц, динамической вязкости и теплоёмкости газа, отношения длины свободного пробега молекул к диаметру частиц (числа Кнудсена для частиц) и др. [15].
Процесс кристаллизации жидких частиц может начаться, если при расширении температура газа становится ниже температуры плавления жидких частиц. Однако равновесный процесс кристаллизации может не осуществляться, так как скорость передачи теплоты от частиц газу и скорость образования кристаллических зародышей конечны [1]. Уравнение, определяющее уменьшение запаса теплоты плавления частицы, имеет вид
¿И„
вЬ
-2п(1р\д(Трс - ТГд(г)), Ь > ¿о,
Нте1 (Ьо) = Дна
Мр
где вр — диаметр частиц, Хд — теплопроводность газа, Трс — температура частиц, Тгд — температура торможения газа относительно частиц. В процессе кристаллизации температура частиц остаётся постоянной.
Быстрое расширение реального газа в струях сопровождается также колебательной релаксацией молекул. Колебательная релаксация молекул влияет как на газодинамические параметры струи, так и на распределение заселенностей колебательных уровней молекул, которое необходимо знать для расчёта радиационных потерь.
При течении высокотемпературной смеси время релаксации колебательных степеней свободы на один-два порядка меньше времени химических реакций [1]. Поэтому при небольших скоростях спутного потока, когда влияние головной ударной волны незначительно, химические реакции в зоне смешения не приводят к нарушению термического равновесия. В условиях интенсивной головной ударной волны время колебательной релаксации может быть сравнимо с характерным временем протекания химических реакций и, строго говоря, процессы колебательной и химической релаксаций в зоне смешения необходимо рассматривать совместно. Как показывают расчёты, для возможных параметров спутного воздушного потока (см. ниже) головная ударная волна становится интенсивной в условиях спутного воздушного потока малой плотности, когда химические реакции практически не идут (хотя модуль химических реакций работает). Поэтому в данной версии модели влияние химических реакций на процесс колебательной релаксации не учитывается. Предполагается, что колебательная релаксация определяется процессами энергообмена между колебательными и поступательно-вращательными степенями свободы при столкновении молекул.
Расчёт колебательной релаксации проводится в рамках модели многотемпературного газа [16]. Колебательная температура связана со средним числом колебательных квантов в расчёте на одну молекулу а. уравнением
ду. ,
ехр | — I - 1
1
где ду. — характеристическая колебательная температура молекулы.
Изменение среднего числа колебательных квантов а в результате колебательно-поступательного ^Т) обмена описывается релаксационным уравнением
ва.
ИГ
Е. = -
а. - а.о (Т)
а.о(Т) =
ехр
туг;
— - 1 Т-
где туг; — время колебательно-поступательной релаксации [16].
Для расчётов радиационных потерь газа и частиц используется формула Стефана-Больцма-на с эффективными коэффициентами излучения, определёнными по результатам численного моделирования излучения потоков газа с частицами.
а
1
т
V. Некоторые результаты моделирования
Сверхзвуковые струи, истекающие на умеренно нерасчётных режимах, то есть при параметре нерасчётности п = ра/рж: 1 < п < 103 (ра — давление на срезе сопла, рж — давление в спут-ном потоке), условно разделяют на три характерных участка: начальный, переходный и основной. Характерные масштабы изменения параметров в начальном, переходном и основном участках значительно различаются. На рис. 3-6 представлены результаты моделирования сверхзвуковой струи при п = 10. Массовая доля частиц 33%. Линейные размеры указаны в нормированных единицах. На рис. 3, 4 показаны поперечные относительные профили основных параметров потока для сечений в начальном, переходном и основном участках струи соответственно. По оси абсцисс отложены нормированные значения поперечной координаты. Значения параметров в относительных единицах откладываются по оси ординат. Значения
температуры показаны в тыс. К, остальные параметры для удобства показа нормированы к различным условным значениям. Сверху рисунков указаны нормированные расстояния сечений от среза сопла. На рис. 5 показаны продольные профили значений основных параметров (в тех же относительных единицах) на оси струи. По оси абсцисс отложено нормированное расстояние от среза сопла. На рис. 6 показана яркостная картина распределения температуры газа в плоскости продольного разреза (ХУ). Градации изменения яркости представлены в правом верхнем углу рисунка.
Полученные распределения параметров струи соответствуют общим физическим закономерностям изменения структуры и параметров течения вдоль струи в условиях относительно плотного спутного потока [1-5]. Начальный участок является относительно коротким, в нём довольно быстро происходит релаксация выхлопной струи от давления на срезе сопла к давлению спутного потока. Взаимодействие со спутным потоком в начальном участке происходит в довольно узком слое смешения.
12 3 4
Рис. 3. Поперечные профили основных параметров в начальном участке струи, п = 10
12 3 4
Рис. 4. Поперечные профили основных параметров в основном участке струи, п = 10
Радиус струи
Радкус потока частиц
Концентрация С02 (%/10) ✓ Температура газа ( ° /10 0 0)
Р
О 100 200 300 400 500 600 700 800 900 1 000 1 100 1 200 1 300 1 400 1 500 1 600 1 700 1 800 1 900 2 000
Рис. 5. Продольные профили значений основных параметров на оси струи, п = 10
0.0 100.0 1ЭЭ.Э 293.9 399.9 499.3 5393 3333 733.3 833.7
212 .23 23.43 1646.72 2363.97
11.271
10.021
8.771
7.521
6 271
5 021
3.771
Рис. 6. Распределение температуры газа в плоскости продольного разреза, п = 10
Переходный участок нерасчётных струй простирается вдоль продольной оси струи от сечения, в котором отражённый висячий скачок пересекает невязкую границу струи, до сечения, в котором внутренняя граница слоя смешения достигает оси струи. Скорость поперечного движения в переходном участке невелика, и внутри области невязкого течения, ограниченной внутренней границей слоя смешения, реализуется режим почти потенциального течения. Поперечные профили продольной скорости, полной энтальпии, концентраций различных составляющих осесимметричного газового потока в слое смешения переходного участка приближённо сохраняют свою форму, то есть ста-
новятся автомодельными. Геометрические и физико-химические параметры слоя смешения определяются главным образом соотношением слабо изменяющихся продольных скоростей в потенциальном ядре и спутном потоке, а также эффективностью процессов физико-химических превращений в реакциях догорания.
Основной участок сверхзвуковых струй при плотном спутном потоке простирается вниз по потоку от точки пересечения условной внутренней границы слоя смешения с осью струи. Течение в этом участке с хорошим приближением можно считать изобарическим. Основным процессом, определяющим изменение параметров струи
вдоль оси, является турбулентное перемешивание с внешним потоком, сопровождаемое процессами догорания продуктов выхлопа. Изменения диаметра турбулентной струи и значений параметров на оси струи подчиняются (приближённо, с учётом влияния частиц) интегральным соотношениям сохранения избыточного импульса, энтальпии и концентраций газовых составляющих [5].
В сверхзвуковых струях при сильно нерасчётных режимах п > 103 значение эффектов тур-
булентного течения уменьшается, взаимодействие выхлопных продуктов со спутным потоком и релаксация проходит в основном на молекулярном уровне. На рис. 7, 8 представлены результаты моделирования сверхзвуковой струи при п = 105. На рис. 7 показаны поперечные относительные профили основных параметров потока для сечения в начальном участке струи. На рис. 8 показана яр-костная картина распределения температуры газа в плоскости продольного разреза.
Рис. 7. Поперечные профили основных параметров в начальном участке струи, п = 105
Рис. 8. Распределение температуры газа в плоскости
В начальном участке (рис. 7) достаточно четко выделяются границы контактного разрыва, висячего скачка, головной ударной волны. В изоэнтро-пийном ядре, внутри области, ограниченной поверхностью висячего скачка, реализуется режим сверхзвукового истечения потока в пустоту. Боль-
40 051
продольного разреза, п = 105
шая часть массового расхода струи приходится на область сжатого слоя, заключённого между поверхностями контактного разрыва и висячего скачка. При сильно нерасчётных режимах имеет место существенное различие поступательной и колебательной температур, обусловленное за-
медленной релаксацией колебательных степеней свободы молекул. Следует отметить также существенное повышение поступательной температуры газа за головной ударной волной и в области прихода сжатого слоя на ось струи (рис. 8). Тем не менее температура в зоне смешения не превосходит 1200 К (рис. 7). Вследствие невысокой температуры и малой плотности газа скорости химических реакций в зоне смешения незначительны.
\
11 1 Темп. (г. 83). К* "3 |
X ; Кол.дол. С02,хЮ
Мол.дол. Н20.х10 !
О 1 2 3 4 5 б 7 8 9 10 11 12 13 14 15
а)
1 Р---^
0 100 200
б)
Рис. 9. Мольные доли частиц и температура газа на оси струи: а) рассчитано предложенным алгоритмом; б) взято из [4]
Для проверки модели проведено сравнение с результатами аналогичных расчётов по другим моделям. На рис. 9, 10 представлены расчётные продольные профили значений основных параметров на оси струи. Профили, показанные на рис. 9а, 10а, рассчитаны на описываемой здесь модели для значений параметров в выходном сечении сопла, приведённых в [18]. Профили, показанные на рис. 9б, 10б, взяты из [4]. На рис. 9 показаны профили мольных концентраций основных газовых составляющих и температуры газа, на рис. 10 — профили температур частиц различных радиусов, а также профили температуры и давления газа. Если учесть неизбежное различие исходных параметров расчётов, то следует отметить, что профили, представленные на рис. 9а, 10а, хорошо согласуются с профилями, представленными на рис. 9б, 10б, как в начальном участке, так и на основном участке струи. Весьма близки и количественные характеристики профилей.
0 1 2 3 4 5 б 7 8 9 10 11 12 13 14 15
а)
б)
Рис. 10. Температуры частиц на оси струи: а) рассчитано предложенным алгоритмом; б) взято из [4]
а
Рис. 11. Яркостная картина распределения температуры газа: а) рассчитано предложенным алгоритмом; б) модель фирмы АегоОр^св
На рис. 11 показаны яркостные картины распределений температуры газа в плоскости продольного разреза, рассчитанные для параметра нерасчётности п ~ 103.
Распределение, показанное на рис. 11а, получено на описываемой здесь модели, распределение, показанное на рис. 11б, — на CFD-модели (www.aero-optics.com/model.htm). Видно, что распределения хорошо согласуются.
VI. Заключение
1. Разработана компьютерная модель для расчёта пространственных распределений параметров сверхзвуковой струи газа с частицами, истекающей в спутный воздушный поток. Модель основана на численном решении методом сеток в сочетании с методом расщепления системы уравнений в частных производных, описывающих течение гетерогенного потока, и на использовании упрощённых моделей релаксационных процессов.
2. Проведено моделирование течения сверхзвукового потока газа с частицами при умеренно нерасчётных и сильно нерасчётных режимах. Полученные результаты расчётов соответствуют общим физическим закономерностям изменения структуры и параметров течения вдоль струи в зависимости от характеристик спутного потока. Сравнение результатов моделирования с результатами аналогичных расчётов по другим моделям показывает их хорошее согласие.
Литература
1. Авдуевский В.С., Ашратов Э.А, Иванов А.В., Пирумов У.Г. Газодинамика сверхзвуковых неизобарических струй. — М.: Машиностроение, 1989.
2. Сверхзвуковые газовые струи / под ред. В.Г. Дулова. — Новосибирск: Наука, 1984.
3. Бондарев Е.Н. Динамика вязких сверхзвуковых нерасчётных струй. — М.: МАИ — 1983.
4. Аэродинамика ракет. Т. 2. / под ред. М. Хем-ша, Дж. Нилсена. — М.: Мир, 1989.
5. Теория турбулентных струй / под ред. Г.Н. Абрамовича. — М.: Наука, 1984.
6. Андреев Е.П., Завелевич Ф.С., Макаров И.П. Сравнение результатов расчёта UK-излучения факела с экспериментальными данными, полученными в вакуумной камере // Оптический журнал. — 1998. — Т. 65, № 11. — С. 34-36.
7. Галицейский К.Б. Моделирование догорания высокоскоростных турбулентных струй // Физика горения и взрыва. — 2006. — Т. 42, № 2. — С. 3--9.
8. Марчук Г.И. Методы вычислительной математики. — М.: Наука, 1977.
9. Годунов С.К., Забродин А.В., Иванов М.Я., Крайко А.Н., Прокопов Г.П. Численное решение многомерных задач газовой динамики. — М.: Наука, 1976.
10. Стулов В.П. Лекции по газовой динамике. — М.: Физматлит, 2004.
11. Лунев В.В. Течение реальных газов с большими скоростями. — М.: Физматлит, 2007.
12. Черняк В.Г, Суетин П.Е. Механика сплошных сред. — М.: Физматлит, 2006.
13. Варнатц Ю, Маас У., Дибл Р. Горение. — М.: Физматлит, 2003.
14. Turns S.R. An Introduction to Combustion: Concepts and Applications. — New York: McGraw-Hill International Editions, 2000.
15. Стасенко А.Л. Физическая механика многофазных потоков. — М.: МФТИ, 2004.
16. Физико-химические процессы в газовой динамике. Т. 1.; Т. 2. / под ред. Г.Г. Черного и С.А. Лосева. — М: МГУ, 1995, 2002.
17. Androulakis I.P. Kinetic Mechanism Reduction Based on an Integer Programming Approach // AIChE Journal. — 2000. — V. 46, N. 2. — P. 361-371.
18. Нельсон Х.Ф. Влияние рассеяния на UK-излучение факелов ракет // Аэрокосмическая техника. — 1986. — № 1. — С. 128-130.
Поступила в редакцию 28.01.2008.