CC BY
66
_______ УЧЕНЫЕ ЗАПИСКИ ЦАГИ
Том XXVIII ' 199 7
М2
УДК 533.6.011.8
ИССЛЕДОВАНИЕ ГИПЕРЗВУКОВОГО ОБТЕКАНИЯ ПЛОСКОЙ ПЛАСТИНЫ НА ОСНОВЕ СПЛОШНОСРЕДНОГО И КИНЕТИЧЕСКОГО ПОДХОДОВ
И. В. Егоров, А. И. Ерофеев
Представлено сопоставление результатов математического моделирования гиперзвукового обтекания плоской пластины, расположенной под нулевым углом атаки, полученных на основе численного решения полных уравнений Навье — Стокса и кинетического уравнения, решаемого методом статистического моделирования Монте-Карло. Рассмотрено влияние вида граничных условий на стенке на решение уравнений Навье — Стокса: Представлен анализ силового и теплового воздействия потока газа на поверхность, а также поведение газодинамических переменных на твердой поверхности и в поле течения.
В настоящее время для исследования течений разреженного газа широко используются методы численного моделирования. Для моделирования течения разреженного газа применяются методы Монте-Карло решения уравнения Больцмана, а для сплошносредного режима (т. е. при достаточно малых числах Кнудсена Кп = 1/Ь, I — средняя длина свободного пробега молекул в невозмущенном потоке, Ь — характерный размер возмущенной области течения) — уравнения Навье — Стокса с различными граничными условиями на стенке, а также другие сплошносредные модели с нелинейными зависимостями переносных свойств от градиентов газодинамических переменных [1] — [9].
Часто сплошносредные модели используются для таких режимов обтекания, когда формально нарушается область их применимости и встает вопрос о достоверности получаемых результатов. Поэтому чрезвычайно важным является прямое сопоставление решений, полученных при кинетическом и сплошносредном подходах. В данной работе такое сопоставление проведено на примере задачи о гиперзвуковом обтекании плоской пластины, расположенной под нулевым углом атаки.
На сегодняшний день одним из основных методов решения уравнений Больцмана является метод статистического моделирования
Монте-Карло [5]. С появлением высокопроизводительных ЭВМ применение метода статистического моделирования все дальше продвигается в область сплошносредных режимов обтекания (Кп«1), что делает реальным исследование границ применимости сплошносредного подхода для исследования течений разреженного газа.
К настоящему времени разработаны эффективные алгоритмы [6] численного интегрирования полных уравнений Навье — Стокса методом сквозного счета, обеспечивающие свойство консервативности. Данные методы основаны на полностью неявных монотонных разностных схемах второго и выше порядка точности и разрешают многие математические трудности, возникающие при численном анализе задач обтекания.
Для повышения эффективности сплошносредных моделей в переходных по числу Кнудсена областях применяются два основных подхода: использование условий скольжения и скачка температуры на твердой поверхности тела для уравнений Навье — Стокса [7] и нелинейная зависимость тензора напряжений от тензора скоростей деформаций [8], [9]. При этом первый из подходов, в частности, устраняет математическую особенность в передней кромке острых тел.
Сопоставлению кинетических и слошносредных моделей посвящено много работ (см., например, [10]—[12]). Однако основным недостатком многих исследований является то, что в них для сравнения используются данные из работ, в которых не приводится полный набор входных параметров численного анализа задачи, поэтому сравнения проводятся в узком наборе определяющих параметров задачи.
Одной из наиболее последовательно изучаемых задач является моделирование течения около плоской пластины, установленной параллельно направлению свободного потока. Несмотря на большие успехи в решении этой проблемы [13], [14], имеется ряд неизученных особенностей структуры поля течения как с кинетической, так и со сплошносредной точки зрения.
В настоящей работе проведено исследование гиперзвукового обтекания плоской пластины, расположенной под нулевым углом атаки. Расчеты проведены на основе кинетического подхода методом Монте-Карло и на основе уравнений Навье — Стокса при значениях числа Ие«, = 10000, температурного фактора ^ = 0,2 и числах Маха М«, = 10 и 23. Для М„ = 23 рассматривался как одноатомный газ, так и газ с вращательными степенями свободы молекул.
1. Постановка и решение задачи для уравнений Навье — Стокса. Математическая постановка задачи обтекания произвольного двумерного тела потоком вязкого теплопроводного газа сводится к системе уравнений Навье — Стокса с соответствующими граничными условиями. Уравнения Навье — Стокса в произвольной криволинейной системе координат (|, -п):
* = *(£,т0, у =
где X, у форме:
декартовы координаты, можно записать в дивергентной
д(2 дЕ дв л
-=- + — + — = 0. дt дг\
Здесь О — вектор консервативных зависимостей переменных задачи, Е и в — векторы потоков в криволинейной системе координат. Векторы О, Е и б связаны с соответствующими векторами Ос, Ес и Сс в декартовой системе координат формулами
й = 10с, Е
,№1+с4>- с=,(£‘5+с‘0>’
в которых 3 = д(х, у) / д(%, т\) — якобиан преобразования.
Декартовы компоненты векторов Ос, Ес и Сс для двумерных уравнений Навье — Стокса имеют вид
р
0с = ри II
ри
е
р и
ри2 + р - т
XX
рии
риН - итж - хПф - X
дТ
дх
Сс =
ри
рт - Тху Р V2 +Р~ Хуу
рьН - их™ - их™ - X
дТ
ду
где р Р
плотность, и, V
декартовы компоненты вектора скорости, давление, е- р(СуТ + (и2 + и2) / 2) — полная энергия на единицу объема, Н = срТ + (и2 + V2) / 2 — полная энтальпия, ср и с„ — теплоемкости при постоянном давлении и объеме, х — тензор напряжений с компонентами
-—сЦуК + 2 3
*уу
1ху — Хух - Ц
ди ду ду дх
йЬ/У + 2
3
ду)>
V— вектор скорости, ц и X, — коэффициенты вязкости и теплопроводности.
Система уравнений для совершенного газа замыкается уравнением состояния р = рКТ/М, где Л — универсальная газовая постоянная, М — молярный вес газа. При определении коэффициента вязкости в
о т
данной работе йспользовалась зависимость ц./цв = (Т/Тх) ’ , что обеспечивало согласование с моделью межмолекулярного взаимодействий, принятой при решении задачи методом Монте-Карло (см. п. 2). Коэффициент теплопроводности определялся из соотношения Рг = цср / X = 0,7.
Задача сверхзвукового обтекания плоской пластины, установленной параллельно направлению свободного потока газа, решалась в прямоугольной области. На входной и верхней границах расчетной области параметры потока задавались равными параметрам в невозмущенном потоке, на выходной хранице использованы «мягкие» условия экстраполяции искомых газодинамических переменных. На нижней 1ранице при х < 0 рассматривалось условие симметрии потока, при х £ 0 использовалось условие прилипания
и = V = 0
и равенство температуры газа температуре стенки Т = Ту,. Наряду с этими граничными условиями рассмотрены также (максвелловские) условия скольжения и скачка температуры [7]:
_ 2 - 9 5л ди Т 2-а 75л Хя дТ
0 16 *1%’ Т^~ + а 128 ~Т~ду’
где X, — средняя длина свободного пробега молекул, определяемая по формуле .
, 16 I п ц
'-¡жТЪрЛёг'
Здесь б и а — коэффициенты аккомодации, у — отношение удельных теплоемкостей совершенного газа (показатель адиабаты).
В работе [15] сделана попытка усовершенствования условий скольжения в переходной по числу Кнудсена области. Основная идея состоит в определении коэффициента пропорциональности между скоростью и ее градиентом из условия согласованности со свободномолекулярным режимом. В [15] показано, что для этого необходимо использовать условие скольжения в виде
2 5л ди и = - —Хс—.
9 16 5 ду
В настоящей работе данное условие (называемое обобщенным) вместе с Максвелловскими условиями, а также условиями прилипания использовано при проведении параметрических исследований гиперзвукового обтекания плоской пластины. 'ь
Для получения безразмерных коэффициентов трения С/ и теплового потока 0 использованы соотношения
ди п( 8и . дТ V
г с 2Г*+Х*)
С'-7Л’ а—& '
Выражение для £? состоит из двух слагаемых. Первое слагаемое связано с вязкой диссипацией, а второе — с теплопроводностью. Следует отметить, что первое из слагаемых в определении 0 равно нулю для условий прилипания и отлично от нуля в случае использования условия скольжения на твердой поверхности тела.
Уравнения Навье — Стокса решались полностью неявным конечно-разностным методом [6] с использованием интегроинтерполяцион-ного метода к аппроксимации. Для вычисления конвективной составляющей вектора потока применена монотонизированная схема второго порядка точности (МШСЬ). При интерполяции консервативных зависимых переменных использован ограничитель вида
/(«,*) =
2 аЬ , п
аЬ £ 0,
а + Ъ ’ 2аЬ(а + Ь)
(а - Ь)
, аЬ < О,
а для вычисления собственных значений и собственных векторов метод Роу [16] для приближенного решения задачи Римана о распаде произвольного разрыва.
Для аппроксимации диффузионной составляющей вектора потока в уравнениях Навье — Стокса применена разностная схема типа центральных разностей. Шаблон, на котором аппроксимируются двумерные уравнения Навье — Стокса, состоит в общем случае из 13 точек.
Для решения нелинейных сеточных уравнений Р (X), где X — вектор искомых сеточных переменных, использован модифицированный метод Ньютона
^*+11 = *[*] _ %к+1 о-1 ґ(хк),
где 2) = дЖ/дХ — матрица Якоби, к — номер итерации. В процессе численного решения параметр т*+ ! определялся по формуле [17]
_ (АЛгМ - АЛГ^1!, дМ _ лК*-1!)
Х*+1 (да^З _ д*!*-1])2 ’
где ДЛГ^ — вектор поправок. По мере сходимости итерационного процесса хк -► 1. Формирование матрицы Якоби осуществлялось при помощи конечных приращений вектора невязки по вектору искомых сеточных переменных. При аппроксимации уравнений Навье — Стокса оператор дґ/дХ имеет разреженную блочную тринадцатидиагональную структуру, а элементарный блок — плотная матрица размером 4x4.
Решение системы линейных алгебраических уравнений, получаемых на итерации по нелинейности, осуществлялось при помощи прямого метода Ы1 — разложения (дЖ/дХ = Ь х и, где Ь — нижняя, и — верхняя треугольные матрицы). Для снижения суммарного числа арифметических операций предварительно анализировалась структура разреженности матриц Ь и и и проводилась перенумерация неизвестных по обобщенному методу вложенных сечений [6]. Эта методика была многократно опробована в численных экспериментах и доказала свою надежность и высокую эффективность [18].
Численное решение уравнений Навье — Стокса осуществлялось с использованием ортогональной неравномерной расчетной сетки, содержащей 101 х 101 точек. Сгущение узлов осуществлялось по двум направлениям: в пристеночной области пластины толщиной
8(5 = 2/БЬ3^) число узлов определялось равным 20% от числа узлов по нормальной координате (101), в области передней кромки толщиной 5Х (5Х = 2/Яе) задавалось такое же количество узлов.
2. Постановка и решение задачи методом прямого статистического моделирования Монте-Карло. Решение задачи при кинетическом подходе проводилось методом прямого статистического моделирования, подробное описание которого дано в [5]. Здесь отметим только, что в этом методе реальное течение потока разреженного газа моделируется движением ансамбля частиц в некоторой расчетной области. Расчетная область разбивается на ячейки, размер которых должен быть меньше местной длины свободного пробега частиц. В начальный момент времени область течения заполняется частицами, поступательные скорости и внутренняя энергия которых определяются по начальной функции распределения, как правило, соответствующей невозмущенному состоянию газа в потоке. Затем последовательно на каждом шаге по времени Л проводятся свободное перемещение частиц и столкновение между ними, причем сталкиваться могут лишь частицы, находящиеся в одной геометрической ячейке. При движении частиц в расчетной области они могут сталкиваться с твердыми поверхностями, а также вылетать за пределы области. Вылетевшие частицы исключаются из дальнейшего рассмотрения, а на каждом временном шаге с границ области проводится вбрасывание частиц в соответствии с граничной функцией распределения. По прошествии некоторого числа шагов в системе устанавливается квазистационарное состояние, с этого момента производится сбор необходимой информации о полях течения, о потоках импульса и энергии на поверхности и других выходных параметрах.
В данной работе решение проводилось одним из вариантов метода Монте-Карло [19], в котором вероятность столкновения частиц в ячейке за время Л определяется соотношением
ад
где — полное сечение столкновения, g — относительная скорость пары частиц, IV — объем ячейки. Пара столкнувшихся частиц опреде-
ляется на основе (2.1) перебором всех возможных пар в данной ячейке. Как показала практика (см., например, [20]), этот метод позволяет проводить расчеты с малым числом частиц, вплоть до числа частиц в невозмущенном потоке Щ = 1. Это обстоятельство позволяет уменьшить требования к оперативной памяти ЭВМ и при ограничениях на объем оперативной памяти проводить расчеты с большим числом ячеек, чем в случае »1, а следовательно, при заданной оперативной памяти проводить расчеты при меньших числах Кнудсена.
Проведенные в данной работе методические исследования (Мзд = 23, Т„/Т0 = 0,2, Яе«, = 1000) показали, что при изменении N0 от
1 до 4 результаты расчетов потоков импульса и энергии не зависят от величины N0 (при одинаковых величинах потоков частиц на поверхность). При проведении систематических расчетов принималось N0 = 2; число ячеек изменялось в пределах (2,24 + 6) • 105; число моделирующих частиц изменялось в пределах (4,7 * 12,5) • 105; число ударов частиц о поверхность для большинства случаев составляло величину
Сечение столкновения определялось на основе модели сфер переменного диаметра [21], [22] для степенного потенциала взаимодействия частиц Л (г) = А/г*. Для этой модели зависимость коэффициента вязкости от температуры имеет вид ц ~ Гш, ю = 0,5 + 2/5, а связь между средней длиной свободного пробега молекул и коэффициентом вязкости описывается соотношением
В дальнейшем будет принято 5 = 10; при этом значении 5 зависимость ц. от Т хорошо аппроксимирует температурную зависимость коэффициента вязкости азота при Т £ 300К. В интервале температур Т = 300 -г- 20000 К погрешность аппроксимации справочных данных [23] не превышает 6%.
Одним из основных требований при применении метода Монте-Карло является выполнение условия И <Х3, где А — размер ячейки, Хв — местная средняя длина свободного пробега частиц. Практика применения метода показала, что результаты расчетов не зависят от размеров ячейки при А < Х5/3; кроме того, выяснилось [24], [25], что такое жесткое требование необходимо выполнять в направлении наибольших 1радиентов функции распределения или ее моментов. Поскольку в различных зонах течения величины Я., могут сильно различаться, то желательно иметь методику, позволяющую в процессе счета «подстраивать» расчетную сетку под местную длину свободного пробега частиц. В данной работе такая «подстройка» проводится на основе расчетов поля плотности (см. также [25], [26]). Действительно, средняя длина свободного пробега частиц для модели сфер переменного диа-
АГуд 2»2-105.
метра в равновесных условиях Х5 - Т2*/п. Если температура в поле течения Т > Т„, то для оценок можно принять, что
Л -Л За.
* Л3 — Л*оо ,
П
и использовать это соотношение для определения размеров ячеек Адаптация размеров ячеек под поле плотности проводилась следующим образом: в начальный момент времени базовый размер ячейки устанавливался равным й (расчеты велись на декартовой сетке) и на текущей итерации определялось поле плотности; на следующей итерации в том случае, если плотность газа в ячейке превышает значение и, (величина и, принималась равной и, = 1,5/1«,), то проводится дробление базовой ячейки в направлении той координатной оси, вдоль которой градиент плотности максимален. Если ячейка примыкает к твердой поверхности, то дробление проводится по той координатной оси, вдоль которой максимальна проекция нормали к поверхности. Такое дробление в определенном смысле обеспечивает и адаптацию сетки к обтекаемым поверхностям. Как правило, двух-трех итераций обычно бывает достаточно для получения сходящихся результатов.
Расчеты проводились для газа, имеющего вращательные степени свободы. Обмен энергией между поступательной и вращательной модами описывался моделью Ларсена — Боргнакке [27] с постоянным параметром <ро, который связан с параметром релаксации = т д/тс (тд — время вращательной релаксации, тс — среднее время свободного пробега молекул) соотношением
2,°6
Фо=_т?—• (2.2)
¿я
Для моделирования одноатомного газа параметр сро полагался равным 10-4, что практически означает «замораживание» обмена энергией между поступательными и вращательными степенями свободы. В этом случае число Маха и значение температурного фактора = Т’н'/Т’о определялись для отношения удельных теплоемкостей у = 5/3. В качестве граничных условий на поверхности принималось диффузное отражение с коэффициентом аккомодации, равным единице.' ■
В процессе расчета определялись параметры газа в поле течения: плотность и, компоненты средней скорости .0}, поступательная Тгг (и вращательная 7д) температура, а также температуры газа в направлениях координатных осей . ,
кТи/2 = ||(|, - . (2.3)
где %i,Vi — компоненты скорости и средней скорости молекул. На поверхности пластины вычислялись потоки импульса и энергии
Рп=гп\%п%пГ{\>Е11)<11(ІЕл,
Рх=т\%£п №,Ея)сії<ІЕл,
Ґ е2 Л
е = ] Щ- + Ек ^„/(ІЕ^ІаЕі
/
и «поверхностные» значения плотности газа, тангенциальной скорости, поступательной и вращательной температур. Здесь Ед — вращательная энергия молекул. Для получения безразмерных коэффициентов нормальной, тангенциальной компонент силы и теплового потока использовались скоростные напоры рю 17^/2 и р» и^/2 соответственно.
Программа расчета двумерного обтекания тел методом Монте-Карло разрабатывалась совместно с В. Д. Перминовым.
3. Сравнение результатов для различных моделей газа. В данной работе проводится подробное сопоставление результатов численных расчетов, полученных на основе кинетического подхода методом Монте-Карло и на основе уравнений Навье — Стокса при значениях числа Ке« = 10 ООО. Расчеты проводились при значении температурного фактора = 0,2 и числах М«, = 10 и 23. Для Мте = 23 рассматривался как одноатомный газ, так и газ с вращательными степенями свободы молекул. В последнем случае при решении задачи методом Монте-Карло параметр сро феноменологической модели Ларсена — Боргнакке принимался равным 0,7, что соответствует параметру релаксации Zr = 3 (см. (2.2)). Такое достаточно малое значение параметра вращательной релаксации было принято для того, чтобы обеспечить лучшее согласование с решением уравнений Навье — Стокса, в которых рассматривался совершенный газ. Соответствующие перечисленным параметрам значения числа Кнудсена приведены в таблице.
Следует отметить, что проведенное сравнение имеет более общий характер, а не только относится к конкретному значению числа в силу свойства «параболичности» решения. Под этим термином понимается то, что результаты расчетов характеристик течения для фиксированного значения числа Яе^ = Ке«, х/Ь не зависят от длины пластины за исключением некоторой зоны вблизи задней кромки (І — длина пластины). В качестве примера на рис. 1
приведены результаты расчетов методом Монте-Карло параметров на поверхности пластины для двух значений числа «4400 и 10000. Из этих данных следует, что влияние задней кромки простирается на расстояние, примерно равное 0,1£, вверх по потоку. Ранее для меньших значений числа Не«, в работе [10] отмечалась параболичность решения для сопротивления трения на пластине, но величина давления существенно зависела от размеров пластины.
м. Т Кп
23 5/3 0,00315
23 7/5 0,00292
10 5/3 0,00138
ЦЙ Лл.,,/0*
Рис. 1. Результаты расчетов методом Монте-Карло параметров на поверхности пластины для двух значений числа
/*/_ =*Ю; гш-о,2; у-5/з Обратимся теперь к прямому сопо-
о о пе^-шо...... ставлению результатов. На рис. 2
приведены результаты расчетов давления таза на поверхности пластины, полученные методом Монте-Карло (МК-результаты) и на основе уравнений Навье — Стокса (НС-результаты) с граничными условиями прилипания (N$1) — кривая 1, максвелловского скольжения (Мз1) — кривая 2 и обобщенного скольжения (СИ) — кривая 3 для Мв = 23 и у = 1,4. Анализ этих данных показывает, что лучшее согласование МК-результатов и НС-результатов имеет место для граничных условий прилипания. При х/Ь >0,5 влияние граничных условий для НС-уравнений уменьшается и МК-расчеты хорошо согласуются с НС-результатами при Мв1. Следует подчеркнуть различие в МК- и НС-результатах, обусловленное наличием внутренних степеней свободы молекул и связанное, по-видимому, с неравновесным характером течения газа, не учитываемого при решении НС-уравнений. Отметим также, что использование обобщенных граничных условий привело к уменьшению давления по всей длине пластины и особенно существенно в окрестности передней кромки. Применение этих граничных условий, как и ожидалось, улучшает согласование МК-и НС-результатов для сопротивления трения вблизи передней кромки (рис. 3). В целом при х/Е >0,1 (Мм = 23) и при х/Ь = 0,05 (Мю = 10) результаты МК-расчетов трения хорошо согласуются с данными НС-расчетов (различие не превышает 10%). Из расчетных данных следует, что влияние граничных условий для НС-уравнений уменьшается с уменьшением числа Маха. Что же касается теплового потока к поверх-
Рис. 2. Давление газа на поверхности пластины:
1 — уравнения Навье — Стокса с граничными условиями прилипания;
2 — скольжения (Мб!); 3 — скольжения (вй); 4 — метод Монте-Карло
Рис. 3. Сопротивление трения:
1 — уравнения Навье — Стокса с граничными условиями прилипания; 2 — скольжения (Ш1); 3 — скольжения (Єві); 4 — метод Монте-Карло
ности пластины, то для этой характеристики лучшее согласование МК-и НС-результатов имеет место для обобщенных граничных условий скольжения (рис. 4). Влияние граничных условий на решение НС-уравнений при х/ь >0,1 не превышает 10%, причем в противоположность давлению и сопротивлению трения; тепловой поток к поверхности на большей части пластины при (М выше, чем при других граничных условиях. Аналогичные результаты получены и для случаев у = 5/3, Мм = 23 и М.*, =10, причем с уменьшением числа М согласование МК- и НС-результатов происходит при меньших числах Ке^.
Следующая серия сравнений результатов относится к характеристикам: плотности, температуре и тангенциальной скорости газа, определенным на поверхности пластины. На рис. 5 приведены результаты расчета плотности, показывающие хорошее согласование МК- и НС-
Рис. 4. Тепловой поток:
1 — уравнения Нанье — Стокса с граничными условиями прилипания; 2— скольжения (М&1); 3 — скольжения (СМ); 4 — метод Монте-Карло
0,2 0,4 0,6 0,8
х
Рис. S. Плотность газа на поверхности пластины:
1 — уравнения Навье — Стокса с граничными условиями прилипания; 2 — скольжения (Msl); 3 — скольжения (Gsl); 4 — метод Монте-Карло
расчетов с максвелловскими граничными условиями. Влияние граничных условий наиболее сильно проявляется при переходе от условий прилипания к условиям со скольжением и температурным скачком, чем при изменении условий скольжения. Расчеты, проведенные при других значениях температурного фактора - 0,081 — 1,0, показали достаточно сильное влияние граничных условий при <0,2, а при ^ >0,5 это влияние существенно уменьшается за исключением зоны течения вблизи передней кромки. На рис. 6 показаны результаты расчетов тангенциальной скорости газа на поверхности, из которых видно согласование МК- и НС-расчетов с Мя1 и сильная зависимость от граничных условий при решении НС-уравнений.
Наиболее сильное различие МК- и НС-результатов проявляется при расчете температуры газа на поверхности (рис. 7). При больших значениях числа М«, (М*, = 23) область достаточно сильного различия
и
y~7/S;M„=2J
0°СО ООО ОО ООО ооб ойоовв"
J------1-----11 ..I______' 1
_________I /— 1
О 0,2
0,4 0,6 Oj х
Рис. 6. Скорость газа на поверхности пластины:
1 — уравнения Навье — Стокса с граничными условиями прилипания; 2 — скольжения (Мй); 3 — скольжения (вй); 4 — метод Монте-Карло
Рис. 7. Температура газа на поверхности пластины:
1 — уравнения Навье — Стокса с граничными условиями прилипания; 2 — скольжения (Мй); 3 — скольжения (СМ); 4 — метод Монте-Карло
результатов простирается почти до середины пластины, причем учет вращательных степеней свободы молекул приводит к более сильному различию результатов. С уменьшением числа Мв различия уменьшаются.
При сопоставлении температуры и скорости газа на поверхности необходимо иметь в вицу, что эти величины определены для уравнений Навье — Стокса как некоторые фиктивные значения на стенке, обеспечивающие вне кнудсеновского слоя совпадение решений НС-уравнений с решением уравнения Больцмана с истинными условиями на границе. Таким образом, на рис. 6, 7 сопоставляются «истинные» значения параметров газа, полученные в результате решения уравнения Больцмана методом Монте-Карло, и некоторые «фиктивные» значения, определяемые из решения сплошносредных уравнений, и в общем случае могут и не совпадать. Тем не менее, как это следует из рис. 6, 7, по мере удаления от передней кромки имеется хорошее количественное согласование параметров газа на поверхности, полученных из МК-и НС-решений. Общая тенденция такова, что лучшее согласование происходит при меньших градиентах газодинамических величин в поле течения и при меньшей степени неравновесности газа, т. е. с уменьшением числа Мдз и по мере удаления от передней кромки.
Наконец, обратимся к сопоставлению результатов расчета полей течения. Результаты расчета профилей плотности р-п(у)/п0 приведены на рис. 8, 9 для у = 5/3. Анализ приведенных данных показывает, что при М«, = 23 вблизи передней кромки МК-скачок более размыт, а положение и максимальное значение р примерно соответствуют навье-стоксовскому скачку при граничных условиях обобщенного скольжения (СМ). По мере удаления от передней кромки толщина и амплитуда МК-скачка приближаются к НС-скачку, а положение максимума смещается к положению НС-скачка с граничными условиями Мв1. Учет внутренней структуры молекул не вносит качественно новых эффектов. При М«, = 10 все характеристики монтекарловского скачка плотности близки к характеристикам навьестоксовского скачка. Результаты расчета поступательной температуры в поле течения, приведенные на рис. 10, показывают наиболее сильные отличия для разных моделей газа вблизи поверхности пластины, что соответствует приведенным выше данным на самой поверхности. По мере отхода от поверхности МК-результаты приближаются к Н С -результатам и располагаются между теми из них, которые получены с граничными условиями Мв1 и СМ
Рис. 8. Профили плотности газа:
1 — уравнения Навье — Стокса с граничными условиями прилипания; 2 — скольжения (Мй); 3 — скольжения (051); 4 — метод Монте-Карло
Рис. 9- Профили плотности газа:
1 — уравнения Навье — СГоксас граничными условиями прилипания; 2 — скольжения (Ш1); 3 — скольжения (СЗД; 4 — метод Монте-Карло
Рис. 10. Профили поступательной температуры:
1 — уравнения Навье — Стокса с граничными условиями пряли- '
пания; 2 — скольжения (М^); 3 — скольжения (вя!); 4 — метод Монте-Карло
С уменьшением числа Мот МК-профили температур приближаются к НС-профилям, полученным для граничных условий Мв1. В отличие от других характеристик профили компоненты скорости их (рис. 11), полученные в МК-расчетах, при Мм = 23 лучше всего согласуются с НС-результатами с граничными условиями обобщенного скольжения (Ой). При М«, = 10 влияние граничных условий уменьшается и МК-результаты располагаются между НС-результатами для разных’ условий скольжения.
4. О границе применимости уравнений Навье — Стокса для описания течений разреженного газа. Выше отмечалось, что применение уравнений Навье — Стокса для описания течений разреженного газа требует обоснования в каждом конкретном случае, поскольку область их при-
Рис. 11. Профили продольной компоненты скорости:
1 — уравнения Навье — Стокса с граничными условиями прилипания; 2 — скольжения (Msl); 3 — скольжения (Gsl); 4 — метод Монте-Карло
менимости определяется асимптотически при Кп -> 0 и предполагается малое отклонение функции распределения от локально равновесной. Приведенные выше результаты расчетов показывают, что в рассматриваемой задаче имеются области течения вблизи поверхности и в зоне ударной волны, где состояние газа является существенно неравновесным и где применение уравнений Навье — Стокса является недостаточно точным. В других областях течения, где градиенты переменных невелики, имеется хорошее согласование результатов для различных моделей описания течения газа.
Если же говорить о распределении силовых и тепловых характеристик по пластине, то здесь положение таково: МК- и НС-результаты хорошо согласуются между собой на большей части пластины за исключением области вблизи передней кромки. Границы применимости уравнений Навье — Стокса оказываются различными для различных областей и различных параметров течения. Если принять в качестве критерия примерно 10% различия в результатах расчета указанных характеристик, то суммарно можно сказать, что согласование молекулярной и сплошносредной моделей газа начинается при x/L = 0,1 для М* = 23 и при x/L = 0,05 (Мда = 10). Это соответствует локальным числам Re^ = 1000 (Mw = 23) и 500 (Мда = 10) или Кпх = 0,03.
Приведем также «граничное» значение критерия подобия Re0^, поскольку применение этого критерия позволяет скоррелировать результаты для разных моделей взаимодействия между молекулами [10], [20]: Мте = 23 — Re0jc = 30, Мм = 10 — Re0x = 45. Эти значения отличаются от полученных ранее в работе [10] оценок, а именно: Re0x > 10, но в этой работе оценка носила качественный характер без привлечения количественного критерия. Приведенные оценки для Re0x согласуются с данными, полученными в работе [15].
Полученные данные позволяют с большей степенью достоверности судить о возможностях применения уравнений Навье — Стокса для описания течений разреженного газа. Анализ полей течения, полученных методом Монте-Карло, показывает, что имеются области с нерав-
новесным характером течения (вблизи передней кромки, вблизи поверхности, в зоне ударной волны), о чем свидетельствует анизотропия поступательной температуры в различных направлениях и различия между поступательной и вращательной температурами для газа с внутренними степенями свободы. В этой связи представляет интерес выявление возможной анизотропии при решении уравнений Навье — Стокса (например, определение компонент тензора напряжений в поле течения).
Результаты настоящей работы показывают достаточно сильное влияние граничных условий на поверхности пластины на распределение газодинамических переменных для уравнений Навье — Стокса. Сопоставление скорости и температуры газа на поверхности, полученных при решении уравнений Навье — Стокса и методом Монте-Карло, показывает, что наибольшее отличие имеет Место для температуры газа на поверхности. Выяснение причин этого расхождения требует, на наш взгляд, дальнейшей теоретической проработки проблемы граничных условий.
ЛИТЕРАТУРА
1. В ird G. A. Low density aerothermodynamics // AIAA Paper 85-0994.
2. Ivanov M. S. Rarefied CFD for reentry problems. Proc. 5-th ISCFD.
Sendai.- 1993. V. II. !
3. Mundt C., Monnoyer F., Hold R. Computational simulation of aerothermodynamic characteristics for the reentry of HERMES//AIAA Paper.—
1993, N 93-5069. : ^
4. Yamamoto Y. Numerical simulation of real gas effects and aerodynamic heating of hypersonic space transportation vehicles'// Computational Methods in Applied Sciences.— 1992.
5. Берд Г. Молекулярная газовая динамика,— М.: Мир.— 1981.
6. Егоров И. В.( Зайцев О. JI. Об одном подходе к численному решению двумерных уравнений Навье — Стокса методом сквозного счета //
Ж. вычисл. матем. и матем. физ,— 1991. Т. 31, № 2.
7. Коган М. Н. Динамика разреженного газа.— М.: Наука.— 1967.
8. Кузнецов М. М., Никольский В. С. Кинетический анализ шперзвуковых вязких течений многоатомного газа в тонком трехмерном ударном слое // Ученые записки ЦАГИ.— 1985. T. XVI, № 3.
9. Cheng Н. К., Wong Е. Y., Dogra V. К. A shock-layer theory based on thirteen-moment equations and DSMC calculations of rarefied hypersonic flows // AIAA Paper 91-0783.- 1991.
10. Гусев В. H., Ерофеев А. И., Климова Е. В., Перепу-хов В. А., Рябов В. В., Толстых А. И. Теоретические и экспериментальные исследования обтекания теп простой формы шперзвуковым потоком разреженного газа // Труды ЦАГИ,— 1977. Вып. 1855.
11. Ботин А. В., Егоров И. В., Николаев К. В. Сравнение сплошносредных и кинетических методов расчета шперзвукового обтекания затупленного тела в промежуточном режиме // Труды X юбилейной научнотехнической конференции по аэродинамике больших скоростей,— 1990.
12. Lengrand J. С. Using experimental results to validate numerical solutions of rarefied flow problems // AIAA Paper 94-2632.— 1994.
13. Pullin D. I., Harvey J. K. A numerical simulation of the rarefied hypersonic flat-plate problem // J. Fluid Mech.— 1976. V. 78, pt 4.
14. Nagamatsu H. T., Messitt D. G., Myrabo L. N., Sheer R. E. Computational, theoretical and experimental investigation of flow over a sharp flat-plate, Mj = 10 — 25 // AIAA Paper 94-2350.— 1994.
15. Gokcen T., MacCormack R. W., Chapman D. R. Computational fluid dynamics near the continuum limit // AIAA 8-th Computational Conference.— 1987, AIAA Paper 87-1115.— 1987.
16. Roe P. L. Approximate Riemann solvers, parameter vectors and difference scheme // J. Comput. Phys.— 1981. V. 43.
17. Каримов Т. X. О некоторых итерационных методах решения нелинейных уравнений в гильбертовом пространстве // ДАН СССР.— 1983. Т. 269, № 5.
18. Yegorov I. V., Zaitsev О. L. Development of efficient algorithms for computational fluid dynamic problems // Proc. 5-th ISCFD. Sendai.— 1993. V. III.
19. Белоцерковский О. М., Яницкий В. Е. Статистический метод частиц в ячейках для решении задач динамики разреженного газа // Ж. вычисл. матем. и матем. физ,— 1975. Т. 15, № 5, б.
20. Белоцерковский О. М., Ерофеев А. И., Яницкий В. Е. О нестационарном методе прямого статистического моделирования течений разреженного газа // Ж. вычисл. матем. и матем. физ.— 1980. Т. 20, № 5.
21. Ерофеев А. И. О моделировании межмолекулярного взаимодействия при решении уравнения Больцмана методом Монте-Карло // Изв. АН СССР, МЖГ,- 1977, № 6.
22. В i г d G. A. Monte-Carlo simulation in an engineering context // AIAA Progress in Astronautics and Aeronautics: Rarefied Gas Dynamics. V. 74, pt. 1, AIAA. New York.- 1981.
23. Гордеев О. А. Потенциалы взаимодействия, упругие сечения, интегралы столкновений компонентов воздуха для температур до 20 ОООК (Методы определения, рекомендуемые данные). Обзоры по теплофизическим свойствам веществ № 5 (55).- М.: ИВТ АН СССР.- 1985.
24. Николаев К. В. Аэродинамические и тепловые характеристики обтекания затупленных тел разреженным газом. Канд. диссертация.— М.: МФТИ.- 1990.
25. Taylor J. С., Moss J. N., Hassan Н. A. Study of hypeisonic flow past sharp cones // AIAA.— 1989, N 1713.
26. Olynic D. P., Hassan H. A., Moss J. N. Grid generation and adaptation for direct simulation Monte-Cado method // AIAA Thermophysics, plasmadynamics and lasers conference, June 27—29 1988, San Antonio, Texas.
27. Larsen P. S., Borgnakke C. Statistical collision model for simulating polyatomic gas with restricted energy exchange.— Proc. 9-th. Intern. Symp. RGD, DFVLR Press, Potz-Wahn. V. 1.- 1974.
Рукопись поступила 9/11996
