МОДЕЛИРОВАНИЕ ПРОЦЕССОВ
УДК 533.6.011.8
А. Л. Железняков а, В. В. Кузенов, А. С. Петрусев, С. Т. Суржиков
ЧИСЛЕННЫЙ АНАЛИЗ КОНВЕКТИВНОГО НАГРЕВА ДВУХ МОДЕЛЕЙ СПУСКАЕМЫХ КОСМИЧЕСКИХ АППАРАТОВ
Приведены результаты численного моделирования обтекания сверхзвуковым потоком совершенного газа двух моделей космических спускаемых аппаратов: КА сегментально-конической формы и разрабатываемого аппарата MSL (Mars Science Laboratory). Проанализированы результаты, полученные при различных условиях в набегающем потоке во всей области течения: от головной ударной волны до дальнего следа. Исходные данные, используемые в расчетах, отвечают условиям экспериментальных исследований на импульсных аэродинамических трубах.
Ключевые слова: численное моделирование, аэротермодинамика, газовая динамика, теплообмен, сверхзвуковые течения, спускаемые космические аппараты.
К одной из наиболее важных проблем компьютерной аэротермодинамики космических аппаратов (КА) относится математическое моделирование сложных (двух- и трехмерных) элементов течения у поверхности спускаемых КА (отрывные течения и течения присоединения, ламинарно-турбулентный переход, сложные ударно-волновые взаимодействия). Отметим, что необходимость в определении параметров газа за головной ударной волной и в зоне отрывных течений газа, которые наблюдаются в донной области спускаемых КА, связана с существенной перестройкой первоначально невозмущенного течения. Эта перестройка течения сопровождается значительным изменением теплового режима и аэродинамических характеристик вблизи поверхности спускаемых КА. Сложность рассматриваемых течений и геометрических форм летательных аппаратов приводит к необходимости разработки новых методов расчета многомерных (включая трехмерный случай) и осесимметричных течений. Использованные в настоящей работе методы позволяют находить приближенное решение усредненных по Рейнольдсу уравнений Навье-Стокса на структурированных и неструктурированных расчетных сетках в широком диапазоне физических параметров набегающего потока. Они дают возможность проводить расчетно-теоретические исследования для обеспечения опытно-конструкторских разработок аэротермодинамики КА, предназначенных для входа в плотные слои атмосфер планет и возвращения на Землю [1, 2].
Постановка задачи. Определение газодинамических параметров в основной области течения и вблизи поверхности спускаемого КА
проводится на основе модели вязкого теплопроводного газа с использованием трехмерных уравнений Навье-Стокса, усредненных по Рей-нольдсу. Расчетная область включает в себя поле течения в невозмущенном потоке, за фронтом ударной волны и в следе за обтекаемым телом. Течение во всей расчетной области предполагается ламинарным. Для решения рассматриваемой системы уравнений в настоящей работе расчетная область представляется в виде регулярной (структурированной) и неструктурированной сеток.
Наиболее часто для решения задач аэротермодинамики применяются регулярные (структурированные) сетки. В этом случае сетка представляет собой упорядоченную по определенным правилам структуру с явно выраженными сеточными направлениями.
При построении регулярных сеток в геометрически сложных физических областях применяют преобразование координат, связанных с поверхностью тела. Такое преобразование, эквивалентное введению криволинейной системы координат общего вида, позволяет построить равномерную расчетную сетку в преобразованном (вычислительном) пространстве. При этом физические границы области совпадают с координатными линиями в вычислительной области, а ячейки сетки являются криволинейными параллелепипедами (случай трехмерного моделирования). Однако при этом в системе уравнений Навье-Стокса, записанных в криволинейных координатах, появляются дополнительные слагаемые — параметры преобразования, определяющие отображение физической области на пространство обобщенных координат.
В случае дискретизации системы уравнений Навье-Стокса с использованием квазиортогональной структурированной сетки близки к нулю некоторые параметры преобразования — компоненты метрического тензора преобразования (матрицы Якоби), находящиеся не на главной диагонали данного тензора. В этом случае наблюдается уменьшение погрешности аппроксимации в сравнении с основным случаем регулярной сетки и, следовательно, повышение точности получаемого решения.
Важная особенность неструктурированных сеток — произвольное расположение узлов сетки в физической области. Эта особенность обусловливает основное преимущество неструктурированных сеток по сравнению с регулярными сетками, которое заключается в большей гибкости при дискретизации физической области сложной формы. При этом процесс построения неструктурированной сетки в области сложной пространственной формы происходит во много раз быстрее по сравнению с регулярным случаем (когда часто необходимо использовать многоблочный вариант метода построения сетки). Произвольность расположения узлов сетки следует понимать в том смысле, что отсутствуют сеточные направления и нет структуры сетки, подобной регулярным сеткам. В случае нерегулярных сеток узлы объединяются в многогранники произвольной формы (тетраэдры и призмы — для трехмерного моделирования).
Более подробно с проблемами построения регулярных и нерегулярных сеток можно познакомиться в работах [3-5].
Для описания течения газа использовались записанные в векторной форме в декартовой системе координат хуг трехмерные нестационарные усредненные по Рейнольдсу уравнения Навье-Стокса, замыкаемые с помощью соотношений модели вихревой вязкости и теплопроводности:
д( дС дН _ о
дЬ дх ду дг
Вектор консервативных переменных (( и векторы потоков Р, СС, Н имеют следующий вид:
( Р \ ( ри \
Q =
pu pv pw \pE
G =
F =
puu + P - Txx
PUV - Txy PUW - Txz
\ (pE + P) u - UTxx - VTxy - WTxz + qx J
( pv
pvu - Tyx pvv + P - Tyy pvw - Tyz
\
H =
\ (pE + P) v - uTyx - vTyy - WTyz + qy )
pw
pwu - Tzx pwv - Tzy pww + P - Tzz
\ (рЕ + Р) т - ит2Х - ит2у - тт22 + qz )
Составляющие вектора теплового потока q^ и компоненты тензора вязких напряжений т находят из следующих соотношений:
дГ дГ дГ Г ди ди дт
qx = -Хтт; qy = -хтт; qz = -хтт; й1уУ = + + тт;
дх ду дг дх ду дг
2 Г ди 2 Г ди
Тхх = - 3 Р шук + 2р—; Туу = - 3 р ¿1у\/ + 2р—;
2 Л- Ъ о dw Tzz = -ё Р dlv v + 2p^—;
3 dz
du dv
Tyx = Txy = p о + О ; Tzx = Txz = Р I о _ +
dy dx
i dv övj\
Tzy = Tyz = + gyj
du dvj\
-
dz dx 1
В приведенных уравнениях введены следующие параметры: Ь — время; р — плотность; и, V, ы — составляющие скорости V в коорди-
. и2 + V2 + ы2
натных направлениях (х, у, г); Р — давление; Е = е +-----
2
полная энергия единицы массы; Л — коэффициент теплопроводности; ц — динамический коэффициент вязкости; е — удельная внутренняя энергия; Т — температура. Газовая среда, натекающая на спускаемый КА, рассматривается как совершенный газ с 7 = 1,4.
Как уже отмечалось ранее, в настоящей работе применяются две различные численные методики решения поставленной задачи.
В первой из них реализуется принцип расщепления по физическим процессам полной системы уравнений динамики вязкого и теплопроводного газа. Первую группу составляют уравнения Навье-Стокса, а вторую — уравнения, описывающие энергетическое состояние газа. Специфика данного расчетного подхода заключается в алгоритмическом решении процедуры расщепления. Для численного интегрирования усредненных по Рейнольдсу уравнений Навье-Стокса используются АИБМ конечно-разностные схемы [6]. Полностью неявные конечно-разностные схемы применяются для численного интегрирования уравнения сохранения энергии, которое формулируется в виде уравнения относительно температуры [7]. Такой подход к решению обсуждаемых задач аэротермодинамики КА имеет ряд существенных преимуществ по сравнению с полностью явными или неявными численными методиками.
Вторая методика основана на решении уравнений Навье-Стокса, усредненных по Рейнольдсу, записанных в декартовой системе координат хуг, методом контрольного объема. В этом случае применяется нерегулярная сетка, и вся расчетная область разбивается на конечное число непересекающихся контрольных объемов — криволинейных параллелепипедов (тетраэдров), которые покрывают всю расчетную область вплоть до ее границы.
Для каждого контрольного объема V данная система уравнений может быть записана в интегральной форме:
дъЩ я^йф + /// {дх + ^ + =°
Применяя для второго слагаемого теорему Остроградского-Гаусса, получаем следующую систему уравнений:
v f + / = »■
пг
1
где Я г = у, J J J Яйх йу йг — среднее значение вектора консерватив-
V
ных переменных Я в г-м контрольном объеме; V — объем замкнутой
области (контрольного объема); П» — поверхность ¿-го контрольного объема; дБ — элемент поверхности; Ф = , б, И^ — обобщенный вектор потоков, пересекающих поверхность П» ¿-го контрольного объема; п — единичный вектор внешней нормали к элементу поверхности дБ.
Вычисляя интеграл ^ ^ Ф • п^ дБ по границе контрольного объема V как сумму произведений значений вектора потока Ф в центрах граней на площади граней Бj, запишем приведенное уравнение в следующем дискретном виде:
мг
V- ef+Е
j=i
Ф^и) Sj
= 0, г = 1,N,
dQ i ■ rr тг „
где —--усредненное по г-му объему Vi значение производной по
времени t консервативной переменной; ^ Ф • n^ Sj — аппроксимация
потока консервативной переменой через j-ю грань i-го контрольного объема; Mi — число граней контрольного объема; N — число контрольных объемов в расчетной области.
В методе контрольного объема для определения газодинамических переменных в произвольный момент времени t необходимо рассчитать потоки массы, импульса и энергии, протекающие через боковые поверхности контрольного объема. Эти потоки находят путем решения автомодельной задачи о распаде произвольного разрыва или путем использования метода AUSM. Для определения положения в пространстве центра тяжести Pi криволинейного параллелепипеда (гексаэдра), вычисления объема V гексаэдра и нахождения единичного вектора внешней нормали n к элементу поверхности dS используются формулы, приведенные в работе [8] (в главе, посвященной методу контрольного объема).
Результаты численного моделирования. Разработанные численные методики были применены для расчета аэротермодинамики экспериментальных моделей КА: MSL (Mars Science Laboratory) и сегментально-конического КА (рис. 1).
В набегающем на экспериментальную модель КА невозмущенном потоке были приняты следующие значения газодинамических параметров.
1. При расчетах обтекания модели MSL задавались исходные данные, отвечающие экспериментальному исследованию [9]: число Маха M^ = 10,3; давление Рж = 2079,3 Па; температура Тх = 49,5 K (расчеты выполнены на регулярных сетках).
2. Для той же модели проводили расчеты при M^ = 10,3; Рж = = 1068,3 Па; Тж = 51,9 K (расчеты выполнены на регулярных сетках).
j
Рис. 1. Расчетные регулярная (а) и нерегулярная (б) сетки
3. Исходные данные для расчетов модели сегментально-конического КА (М^ = 19,1; Р^ = 40,2 Па; Т^ = 31,1 К) соответствовали условиям в экспериментах [10].
Расчеты выполнены на регулярной сетке (рис. 1, а) с числом узлов по трем пространственным переменным 150x150x300 и на нерегулярной сетке с 2 млн расчетных ячеек (рис. 1, б).
На рис. 2, а и б (см. 3-ю полосу обложки) показаны распределения чисел Маха для двух типов спускаемых КА (расчеты выполнялись на неструктурированных сетках). В рассматриваемом случае реализуется осесимметричная картина течения: вблизи лобовой поверхности спускаемого аппарата формируется головной скачок уплотнения и пограничный слой, между ними располагается область слабовязкого течения, в области перехода от лобовой поверхности к задней части спускаемого аппарата формируется область волн разрежения, далее по потоку наблюдается "горло" следа с примыкающими к нему хвостовыми скачками, замыкается течение ближним и дальним следами. Вблизи задней поверхности аппаратов наблюдается торообразная зона возвратно-вихревого течения. За головным скачком уплотнения вблизи поверхности спускаемого аппарата значение температуры ударного сжатого газа достигает нескольких тысяч градусов. Отметим, что в реальных условиях обтекания, когда в ударном слое достигаются более высокие температуры, для правильного описания картины течения газа, его необходимо рассматривать как многокомпонентную химически реагирующую смесь. Однако для указанных исходных данных используемое приближение совершенного газа правомерно, что было проверено также в специальных вычислительных экспериментах.
При изменении угла натекания набегающего потока (угла атаки) картина течения заметно меняется (рис. 2, в... е, см. 3-ю полосу обложки; неструктурированная сетка): форма головной ударной волны теперь становится симметричной по отношению к плоскости симметрии, проходящей через ось симметрии спускаемого аппарата, критическая точка при соответствующем изменении угла атаки а > 0 опускается в нижнюю часть (см. рис. 2, в. . . е, 3-я полоса обложки) горло следа слабо выражено, возвратно-вихревое течение в донной области перестает быть осесимметричным. Из распределений газодинамических параметров следует, что при увеличении угла атаки головная ударная волна прижимается к нижней (по отношению к оси симметрии) части спускаемого КА. Это приводит к тому, что ударно сжатый и нагретый газ все в большей степени попадает на заднюю поверхность спускаемого аппарата и сильно меняет характер теплообмена не только на лобовой поверхности, но и в донной области (рис. 3, см. 4-ю полосу обложки).
Важной особенностью течения при ненулевых углах атаки является появление за задней частью спускаемого КА в ближнем следе малых неосесимметричных областей, в которых сходятся линии тока (сингулярная точка течения, рис. 4). Каждая из этих областей в пространстве
Иллюстрации к статье А.Л. Железняковой, В.В. Кузенова, А.С. Петрусева, С.Т. Суржикова «Численный анализ конвективного нагрева двух моделей спускаемых космических аппаратов»
Рис. 2. Распределение числа Маха для MSL-аппарата (а, в, г) и сегментально-конического КА (б, д, е) при (а, б), (в, д) и (г, е)
соответственно
Иллюстрации к статье А.Л. Железняковой, В.В. Кузенова, А.С. Петрусева, С.Т. Суржикова «Численный анализ конвективного нагрева двух моделей спускаемых космических аппаратов»
Рис. 3. Распределение конвективных тепловых потоков вблизи лобовой (а) и задней (б) поверхностей MSL-аппарата при угле атаки (структурированная сетка)
Рис. 4. Распределения линий тока на задней поверхности сегментально-конического спускаемого КА (а, б) и аппарата MSL (в, г) при а = 10° (а, в)
и 20° (б, г)
представляет собой вихревую линию или поверхность. Топология такого сложного трехмерного течения на качественном уровне описана в работе [11]. Сравнение распределений числа Маха, приведенных на рис.2,в, г (см. 3-ю полосу обложки) показывает, что возвратно-вихревое течение у разных по форме аппаратов наблюдается в разных пространственных областях: в нижней части (по отношению к оси симметрии) донной области для МБЬ-аппарата (см. рис. 2, в, г, 3-я полоса обложки) и в верхней части для сегментально-конического аппарата (см. рис. 2, д, е, 3-я полоса обложки).
Таким образом, можно сделать вывод, что пространственное расположение сингулярных точек и поверхностей зависит от угла атаки, теплофизических условий в натекающем потоке и пространственной формы спускаемого аппарата.
На рис. 5 показаны распределения давления, продольной скорости и температуры в плоскости симметрии течения для модели МБЬ (регулярная сетка [12]) при следующих значениях аэротермодинамических параметров в набегающем потоке: а = 5°, М^ = 10,31; Р^ = 2079,3 Па; Тх = 49,5 К. Поверхностное распределение конвективного теплового потока для аналогичной расчетной сетки и условий в набегающем потоке приведены на рис. 3, а, б (см. 3-ю полосу обложки).
Из этих распределений следует, что максимальные значения давления Р & 210 Па, температуры Тш & 1100 К и конвективного теплового потока Цщ & 10 Вт/см2 реализуются в точке торможения потока на лобовой поверхности спускаемого КА. На задней поверхности спускаемого аппарата наибольшие значения Т,ш & 700 К, дш & 0,5 Вт/см2 достигаются в зоне сингулярности потока и в нижней (по отношению к оси симметрии) части спускаемого аппарата, где максимально влияние головной ударной волны. Сравнение интенсивности конвективного нагрева, предсказываемого расчетами (регулярной сетки, сегментально-конического спускаемого КА), с экспериментальными данными, приведенными в работе [10], показано на рис. 6. В этом случае максимальные значения конвективного теплового потока также наблюдаются в точке торможения, но уровень их значений существенно ни^ке и составляет & 1 Вт/см2.
Выводы. С использованием разработанных трехмерных кодов, предназначенных для расчета аэротермодинамики спускаемых КА выполнены расчеты аэротермодинамики двух моделей КА, имеющих разную геометрическую форму. Выполненные расчеты дают представление о топологии поля течения вблизи КА, обтекаемого под различными углами атаки. Расчетные значения конвективных тепловых потоков удовлетворительно согласуются с имеющимися экспериментальными данными. Вычислительные эксперименты, проведенные при подготовке настоящей статьи (использование разных методов, структурированных и неструктурированных сеток, разных расчетных областей и т.д.), приводят к выводу, что в исследовании закономерностей пространственного обтекания спускаемых КА вязким
то Ч"
q, Вт/см2
1,2 -,
0,3 -
0,9 -
0,6 -
0
2
3
S/R
Рис. 6. Конвективный поток на подветренной поверхности сегментально-конического спускаемого аппарата. Продольная координата поверхности (от передней критической точки) отнесена к радиусу Миделя R:
M^ = 19,8, а = 0, 1 — эксперимент, 2 — расчет на регулярной сетке
теплопроводным газом сделаны лишь первые шаги. Для решения указанных задач в будущем понадобятся комплексные решения ряда задач аэротермодинамики вычислительной аэрофизики, теории численных методов, использование новых компьютерных технологий.
Работа выполнена в рамках проекта РФФИ № 07-01-00133, международного проекта РФФИ № 09-08-92422-КЭ_а, программы фундаментальных исследований РАН и программы Минобразования РФ РНПВШ2.1.1/4693.
СПИСОК ЛИТЕРАТУРЫ
1. Surzhikov S. T. 2D CFD/RGD model of space vehicles // Proc. of the Int. Workshop on Radiation of High Temperature Gases in Atmospheric Entry. October 2003, Lisbon, Portugal, European Space Agency, SP-533, 2003. - P. 95-102.
2. Surzhikov S. T. Numerical simulation of heat radiation generated by entry space vehicle // AIAA 2004-2379, 2004. - 11 p.
3. L i s e i k i n V. D. Grid generation methods. - Berlin: Springer, 1999.
4. Скворцов А. В. Обзор алгоритмов построения триангуляции Делоне // Вычислительные методы и программирование. - 2002. - № 3. - C. 14-39.
5. ГаланинМ. П., Щ е г л о в И. А. Разработка и реализация алгоритмов трехмерной триангуляции сложных пространственных областей: прямые методы: Препринт ИПМ им. М.В. Келдыша РАН, 2006.
6. Edwards J. R., Lion M. -S. Low-diffusion flux-splitting methods for flow at all speeds // AIAA Journal. - 1998. - Vol. 36. No. 9. - P. 1610-1617.
7. Суржиков С. Т. Физическая механика газовых разрядов. - М.: Изд-во МГТУ им. Н.Э. Баумана. - 639 c.
8. Суржиков С. Т. Тепловое излучение газов и плазмы. - М.: Изд-во МГТУ им. Н.Э. Баумана. - 543 c.
9. Brian R. Hollis, Arnold S. Collier. Turbulent aeroheating testing of Mars science laboratory entry vehicle in Perfect-Gas Nitrogen // AIAA Paper 2007-1208.
10. B o r o v o i V. Y a., S k u r a t o v A. S., S u r z h i k o v S. T. Study of convective heating of segmental-conical Martian descent vehicle in shock wind tunnel // AIAA Paper 2004-2634. - Reno, NV, 2004.
11. T o b a k M. and P e a k e D.J. Topology of three-dimensional separated flows // Annual Review of Fluid Mechanics, 14: 61-85, 1982.
12. Суржиков С. Т. Аналитические методы построения конечно-разностных сеток для расчета аэротермодинамики спускаемых космических аппаратов // Вестник МГТУ Серия "Машиностроения". - 2004. - № 2. - C. 24-50.
Статья поступила в редакцию 12.03.2009
Александра Львовна Железнякова окончила МГТУ им. Н.Э. Баумана в 2008 г. Сотрудник лаборатории "Радиационная газовая динамика" Института проблем механики им. А.Ю. Ишлинского РАН. Автор 6 научных работ в области теплофизики и вычислительной газовой динамики.
A.L. Zheleznyakova graduated from the Bauman Moscow State Technical University in 2008. Researcher of Laboratory for Radiative Gas Dynamics of the Institute for Problems in Mechanics n.a. A.Yu. Ishlinsky, RAS. Author of 6 publications in theory of heat and mass transfer and computational fluid dynamics.
Виктор Витальевич Кузенов окончил в 1983 г. Московский государственный университет им. М.В. Ломоносова. Канд. техн. наук, старший научный сотрудник лаборатории "Радиационная газовая динамика" Института проблем механики им. А.Ю. Ишлинского РАН, доцент кафедры "Теплофизика" МГТУ им. Н.Э. Баумана. Автор более 70 научных работ в области теплофизики и радиационной газовой динамики.
V.V. Kuzenov graduated from the Lomonosov Moscow State University in 1983. Ph. D. (Eng.), senior researcher of "Radiation Gas Dynamics" laboratory of A.Yu. Ishlinskii Institute for PGroblems of Mechanics, RAS, assoc. professor of "Thermal Physics" department of the Bauman Moscow State Technical University. Author of more than 70 publications in the field of thermal physics and radiation gas dynamics. Андрей Сергеевич Петрусев окончил в 1989 г. МФТИ. Канд. физ.-мат. наук, старший научный сотрудник лаборатории "Радиационная газовая динамика" Института проблем механики им. А.Ю. Ишлинского РАН, доцент кафедры "Физическая и химическая механика" МФТИ. Автор более 40 научных работ в области численного моделирования физико-химических процессов в газовых потоках. A.S. Petrusev graduated from the Moscow Physics and Technology Institute in 1989. Ph. D. (Phys.-Math.), senior researcher of "Radiation Gas Dynamics" laboratory of A.Yu. Ishlinskii Institute for Problems of Mechanics, RAS, assoc. professor of "Physical and Chemical Mechanics" department of the Moscow Physics and Technology Institute. Author of more than 40 publications in the field of numerical simulation of physical and chemical processes in gas flows.
Сеpгей Тимофеевич Суржиков окончил МВТУ им. Н.Э.Баумана в 1975г., МГУ им. М.В. Ломоносова в 1980 г. Чл.-кор. РАН, д-p физ.-мат. наук, заведующий лабоpатоpией "Радиационная газодинамика" Института пpоблем механики им. А.Ю. Ишлинского РАН, заведующий кафедрой "Физическая и химическая механика" МФТИ, пpофессоp кафедpы "Теплофизика" МГТУ им. Н.Э. Баумана. Автоp более 400 научных pабот в области теплофизики и pадиационной газодинамики. S.T. Surzhikov graduated from the Bauman Moscow Higher Technical School in 1975 and Lomonosov Moscow State University in 1980. Dr. Sc. (Phis.), Head of the Radiative Gas Dynamics Laboratory of the Institute for Problems in Mechanics n.a. A.Yu. Ishlinsky of Russian Academy of Sciences. Head of the Chair of Physical and Chemical Mechanics of Moscow Physics and Technology Institute. Author of more than 400 publications in radiative gas dynamics and theory of heat and mass transfer.