Известия РАН. Механика жидкости и газа, 2021, № 6, стр. 109-124

ТЕПЛООБМЕН И ИОНИЗАЦИЯ ПРИ НЕРАВНОВЕСНОМ ОБТЕКАНИИ ЗАТУПЛЕННОЙ ПЛАСТИНЫ ГИПЕРЗВУКОВЫМ ПОТОКОМ

С. Т. Суржиков ab*

a Институт проблем механики им. А.Ю. Ишлинского РАН
Москва, Россия

b Всероссийский научно-исследовательский институт автоматики им. Н.Л. Духова
Москва, Россия

* E-mail: surg@ipmnet.ru

Поступила в редакцию 16.06.2021
После доработки 22.06.2021
Принята к публикации 22.06.2021

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

Аннотация

Выполнено расчетное исследование процессов теплообмена и ионизации воздуха у поверхности затупленных пластин с радиусами затупления Rn = 0.2 и 0.66 см, обтекаемых газом в широком диапазоне параметров набегающего потока. Построены поля температур поступательных и колебательных степеней свободы молекулярных компонент, концентраций атомарных и молекулярных компонент, а также электронных концентраций вблизи поверхности. Рассчитанные плотности конвективных тепловых потоков сравниваются с данными, предсказываемыми корреляционными соотношениями Фея–Риддела, Детры–Кэмпа–Риддела и Фенстера. В рассмотренных случаях объемная доля ионизованных компонент не превосходит 1%. Однако выделены режимы слабой и сильной ионизации. Граница между этими двумя режимами определена в работе по величине критической электронной концентрации, определяемой для частот 3.3–35 ГГц радиоволн, использованной с целью измерения электронных концентраций в летном эксперименте RAMC-2.

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

Задача обтекания затупленной пластины сверхзвуковым потоком газа является одной из классических задач, привлекающих внимание специалистов в области теоретической и экспериментальной аэродинамики высоких скоростей [13]. Эта задача изучалась в разнообразных постановках фундаментальной теоретической аэротермодинамики [4, 5] и имеет большое практическое значение применительно к гиперзвуковому обтеканию элементов конструкции c затупленными кромками (крылья, органы рулевого управления, воздухозаборники) [68]. Основы теории обтекания и теплообмена плоских затупленных пластин изложены в фундаментальных монографиях [912] как итог многолетних исследований.

Среди значительного числа расчетно-теоретических работ, посвященных указанной задаче, выделяются научные направления, связанные с вязким и невязким обтеканием затупленных пластин под разными углами атаки, обтекание плоских и слегка выпуклых крыльев. Подавляющее число работ по изучению структуры поля течения вблизи затупленных пластин основано на модели совершенного газа [1, 5], но достаточно подробные исследования выполнены также с использованием моделей химически равновесного и неравновесного газа [3, 4]. В этих работах изучалось распределение давления и конвективных тепловых потоков вдоль поверхности пластин при разных числах Маха под разными углами атаки.

С использованием асимптотической теории сильного и слабого вязко-невязкого взаимодействия и численного моделирования полных уравнений Навье–Стокса получены решения для полубесконечных крыльев и крыльев конечного размера, в том числе с учетом излома образующей передней кромки [5], постановка которых является развитием задач обтекания затупленных пластин.

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

В работе [13] изучался вопрос о достоверности расчетов абсолютных величин конвективных тепловых потоков для условий экспериментов [14]. Рассматривалось обтекание сферического затупления радиусом Rn = 0.66, 1.27, 2.54 и 5.04 см сверхзвуковым потоком воздуха в диапазоне скоростей 2–5.5 км/с при условиях набегающего потока, отвечающих высотам полета H = 22 и 27 км. Указанное экспериментальное исследование было выполнено в ряду пионерских работ [1518], в которых подробно обсуждаются результаты экспериментальных исследований и рекомендуются корреляционные расчетные соотношения, которые также были многократно подтверждены.

Мотивация выполненного исследования состояла в том, что в большинстве работ по конвективному нагреву элементов конструкций проектируемых гиперзвуковых летательных аппаратов изучается распределение вдоль поверхности относительных величин, отнесенных к величине плотности конвективного теплового потока в критической точке обтекаемого затупления (см., например, [8]). При этом абсолютная величина выбирается по корреляционным соотношениям типа формулы Фея–Риддела [19].

Доводы в надежности инженерной оценки по формуле Фея–Риддела являются вполне обоснованными и оправданными. Однако вопрос о корректности расчета абсолютных величин остается актуальным, поскольку хорошо известно, что для удовлетворительного совпадения с экспериментальными данными требуется специальный выбор параметров в указанной формуле, которые не всегда хорошо известны. Это, в первую очередь, каталитические свойства поверхности, адекватное описание процессов физико-химической релаксации в пограничном слое и вблизи фронта ударной волны. В то же время проектирование перспективных летательных аппаратов требует повышения точности предсказания плотностей конвективных тепловых потоков, когда разброс в данных 20–30% может оказаться не удовлетворительным.

В работе [13] показано, что принципиальное значение также имеет правильный выбор топологии и подробности используемых конечно-разностных сеток. Примечательно, что аналогичные требования предъявляются к расчетным моделям, используемым для определения электронных концентраций вблизи обтекаемой поверхности.

В данной работе представлено дальнейшее развитие расчетной модели [13] в части определения плотностей конвективных потоков на затупленной кромке с радиусами затупления Rn = 0.2 и 0.66 см. Представлено сравнение результатов расчетов с данными корреляционных соотношений Фея–Риддела (F-R) [19], Детры–Кемпа–Риддела (DKR) [18] и Фенстера (Fe) [20] с поправкой на кривизну обтекаемой сферы или цилиндрической кромки. В качестве тестирования вычислительной модели выполнены расчеты обтекания цилиндра, затупленного по сфере, как в экспериментах [14].

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

1) Формула Дэтра–Кэмпа–Риделла (цитируется по [6])

(1)
${{q}_{w}} = 5.7 \times {{10}^{{ - 7}}}\sqrt {\frac{{(1 + \alpha ){{\rho }_{\infty }}}}{{{{R}_{n}}}}} \cdot {{\left( {\frac{{{{V}_{\infty }}}}{{100}}} \right)}^{{3.25}}}\left[ {1 - \frac{{0.00015{{C}_{{p,\infty }}}{{T}_{w}}}}{{V_{\infty }^{2}}}} \right],\,\,{\text{Вт/с}}{{{\text{м}}}^{2}},~~~~~~~~~$

2) Упрощенная формула Фея и Риддела (цитируется по [21]

(2)
${{Q}_{w}} = 1.63 \times {{10}^{5}}\sqrt {\frac{{1 + \alpha }}{{{{R}_{n}}}}} {{\left( {\frac{{{{V}_{\infty }}}}{{{{{10}}^{6}}}}} \right)}^{{3.15}}}\sqrt {\frac{{{{\rho }_{\infty }}}}{{{{\rho }_{0}}}}} ,{\text{ Вт/с}}{{{\text{м}}}^{2}},~$

3) Формула Фенстера [20]

(3)
${{q}_{w}} = 0.454 \times {{10}^{{ - 6}}}\sqrt {\frac{{1 + \alpha }}{{{{R}_{n}}}}} \sqrt {\frac{{{{\rho }_{\infty }}}}{{{{\rho }_{0}}}}} {{\left( {\frac{{{{V}_{\infty }}}}{{100}}} \right)}^{{2.862}}},{\text{ Вт/с}}{{{\text{м}}}^{2}}$
где: ${{\rho }_{\infty }}$ – плотность газа в набегающем потоке, г/см3, ${{R}_{n}}$ – радиус затупления, см, ${{C}_{{p,\infty }}}$ – теплоемкость газа при постоянном давлении, эрг/(г⋅K), ${{V}_{\infty }}$ – скорость набегающего потока газа, см/с, ${{\rho }_{0}} = {\text{1}}{\text{.23}} \times {\text{1}}{{{\text{0}}}^{{ - 3}}}$, г/см3; α = 0 при обтекании кромки, α = 1 – для сферы (геометрические факторы добавлены в данной работе).

Отметим, что явная коренная зависимость от Rn в (1)–(3) позволяет оценить плотность конвективных тепловых потоков для произвольных радиусов скругления (по крайней мере для высот H = 22 и 37 км). Однако здесь следует проявлять определенную осторожность, связанную с тем, что в разреженной среде для малых радиусов могут нарушаться условия допустимости использования модели сплошной среды. О некоторых особенностях таких условий см. в работе [22]. Для справки, в табл. 1 приведены условия в набегающем потоке на двух высотах. Видно, что на высоте 37 км длина свободного пробега молекул составляет 0.130 × 10–2 см, так что на больших высотах указанное замечание становится актуальным.

Таблица 1.

Параметры набегающего потока в экспериментах [14], соответствующие двум высотам воздушной атмосферы

Высота, км 22 37
Температура, К 219 242
Давление, Па (Тор) 4050 (30.4) 433 (3.25)
Плотность, г/см3 0.645 · 10–4 0.624 · 10–5
Скорость звука, м/с 296 312
Длина свободного пробега, см 0.126 · 10–3 0.130 · 10–2

Во второй части представлены результаты расчетов распределений электронных концентраций вблизи поверхности обтекаемой в продольном направлении затупленной пластины. Условия в набегающем потоке взяты из экспериментальной работы [14]. Использованы три модели химической кинетики высокотемпературного воздуха:

– 5-компонентная модель диссоциации (Модель Ns = 5: N, O, N2, O2, NO);

– 7-компонентная модель ионизации (Модель Ns = 7: N, O, e, N2, O2, NO, NO+);

–11-компонентная модель ионизации (Модель Ns = 11: N, O, e, N2, O2, NO, N$_{2}^{ + }$, O$_{2}^{ + }$, NO+, N+, O+).

Совместно с каждой моделью химической кинетики исследовались две модели физической кинетики: эвристическая модель Парка неравновесной диссоциации [23] и первая кинетическая модель Тринора–Мэрроуна (TM1) равновероятной диссоциации из различных возбужденных молекулярных состояний [24]. Применялась также модель локального термодинамического равновесия (ЛТР) – модель больцмановского распределения молекул по внутренним степеням свободы.

Использование различных кинетических моделей обусловлено следующим. При анализе влияния кинетической модели на рассчитываемые величины плотностей конвективных тепловых потоков в [13] и в более ранних, пионерских работах по теплопередаче в пограничном слое на плоской пластине [2528], была установлена заметная, но не очень сильная зависимость от использованных кинетических моделей диссоциации и даже ионизации. Вопрос о таком влиянии на распределения электронных концентраций остался недостаточно изученным, особенно для практически важных условий, как в работе [14].

Разница в трех кинетических моделях, упомянутых выше, состоит в следующем. В простейшей модели (Ns = 5 – число химических компонент) учитываются только процессы диссоциации. Заряженные компоненты в данной модели отсутствуют.

В модели Ns = 7 учитываются только два типа заряженных частиц: NO+ и e. Предполагается, что ионизация протекает посредством механизма ассоциативной ионизации

${\text{N}} + {\text{O}} \to {\text{NO*}} \to {\text{N}}{{{\text{O}}}^{ + }} + {{e}^{ - }}$

Заметим, что это многостадийный и недостаточно изученный процесс, поскольку предполагает предварительную диссоциацию молекулярного азота и кислорода с образованием возбужденных молекул NO*, которые являются неустойчивыми и при распаде приводят к образованию заряженных частиц. Сложность указанного кинетического процесса объясняет большую неопределенность в константах скорости этого кинетического процесса. В данной работе используется кинетическая модель [29].

В модели Ns = 11 рассматриваются 6 заряженных компонент, что означает, в частности, учет прямой столкновительной ионизации атомов и молекул. В работе [30] со ссылкой на более ранние работы этой научной группы отмечено, что в процессах ионизации столкновительные процессы играют важную роль при высоких скоростях (${{V}_{\infty }}$ > 10 км/с).

1. СИСТЕМА ИНТЕГРИРУЕМЫХ УРАВНЕНИЙ И ГРАНИЧНЫЕ УСЛОВИЯ

Рассматривается двухмерная задача обтекания затупленного по сфере цилиндра (α = 1) и затупленной пластины (α = 0) потоком химически и термически неравновесного воздуха без выделения фронта ударной волны гиперзвуковыми потоками воздуха. Используемая вычислительная модель включает в себя уравнения неразрывности, Навье–Стокса, сохранения энергии поступательного движения частиц в форме уравнения теплопроводности Фурье–Кирхгоффа, сохранения массы химических компонент и сохранения колебательной энергии в отдельных модах

(4)
$\frac{{\partial \rho }}{{\partial t}} + {\text{div}}\left( {\rho {\mathbf{V}}} \right) = 0,$
(5)
$\frac{{\partial \rho u}}{{\partial t}} + {\text{div}}\left( {\rho u{\mathbf{V}}} \right) = - \frac{{\partial p}}{{\partial x}} + {{S}_{x}}$
(6)
$\frac{{\partial \rho {v}}}{{\partial t}} + {\text{div}}\left( {\rho {v}{\mathbf{V}}} \right) = - \frac{{\partial p}}{{\partial y}} + {{S}_{y}}$
(7)
$\rho {{c}_{p}}\frac{{\partial T}}{{\partial t}} + \rho {{c}_{p}}{\mathbf{V}}\operatorname{grad} T = {{S}_{T}}$
(8)
$\frac{{\partial {{\rho }_{i}}}}{{\partial t}} + \operatorname{div} {{\rho }_{i}}{\mathbf{V}} = - \operatorname{div} {{{\mathbf{J}}}_{i}} + {{\dot {w}}_{i}},\quad i = 1,\;2, \ldots ,\;{{N}_{s}},$
(9)
$\frac{{\partial {{\rho }_{{i(m)}}}{{e}_{{{\text{V}}m}}}}}{{\partial t}} + \operatorname{div} \left( {{{\rho }_{{i(m)}}}{{e}_{{{\text{V}}m}}}{\mathbf{V}}} \right) + \operatorname{div} \left( {{{e}_{{{\text{V}}m}}}{{{\mathbf{J}}}_{{i(m)}}}} \right) = {{\dot {e}}_{{{\text{V}}m}}},\quad m = 1,\;2,\; \ldots ,\;{{N}_{{\text{V}}}}$
где диссипативная функция и правые части уравнений сохранения импульсов и энергии содержат динамические и диссипативные компоненты, а также слагаемые, обусловленные химическим реакциями, диффузией и колебательной релаксацией
${{\Phi }_{\mu }} = \mu \left[ {2\alpha {{{\left( {\frac{{v}}{y}} \right)}}^{2}} + 2{{{\left( {\frac{{\partial {v}}}{{\partial y}}} \right)}}^{2}} + 2{{{\left( {\frac{{\partial u}}{{\partial x}}} \right)}}^{2}}} \right.\left. { + {{{\left( {\frac{{\partial {v}}}{{\partial x}} + \frac{{\partial u}}{{\partial y}}} \right)}}^{2}} - \frac{2}{3}{{{\left( {\frac{{\partial u}}{{\partial x}} + \frac{{\partial {v}}}{{\partial y}} + \alpha \frac{{v}}{y}} \right)}}^{2}}} \right]$
${{S}_{x}} = - \frac{2}{3}\frac{\partial }{{\partial x}}\left( {\mu {\text{div}}{\mathbf{V}}} \right) + \frac{1}{{{{y}^{\alpha }}}}\frac{\partial }{{\partial y}}\left[ {{{y}^{\alpha }}\mu \left( {\frac{{\partial {v}}}{{\partial x}} + \frac{{\partial u}}{{\partial y}}} \right)} \right] + 2\frac{\partial }{{\partial x}}\left( {\mu \frac{{\partial u}}{{\partial x}}} \right)$
${{S}_{y}} = - \frac{2}{3}\frac{\partial }{{\partial y}}\left( {\mu {\text{div}}{\mathbf{V}}} \right) + \frac{\partial }{{\partial x}}\left[ {\mu \left( {\frac{{\partial {v}}}{{\partial x}} + \frac{{\partial u}}{{\partial y}}} \right)} \right] + 2\frac{\partial }{{\partial y}}\left( {\mu \frac{{\partial {v}}}{{\partial y}}} \right) + 2\mu \alpha \frac{\partial }{{\partial y}}\left( {\frac{{v}}{y}} \right)$
${{S}_{T}} = \operatorname{div} \left( {\lambda \operatorname{grad} T} \right) + \frac{{\partial p}}{{\partial t}} + {\mathbf{V}}\operatorname{grad} p + {{\Phi }_{\mu }} + {{Q}_{{\text{V}}}} - \sum\limits_{i = 1}^{{{N}_{s}}} {{{h}_{i}}{{{\dot {w}}}_{i}}} + \sum\limits_{i = 1}^{{{N}_{s}}} {\rho {{c}_{{p,i}}}{{D}_{i}}} \left( {\operatorname{grad} {{Y}_{i}} \cdot \operatorname{grad} T} \right)$
t – время; $x$, $y$ – ортогональные цилиндрические координаты; u, ${v}$ – проекции скорости V на оси координат x и y; $p$, $\rho $ – давление и плотность; T – температура поступательного движения частиц; $\mu $, $\lambda $ – динамический коэффициент вязкости и коэффициент теплопроводности; cp – удельная теплоемкость смеси при постоянном давлении; cp = $\sum\nolimits_i^{{{N}_{s}}} {{{Y}_{i}}{{c}_{{p,i}}}} $; ${{Y}_{i}}$ – массовая доля i-го компонента смеси; ${{c}_{{p,i}}}$, ${{h}_{i}}$ – удельная теплоемкость при постоянном давлении, связанная с поступательными и вращательными степенями свободы, и энтальпия i-го компонента смеси; ${{\rho }_{{i(m)}}}$ – плотность i-го компонента газовой смеси, обладающего m-й модой колебательного движения; ${{\dot {w}}_{i}}$ – массовая скорость химических превращений для i-го компонента смеси; ${{D}_{i}}$ – эффективный коэффициент диффузии i-го компонента смеси; ${{{\mathbf{J}}}_{i}}$ – плотность диффузионного потока i-го компонента; ${{{\mathbf{J}}}_{i}} = - \rho {{D}_{i}}\operatorname{grad} {{Y}_{i}}$; ${{N}_{s}}$ – число химических компонентов смеси газов; ${{Q}_{{\text{V}}}}$ – объемная мощность тепловыделения, обусловленная процессами VT колебательной релаксации в газовой смеси, ${{Q}_{{\text{V}}}} = - \sum\nolimits_{m = 1}^{{{N}_{{\text{V}}}}} {{{{\dot {e}}}_{{{\text{V}}m}}}} $; ${{N}_{V}}$ – число колебательных мод (в рассматриваемом случае ${{N}_{{\text{V}}}} = 3$: $m = 1$ для колебательной энергии N2, $m = 2$ для колебательной энергии O2, $m = 3$ для колебательной энергии NO); ${{\dot {e}}_{{{\text{V}}m}}}$ – источник колебательной энергии в m-й моде;

(10)
${{e}_{{{\text{V}},m}}} = \frac{{{{R}_{{i(m)}}}{{\theta }_{m}}}}{{\exp \left( {{{{{\theta }_{m}}} \mathord{\left/ {\vphantom {{{{\theta }_{m}}} {{{T}_{{{\text{V}},m}}}}}} \right. \kern-0em} {{{T}_{{{\text{V}},m}}}}}} \right) - 1}},$

${{e}_{{{\text{V}}m}}}$ – удельная энергия колебательного движения в m-й колебательной моде i-го компонента газовой смеси; ${{R}_{{i(m)}}} = {{{{R}_{0}}} \mathord{\left/ {\vphantom {{{{R}_{0}}} {{{M}_{{i(m)}}}}}} \right. \kern-0em} {{{M}_{{i(m)}}}}}$, ${{R}_{0}} = 8.314 \times {{10}^{7}}$ эрг/(K ⋅ моль) – универсальная газовая постоянная; ${{T}_{{V,m}}}$ – колебательная температура, соответствующая m-й колебательной моде; ${{D}_{{i(m)}}}$, ${{{\mathbf{J}}}_{{i(m)}}}$, ${{M}_{{i(m)}}}$ – эффективный коэффициент диффузии в многокомпонентной газовой смеси, вектор плотности диффузионного потока и молекулярный вес i-го компонента газовой смеси, обладающего m-й модой колебательного движения; ${{\theta }_{m}}$ – характеристическая колебательная температура (${{\theta }_{{m = 1}}}({{{\text{N}}}_{{\text{2}}}})$ = = 3396 K, ${{\theta }_{{m = 2}}}({{{\text{O}}}_{{\text{2}}}})$ = 2275 K, ${{\theta }_{{m = 3}}}({\text{NO}})$ = 2742 K).

Замыкающие соотношения для решаемой системы уравнений включают в себя термическое уравнение состояния идеального газа

(11)
$\frac{p}{\rho } = \frac{{{{R}_{0}}}}{{{{M}_{\Sigma }}}}T,\quad \frac{1}{{{{M}_{\Sigma }}}} = \sum\limits_i^{{{N}_{s}}} {\frac{{{{Y}_{i}}}}{{{{M}_{i}}}}} $
и калорическое уравнение состояния
(12)
$e = \sum\limits_i^{{{N}_{s}}} {{{Y}_{i}}{{e}_{i}}} ,\quad {{e}_{i}} = \int\limits_{{{T}_{0}}}^T {{{c}_{{{\text{V}},i}}}dT + {{e}_{{i,0}}}} ,$
${{e}_{{i,0}}}$ – внутренняя энергия при ${{T}_{0}}$.

С учетом того, что колебательное возбуждение молекул может отличаться от равновесного, удельная внутренняя энергия может быть записана в виде

(13)
${{e}_{i}} = \frac{3}{2}{{R}_{i}}T + {{R}_{i}}T + {{R}_{i}}\frac{{{{\theta }_{i}}}}{{\exp \left( { - \frac{{{{\theta }_{i}}}}{{{{T}_{{V,i}}}}}} \right) - 1}}$
где, однако, предполагается, что используется межъядерный потенциал вида гармонического осциллятора. Допустимость этого для рассматриваемых условий подтверждается тем, что заданная кинетическая энергия соответствует первым ~ 10 уровням возбуждения колебаний: ${{E}_{{{\text{10V}}}}}({{{\text{N}}}_{2}})$ ∼ 9.1 × 10–12 эрг, ${{E}_{{{\text{10V}}}}}({{{\text{O}}}_{2}})$ = 3.13 × 10–12 эрг, которые можно с удовлетворительной точностью описать моделью гармонического осциллятора [31].

Кроме термического и калорического уравнений состояния важным для решаемой задачи является способ расчета переносных свойств многокомпонентного газа: вязкости, теплопроводности и коэффициентов диффузии. Подробное изложение использованной модели дано в работе [13].

В набегающем со скоростью ${{V}_{\infty }}$ потоке задавались условия невозмущенного газа (${{p}_{\infty }}$, ${{\rho }_{\infty }}$), связанные термическим уравнением состояния (11) с температурой ${{T}_{\infty }}$ = 450 K и молекулярным весом ${{M}_{\Sigma }}$ = 29 г/моль. Эти граничные условия задавались на всей внешней поверхности расчетной области при $x < {{x}_{{\max }}}$, поскольку расчетная область выбиралась так, чтобы возмущения в поле течения не достигали этой границы. Эти исходные данные соответствуют экспериментам [14]. Предполагается, что в набегающем потоке газ находится в термическом равновесии.

На поверхности затупленного по сфере цилиндра задавались условия прилипания.

В соответствии с рекомендациями [13, 14] температура поверхности задавалась фиксированной (${{T}_{w}}$ = 550 K). Поверхность полагалась некаталитической (кроме некоторых специально оговоренных ниже случаев). Считалось, что на обтекаемой поверхности достигается термическое равновесие (${{T}_{V}} = {{T}_{w}}$).

Граничные условия в выходном сечении $x = {{x}_{{\max }}}$ задавались в виде $\frac{{\partial \psi }}{{\partial \xi }} = 0$, где ψ = {u, ${v},T,\rho ,{{Y}_{i}},{{e}_{{V,m}}}\} $, $\xi $ – координатная линия, подстраиваемая к линиям тока. Кроме пограничного слоя у поверхности, здесь всегда наблюдалось сверхзвуковое течение.

Подробности компьютерной модели интегрирования системы уравнений, реализованной в авторском компьютерном коде NERAT-2D, приведены в [32].

2. РЕЗУЛЬТАТЫ ЧИСЛЕННОГО МОДЕЛИРОВАНИЯ

Первая серия расчетов с использованием изложенной компьютерной модели проведена для условий экспериментов [14], приведенных в табл. 1 . По сравнению с работой [13] расширен диапазон исследованных скоростей набегающего потока до ${{V}_{\infty }}$ = 8 км/с. Результаты численного моделирования плотностей конвективных тепловых потоков Qw в критической точке сферического затупления радиусом Rn = 0.66 см показаны на рис. 1а совместно с экспериментальными данными [14] и расчетными данными, полученными по корреляционным соотношениям (1)–(3). Заметим, что разброс экспериментальных значений доходит до 50%, и значения, полученные по корреляционным соотношениям, практически охватывают облако экспериментальных данных.

Рис. 1.

Плотность конвективных тепловых потоков к сферическому (а) и цилиндрическому (б, в) скруглению на разных высотах.

Разброс значений Qw, полученных при численном моделировании с использованием трех кинетических моделей химических превращений и ионизации, а также трех моделей неравновесной диссоциации (модели ЛТР, Парка и Тринора–Мэрроуна) достигал 30%. Эти данные лежат внутри коридора, обозначенного на рис. 1а. На рис. 1а приведены расчетные данные, полученные с использованием кинетической модели Ns = 11 и модели неравновесной диссоциации ТМ1. Как и в [13], отметим большое влияние на расчетные данные подробности и топологии расчетных сеток.

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

Расчетные и экспериментальные данные на рис. 1а подтверждают надежность численного моделирования плотностей конвективных тепловых потоков с использованием компьютерной модели NERAT-2D.

На рис. 1б,в представлены результаты расчетов плотностей конвективных тепловых потоков в критической точке обтекания затупленных кромок с радиусами скругления Rn = 0.66 и 0.2 см. Как и прежде, численные данные находятся внутри коридора расчетных величин, даваемых корреляционными соотношениями (1)–(3) при α = 0. Приведенные данные также дают основания рекомендовать (1)–(3) для оценок тепловых потоков к затупленным кромкам в широком диапазоне скоростей.

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

Во-вторых, плотности конвективных тепловых потоков возрастают с увеличением плотности набегающего потока газа (на рис. 1б,в – соответствует уменьшению высоты полета от 37 до 22 км).

В-третьих, уменьшение радиуса скругления приводит к заметному возрастанию плотности конвективных тепловых потоков (сравните рис. 1б и в).

Таким образом, первая серия расчетов позволила выполнить валидацию компьютерного кода NERAT-2D во всем диапазоне экспериментально изученных скоростей [14] и выработать рекомендации по инженерной оценке тепловой нагрузки на элементы конструкции летательных аппаратов.

Вторая серия расчетов была направлена на изучение закономерностей изменения ионизационной структуры течения у затупленных пластин при изменении радиусов затупления, скорости и плотности набегающего потока. В качестве базовых условий обтекания выбраны также экспериментальные условия [14].

Общие закономерности гиперзвукового течения и ионизации воздуха у затупленной пластины рассмотрим на примере обтекания кромки с радиусом затупления Rn = 0.66 см на высоте 37 км при скорости ${{V}_{\infty }}$ = 5 и 8 км/с. На рис. 2 показаны поля температуры и давления, а на рис. 3 – плотности и продольной скорости. Структурные особенности течения являются характерными для обтекания затупленных пластин. Отметим основные из них.

Рис. 2.

Поле температуры поступательных степеней свободы (в K) и давления (Pres = ${p \mathord{\left/ {\vphantom {p {{{\rho }_{\infty }}V_{\infty }^{2}}}} \right. \kern-0em} {{{\rho }_{\infty }}V_{\infty }^{2}}}$) при обтекании кромки с радиусом скругления Rn = 0.66 см со скоростью ${{V}_{\infty }}$ = 8 км/с на высоте Н = 37 км.

Рис. 3.

Поле плотности (Rho = ${\rho \mathord{\left/ {\vphantom {\rho {{{\rho }_{\infty }}}}} \right. \kern-0em} {{{\rho }_{\infty }}}}$) и продольной скорости (Vx = ${u \mathord{\left/ {\vphantom {u {{{V}_{\infty }}}}} \right. \kern-0em} {{{V}_{\infty }}}}$) при обтекании кромки с радиусом скругления Rn = 0.66 см со скоростью ${{V}_{\infty }}$ = 8 км/с на высоте Н = 37 км.

Вблизи критической точки давление достигает давления торможения, а затем вдоль поверхности падает весьма быстро. В рассматриваемом случае – более, чем в 50 раз у самой поверхности и в 2–3 раза у фронта ударной волны.

Температура в сжатом слое у критической точки имеет характерное распределение: резкий скачок поступательной температуры во фронте отошедшей ударной волны с возникновением непротяженной релаксационной зоны, где температура поступательных степеней свободы быстро падает до температуры в сжатом слое, а температуры колебательных степеней свободы молекул N2, O2 и NO возрастают от температуры в набегающем потоке до температуры сжатого слоя (см. рис. 4). Причиной резкого падения температуры сразу за фронтом ударной волны является то, что значительная энергия поступательных степеней свободы молекул N2 и O2 тратится на их диссоциацию. Именно в этой релаксационной зоне рождается основная масса молекул NO и в результате реакции ассоциативной ионизации – рождаются электроны и ионы NO+.

Рис. 4.

Распределение температуры поступательных и колебательных степеней свободы вдоль критической линии тока при обтекании кромки с радиусом затупления 0.66 см со скоростью 8 (1) и 5 (2) км/с на высоте 37 км.

При увеличении скорости потока от 5 до 8 км/с температура в сжатом слое возрастает на 2000 K. При этом давление торможения увеличивается от 1.6 × 106 эрг/см3 до 4.0 × 106 эрг/см–1.

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

На рис. 5 показано распределение концентраций компонент частично ионизованного воздуха. При скорости ${{V}_{\infty }}$ = 8 км/с объемные доли электронов и ионов NO+ в сжатом слое достигают величин 10–3. Это означает, что степень ионизации воздуха составляет примерно 0.1%. Это весьма высокая концентрация электронов. Учитывая связь между массовыми, объемными и числовыми концентрациями

${{Y}_{i}} = {{x}_{i}}\frac{{{{M}_{\Sigma }}}}{{{{M}_{i}}}},\quad {{N}_{i}} = {{x}_{i}}\frac{p}{{1.38 \times {{{10}}^{{ - 6}}}T}},\,\,{\text{с}}{{{\text{м}}}^{{--3}}},$
находим, что при ${{V}_{\infty }}$ = 8 км/с в сжатом слое ${{N}_{e}}$ ∼ 4 × 1015 см–3, а при ${{V}_{\infty }}$ = 5 км/с – ${{N}_{e}}$ ∼ 3 × × 1014 см–3.

Рис. 5.

Распределение относительных мольных концентраций вдоль критической линии тока при обтекании кромки с радиусом затупления 0.66 см со скоростью 8 (1) и 5 (2) км/с на высоте 37 км.

Напомним, что в экспериментальном исследовании [33] использовались СВЧ-датчики, регистрирующие критические концентрации электронов на четырех частотах электромагнитного излучения: f1 = 1.116 ГГц (принадлежит так называемому спектральному диапазону “L-band”), f2 = 3.348 ГГц (“S-band”),  f3 = 10.044 ГГц (“X-band”),  f4 = 35.000 ГГц (“Ka-band”). Учитывая связь между критической частотой отражения электромагнитного сигнала (в ГГц) от электронного облака [35]

$f = 9 \times {{10}^{{ - 6}}}\sqrt {{{N}_{e}}} $
становится понятным выбор указанных частот в [14] при обнаружении критических электронных концентраций Ne, кр1 = 1.6 × 1010 см–3, Ne, кр2 = 1.4 × 1011см–3, Ne, кр3 = 1.23 × 1012 см–3, Ne, кр1 = = 1.7 × 1012 см–3.

Расположение этих датчиков вдоль поверхности экспериментального гиперзвукового аппарата RAMC позволило определить распределение максимальных значений концентраций электронов вдоль обтекаемой поверхности. Использованная в работе [32] модель аэрофизики NERAT-2D показала хорошее согласие с этими экспериментальными данными.

На рис. 6 показано распределение числовых концентраций электронов и ионов NO+ при скорости набегающего потока ${{V}_{\infty }}$ = 5 и 8 км/с. Наибольшая концентрация заряженных частиц формируется вблизи критической линии тока, а ниже по течению превалируют рекомбинационные процессы и концентрации заряженных частиц снижаются на порядки величин. Тем не менее даже при скорости ${{V}_{\infty }}$ = 5 км/с достаточно высокая концентрация электронов наблюдается далеко вниз по течению. Очевидно, что при увеличении скорости набегающего потока заряженные частицы не успевают рекомбинировать, поскольку их скорость рекомбинации становится соизмеримой со скоростью сноса вдоль поверхности, и их высокая концентрация сохраняется дольше.

Рис. 6.

Поля числовых концентраций электронов и ионов NO+ (в см–3) при обтекании кромки с радиусом затупления Rn = 0.66 см со скоростью ${{V}_{\infty }}$ = 5 (а); 8 (б) км/с на высоте Н = 37 км.

На рис. 7 показаны распределения электронных концентраций по нормали к поверхности на разных расстояниях от затупленной кромки. Очевидна корреляция распределений указанных электронных концентраций с конфигурацией поля течения у затупленной пластины. Толщина зоны заметной ионизации увеличивается по мере отхода от передней кромки. Концентрация электронов в сжатом слое заметно падает, но остается весьма высокой. С увеличением скорости от ${{V}_{\infty }}$ = 5 до 8 км/с толщина сжатого слоя уменьшается, но в среднем концентрация электронов оказывается на порядок выше.

Рис. 7.

Распределения числовых электронных концентраций по нормали к поверхности пластины с радиусом скругления 0.66 см на разных расстояниях от кромки: 1x1 = 1 см, 2x2 = 2 см, 3x3 = 4 см, 4 x4 = 6 см, 5 – x5 = 8 см, 6x6 = 10 см при скоростях 5 (а) и 8 (б) км/с. Толстые линии – высота 37 км, тонкие линии – 22 км.

Отмеченные выше закономерности в развитии ионизованной оболочки вблизи обтекаемой затупленной пластины являются весьма общими. Отметим особенности протекания процесса ионизации при увеличении плотности набегающего потока (в рассматриваемом случае – уменьшении высоты полета) и уменьшении толщины обтекаемой пластины.

На рис. 7 приведено сравнение электронных концентраций на высотах 37 и 22 км. Видно, что при скорости 5 км/с различия незначительны, а при ${{V}_{\infty }}$ = 8 км/с на малых расстояниях от кромки электронные концентрации выше примерно на порядок, а ниже по потоку – наоборот, концентрации электронов меньше в более плотном потоке.

Значительные различия наблюдаются в распределениях поступательных и колебательных степеней свободы. Показанные на рис. 8 и 9 температурные распределения свидетельствуют об увеличении степени термической неравновесности с увеличением высоты. При этом с увеличением скорости обтекания наблюдается тенденция к большей термализации, хотя абсолютные значения температур возрастают.

Рис. 8.

Распределения температур поступательных и колебательных степеней свободы молекул N2 по нормали к поверхности пластины с радиусом скругления 0.66 см на разных расстояниях от кромки: 1x1 = 1 см, 2 – x2 = 2 см, 3x3 = 4 см, 4x4 = 6 см, 5x5 = 8 см, 6x6 = 10 см при скорости 5 км/с на высоте 37 (а) и 22 (б) км.

Рис. 9.

Распределения температур поступательных и колебательных степеней свободы молекул N2 по нормали к поверхности пластины с радиусом скругления 0.66 см на разных расстояниях от кромки: 1x1 = 1 см, 2x2 = = 2 см, 3 x3 = 4 см, 4x4 = 6 см, 5x5 = 8 см, 6x6 = 10 см при скорости 8 км/с на высоте 37 (а) и 22 (б) км.

С уменьшением радиуса затупления закономерно уменьшается толщина сжатого слоя у пластины и, как следствие, – толщина ионизованного слоя. Сравнивая распределения электронных концентраций для пластин двух толщин (рис. 7 и 10) отметим, что при ${{V}_{\infty }}$ = 5 км/с у тонкой пластины степень ионизации заметно ниже, хотя при ${{V}_{\infty }}$ = 8 км/с они становятся примерно равными.

Рис. 10.

Распределения числовых электронных концентраций по нормали к поверхности пластины с радиусом затупления 0.2 см на разных расстояниях от кромки: 1x1 = 1 см, 2x2 = 2 см, 3x3 = 4 см, 4x4 = 6 см, 5x5 = 8 см, 6x6 = 10 см при скоростях 5 (а) и 8 (б) км/с на высоте 37 км.

В заключение рассмотрим типичные результаты численного моделирования при низких скоростях потока, при которых еще наблюдается заметная ионизация. На рис. 11 даны распределения концентраций электронов в сжатом слое у пластины, затупленной с радиусом 0.66 см при скорости ${{V}_{\infty }}$ = 3 км/с на высотах 22 и 37 км. Видно, что если на высоте 22 км достигаются концентрации электронов порядка 1011 см–3, то на высоте 37 км уровень максимальных концентраций заметно ниже 1010 см–3.

Рис. 11.

Распределения числовых электронных концентраций по нормали к поверхности пластины с радиусом затупления 0.66 см на разных расстояниях от кромки: 1x1 = 1 см, 2x2 = 2 см, 3x3 = 4 см, 4x4 = 6 см, 5 – x5 = 8 см, 6 – x6 = 10 см при скорости 3 км/с на высоте 22 и 37 км.

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

ЗАКЛЮЧЕНИЕ

Выполненная валидация компьютерного кода NERAT-2D на примере экспериментальных данных по конвективному нагреву сферического затупления радиусом 0.66 см при гиперзвуковом обтекании в диапазоне скоростей ${{V}_{\infty }}$ = 2–8 км/с позволила выработать рекомендации по инженерным оценкам конвективного нагрева сферических и цилиндрических затуплений с использованием трех корреляционных зависимостей.

Численный анализ процессов ионизации воздуха при обтекании затупленных пластин с радиусами скругления 0.2 и 0.66 см на высотах 22 и 37 км показал на образовании вблизи поверхности заметной ионизации сжатого слоя с уровнем концентраций электронов Ne ∼ 1011–1014 см–3.

Показано, что для пластин толщиной 0.4 и 1.32 см ионизованный слой может достигать высоты 3–4 см. По всей толщине сжатого слоя наблюдается высокая степень термической неравновесности, которая при уменьшении высоты имеет тенденцию к снижению.

Установлены наименьшие скорости обтекания, ${{V}_{\infty }}$ ∼ 3 км/с, для которых еще наблюдается заметная ионизация. Однако следует иметь в виду достаточно резкую зависимость всех параметров сжатого слоя от высоты и скорости полета.

Работа выполнена по теме государственного задания (№ госрегистрации АААА-А20-120011690135-5) и при поддержке гранта РФФИ 19-01-00515.

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

  1. Черный Г.Г. Течение газа с большой сверхзвуковой скоростью М.: Физматгиз, 1959. 220 с.

  2. Боровой В.Я., Егоров И.В., Мошаров В.Е., Скуратов А.С., Радченко В.Н. Экстремальный нагрев тел в гиперзвуковом потоке. М.: Наука, 2018. 390 с.

  3. Землянский Б.А., Лунев В.В., Власов В.И., Горшков А.Б., Залогин Г.Н., Ковалев Р.К., Маринин В.П., Мурзинов И.Н. Конвективный теплообмен летательных аппаратов. М.: Физматлит, 2014. 380 с.

  4. Лунев В.В. Течение реальных газов с большими скоростями. М.: Физматлит, 2007. 760 с.

  5. Башкин В.А., Дудин Г.Н. Пространственные гиперзвуковые течения вязкого газа. М.: Наука, Физматлит, 2000. 288 с.

  6. Агафонов В.П., Вертушкин В.К., Гладков А.А., Поляков О.Ю. Неравновесные физико-химические процессы в аэродинамике. М.: Машиностроение, 1972. 226 с.

  7. Власов В.И., Горшков А.Б., Ковалев Р.В., Лунев В.В. Тонкая треугольная пластина с притупленным носком в вязком гиперзвуковом потоке // Изв. РАН. МЖГ. 2009. № 4. С. 134–145.

  8. Ваганов А.В., Дмитриев В.Г., Задонский С.М., Киреев А.Ю., Скуратов А.С., Степанов Э.А. Оценки теплового режима малоразмерного крылатого возвращаемого аппарата на этапе его проектирования // Физико-химическая кинетика в газовой динамике www.chemphys.edu.ru/pdf/2006-11-20-002.pdf

  9. Шлихтинг Г. Теория пограничного слоя. М.: Наука, Главная редакция физико-математической литературы, 1974. 711 с.

  10. Лойцянский Л.Г. Механика жидкости и газа: Учеб. для вузов. М.: Дрофа, 2003. 840 с.

  11. Авдуевский В.С., Галицейский Б.М., Глебов Г.А., Данилов Ю.И., Калинин Э.К., Кошкин В.К., Кошмаров Ю.А., Михайлова М.М., Михайлова Т.В., Михеев Ю.С., Рыжов Ю.А., Солнцев В.П. Основы теплопередачи в авиационной и ракетно-космической технике. М: Машиностроение, 1975. 624 с.

  12. Фабрикант Н.Я. Аэродинамика. Общий курс. М.: Наука, 1964. 814 с.

  13. Суржиков С.Т. Конвективный нагрев сферического затупления малого радиуса при относительно малых гиперзвуковых скоростях // Теплофизика высоких температур. 2013. Т. 51. № 2. С. 261–276.

  14. Rose P., Stark W. Stagnation Point Heat Transfer Measurements in Dissociated Air // JAS. 1958. № 2. P. 86–97.

  15. Lees L. Laminar Heat Transfer over Blunt-Nosed Bodies at Hypersonic Flight Speeds// Jet Propulsion. 1956. № 4. P. 259–274.

  16. Лиз Л. Конвективный теплообмен при наличии подвода вещества и химических реакций / В кн. Газодинамика и теплообмен при наличии химических реакций / Под ред. В.П. Мотулевича и В.П. Ионова. М.: Изд-во Иностранной литературы, 1962. С. 13–69.

  17. Кемп Н., Роуз П., Детра Р. Ламинарный теплообмен тупых тел с потоком диссоциированного воздуха / В кн. Газодинамика и теплообмен при наличии химических реакций / Под ред. В.П. Мотулевича и В.П. Ионова. М.: Изд-во Иностранной литературы, 1962. С. 229–255. (Kemp N.H., Rose P.H., Detra R.W. 1959. JASS. V. 26. № 6. P. 421–430).

  18. Detra R.W., Kemp N.H., Riddell F.R. Addendum to Heat Transferto Satellite Vehicles Re-Entering the Atmosphere // Jet Propulsion. 1957. № 12. P. 1256–1257.

  19. Fay J.A., Riddel F. Theory of Stagnation Point Heat Transfer in Dissociated Air // JAS. 1958. № 2. (см. также в кн. Газодинамика и теплообмен при наличии химических реакций / Под ред. В.П. Мотулевича и В.П. Ионова. М.: Изд-во Иностранной литературы, 1962. С. 190–228.

  20. Fenster S.J. Stagnation-Point Heat Transfer for a New Binary Air Model Including Dissociation and Ionization // AIAA J. 1965. V. 3. № 12. P. 2189–2196.

  21. Мартин Дж. Вход в атмосферу. Введение в теорию и практику. М.: Мир, 1969. 320 с. (John J. Martin Atmospheric Reentry. An Introduction to its Science and Engineering. Prentice-Hall, Inc., Englewood Gliffs N.J.).

  22. Иванов М.С., Хотяновский Д.В., Шершнев А.А., Кудрявцев А.Н., Шевырин А.А., Ёнемура Ш., Бондарь Е.А. Эффекты разреженности при обтекании затупленной передней кромки гиперзвуковым потоком // Теплофизика и аэромеханика. 2011. Т. 18. № 4. С. 543–554.

  23. Park C. Nonequilibrium Hypersonic Aerothermodynamics. N.Y.: Wiley-Intern. Publ.б 1990. P. 358.

  24. Treanor C.E., Marrone P.V. Effect of Dissociation on the Rate of Vibrational Relaxation // Physics of Fluids. 1962. V. 5. № 9. P. 1022–1026.

  25. Стулов В.П. Пограничный слой на пластине с учетом неравновесной диссоциации // Изв. АН СССР. Механика и машиностроение. 1961. № 3. С. 5–12.

  26. Стулов В.П. Теплопередача в ламинарном пограничном слое на пластине с учетом химической неравновесности // Изв. АН СССР. Механика и машиностроение. 1961. № 6. С. 11–14.

  27. Лунькин Ю.П., Колешко С.Б. Ламинарный пограничный слой на пластине при наличии колебательной и диссоционной релаксации / В кн. Гидрогазодинамика. М.–Л.: Машиностроение, 1966.

  28. Провоторов В.П. О вязком взаимодействии на пластине в гиперзвуковом потоке неравновесно диссоциирующего газа // Ученые записки ЦАГИ. 1971. Т. 11. № 1. С. 25–32.

  29. Park C. Review of Chemical-Kinetics Problems of Future NASA Missions, I: Earth Entries // Journal of Thermophysics and Heat Transfer. 1993. V. 7. № 3. P. 385–398.

  30. Горелов В.А., Киреев А.Ю., Шиленков С.В. Неравновесные ионизационные процессы за сильными ударными волнами при высоких скоростях распространения в воздухе // Физико-химическая кинетика в газовой динамике. 2007. Т. 5. 12 с. http://chemphys.edu.ru/issues/2007-5/articles/55/

  31. Ельяшевич М.А. Атомная и молекулярная спектроскопия. М.: Эдиториал УРС, 2000. 896 с.

  32. Суржиков С.Т. Компьютерная аэрофизика спускаемых космических аппаратов. Двухмерные модели. М.: Наука, 2018. 543 с.

  33. Akey N.D., Cross A.E. Radio Blackout Alleviation and Plasma Diagnostic Results from a 25000 Foot per Second Blunt Body Reentry. NASA TN D-5615. February, 1970.

  34. Суржиков С.Т. Двумерный численный анализ ионизации потока в летном эксперименте RAM-C-II // Химическая физика. 2015. Т. 34. № 2. С. 24–42.

  35. Гинзбург В.Л. Распространение электромагнитных волн в плазме. М.: Изд-во URSS, 2015. 688 с.

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