Численное моделирование сверхзвукового течения газа около плоской пластины на основе кинетических и континуальных моделей Текст научной статьи по специальности «Математика»

Вычислительные технологии

Том 16, № 6, 2011

Численное моделирование сверхзвукового течения газа около плоской пластины на основе

•J» >к

кинетических и континуальных моделей*

А. А. Шершнёв, А. Н. Кудрявцев, Е. А. Бондарь Институт теоретической и прикладной механики им. С.А. Христиаповича СО РАН, Новосибирск, Россия e-mail: [email protected], [email protected], [email protected]

Рассмотрено сверх- и гиперзвуковое обтекание конечной плоской пластины под нулевым углом атаки разреженным потоком одноатомного газа. Численное моделирование задачи проведено с помощью метода прямого статистического моделирования, а также решения модельных кинетических уравнений релаксационного типа и уравнений Навье — Стокса с граничными условиями, учитывающими скольжение и скачок температуры на поверхности. Исследованы границы применимости континуального подхода и подхода, основанного на решении модельных кинетических уравнений. Показано, что сильная неравновесность течения при высоких гиперзвуковых скоростях приводит к тому, что уравнения Навье — Стокса не позволяют даже качественно верно предсказать поведение основных газодинамических параметров вблизи пластины. В то же время подход, основанный на решении кинетических уравнений релаксационного типа, демонстрирует хорошее согласие с данными прямого статистического моделирования при всех исследованных условиях.

Ключевые слова: кинетические модели, конечно-разностные схемы, классическая задача.

Введение

Хорошо известно, что континуальный подход к описанию движения жидкостей и газов теряет применимость, когда средняя длина свободного пробега молекул Л становится сравнимой с характерным размером задачи Ь, Обычно считается, что уравнения Навье — Стокса с граничными условиями прилипания могут быть использованы до чисел Кнудсена Кп = Л/Ь > 0.01, а с условиями, учитывающими скольжение и скачок температуры на поверхности обтекаемого тела, — до Кп = Л/Ь > 0.1. В случае тела с острой передней кромкой характерным размером течения вблизи кромки фактически является радиус ее закругления. Это означает, что эффекты разреженности могут играть существенную роль на передней части тела даже в случаях, когда число Кнудсена,

*Работа выполнена при финансовой поддержке РФФИ (грант № 10-08-01203-а), Программы Президиума РАН "Фундаментальные проблемы механики взаимодействий в технических и природных системах, материалах и средах" и Междисциплинарного интеграционного проекта СО РАН № 40 (20092011).

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

Численное моделирование течений разреженного газа должно основываться на моле-кулярно-кинетичееком подходе, т, е, на решении уравнения Больцмана для функции распределения молекул, В настоящее время наиболее распространенным численным методом решения данного уравнения является метод прямого статистического моделирования (ПСМ) [1]. Однако присущий методу ПСМ статистический разброс данных снижает его эффективность при расчете нестационарных и/или низкоскоростных течений, Альтернативный подход, основанный на детерминистическом решении уравнения Больцмана, из-за необходимости вычисления сложных многомерных интегралов при определении столкповительпого члена требует очень больших вычислительных ресурсов, Фактически это делает очень дорогостоящим его применение при решении практических задач даже с помощью мощных современных многопроцессорных ЭВМ,

Вычислительные затраты можно существенно уменьшить, если вместо уравнения Больцмана использовать модельное кинетическое уравнение, в котором интеграл столкновений заменен на простой член релаксационного типа. Такой подход в наиболее систематическом виде представлен в монографии Е.М, Шахова [2]. Были предложены различные аппроксимации интегрального члена, наиболее известными являются модель Бхатнагара — Гросса — Крука (БГК) [3] и ее модификации: эллипсоидальная статистическая (ЭС) модель [4] и модель Шахова [5].

Несмотря на то что данные модели применялись для исследования широкого набора задач динамики разреженного газа, априори их точность для сильно неравновесных течений не очевидна. Основной недостаток этих модельных уравнений состоит в том, что в них в отличие от исходного уравнения Больцмана вероятность столкновения между двумя молекулами не зависит от относительной скорости последних.

Целью настоящей статьи является исследование точности модельных кинетических уравнений релаксационного типа на примере сильно неравновесного течения — еверх-и гиперзвукового обтекания конечной плоской пластины с острой передней и задней кромками.

Обтекание конечной плоской пластины потоком разреженного газа под нулевым углом атаки — одна из классических задач динамики разреженного газа. При различных параметрах набегающего потока и температурах пластины эта задача изучалась многими исследователями, использовавшими различные подходы, В частности, задача моделировалась на основе уравнения БГК [6], эллипсоидальной статистики (ЭС) модели [7] и модели Шахова [8]. Она решалась также с помощью прямого статистического моделирования [9]. Из наиболее поздних работ следует отметить [10-12], в которых сравниваются результаты, полученные при помощи континуального (уравнения Навье —Стокса и Барнетта) и кинетического (метод ПСМ) подходов. Предпринимались также попытки приближенного асимптотического анализа проблемы [13]. Несмотря на такое обилие работ, в литературе отсутствует достаточно подробное сравнение решений, полученных при одних и тех же условиях из различных модельных кинетических уравнений, уравнений Навье —Стокса и метода ПСМ, который в данном случае может быть принят за эталонный. Кроме того, многие исследования, где течение над плоской пластиной моделировалось с использованием модельных кинетических уравнений, были выполнены достаточно давно, когда возможности ЭВМ еще не позволяли получить решение с высокой точностью, В [14] данная задача рассматривалась только с помощью уравне-

ния БГК, к тому же основное внимание в ней было уделено возникновению разрывов функции распределения на передней и задней кромках,

В представленной работе различные модельные уравнения релаксационного типа используются для моделирования обтекания бесконечно тонкой пластины единичной длины, расположенной вдоль сверхзвукового потока. Все расчеты проведены для газа из твердых сфер (вязкость ^ ~ T0'5), Число Кнудеена, вычисленное по длине пластины, Кп = 0.01 а число Маха набегающего потока М = 2 и 10, Температура пластины принята равной температуре набегающего потока — это приблизительно соответствует условиям реального полета, когда поверхность летательного аппарата активно охлаждается, Результаты моделирования сравниваются с данными ПСМ и с решениями уравнений Навье — Стокса,

1. Основные уравнения

1.1. Модельные кинетические уравнения

Двумерное кинетическое уравнение релаксационного типа можно записать в следующем виде:

Здесь f (t, x, £) — функция распределения молекул по скоростям, x = (x,y) — вектор пространственных координат, £ = ) = (£ъ£2,£3) _ вектор молекулярных ско-

ростей, v — частота столкновений, fN — равновесная функция, вид которой зависит от выбранной модели.

Макроскопические величины, такие как плотность р, средняя скор ость u, температура T и др., могут быть найдены из соответствующих моментов функции распределения:

+оо +оо +оо

Р = UI fdíx diy , pu = jjJ £fdíx diy , 3pRT = jjJ c2fdix diy d^z,

—о — о —о

где c = £ — u, R — газовая постоянная, 1.1.1. Модели

fN

fN = fM = _Р_

1 1 (2тг ДТ)3/2 Р

1 3

2ВТ ^

i,j=i

Основной недостаток этой модели — неправильное значение числа Прандтля, Это шн-чит, что разложение Чепмена —Энскога уравнения БГК приводит к уравнениям Навье — Стокса с числом Прандтля Рг = 1, в то время как правильное значение для одноатомного газа Рг = 2/3, Желание исправить данный недостаток привело к построению двух следующих моделей.

Модель Шахова. В модели Шахова /м определяется следующим образом:

1 , „ , /с2 "

N

1 +-(1 — Рг) с • q--5

bpRT J \RT

f

Здесь q — вектор теплового потока, p = pRT — давление

/м ч= ///

Эллипсоидальная статистическая модель. В ЭС-модели fм задается анизотропным гауссианом

f

N _ Р

п

3/2

лАыа

где

4 - (ШТ * ~

11 V Рг рРг

ехр

-1

г,] = 1

г] СгС]

+оо

Рг]

Сг С] fd£x d£у ^

В уравнении БГК и модели Шахова частота столкновений вычисляется по формуле V = р/^, где ^ — вязкость газа, в то время как в ЭС-модели V = Ргр/^.

3

1.1.2. Укороченные функции распределения

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

Функций д и К достаточно для вычисления всех основных макроскопических газодинамических переменных.

Уравнение (1) тогда может быть заменено системой из двух уравнений для укороченных функций распределения

01 Ч-Аг Ч-У// '

которая решается численно.

1.1.3. Численный метод

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

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

предпочтительно применять метод высокого порядка точности, В настоящей работе использовалась одна из современных схем сквозного счета — WENO-ехема пятого порядка точности [15], успешно примененная ранее при решении многих задач вычислительной аэродинамики [16],

На твердой стенке ставились условия диффузного отражения, в соответствии с которыми функция распределения частиц, отраженных от стенки, является максвелловской с температурой, равной температуре стенки, а их плотность определяется из условия равенства потоков падающих и отраженных молекул.

Вычислительный код был распараллелен путем декомпозиции расчетной области в физическом пространстве на подобласти, каждая из которых считается на отдельном процессоре. Для пересылки сообщений между процессорами используется библиотека MPI.

1.2. Уравнения Навье — Стокса

При решении сжимаемых уравнений Навье — Стокса для вычисления конвективных потоков применена упомянутая WENO-ехема, диффузионные члены аппроксимируются центральными разностями четвертого порядка.

Для учета эффектов разреженности на твердой стенке ставятся граничные условия скольжения и скачка температуры [17]

Здесь индекс s означает газодинамические параметры вблизи твердой стенки, u — компонента скорости, касательная к стенке, Tw — температура стенки. Граничные условия (3) выводятся из приближенного решения кинетического уравнения в так называемом кнудееновеком слое, т, е, в узком слое вблизи поверхности, размеры которого сравнимы с локальной длиной свободного пробега молекул. Это приводит к следующим значениям коэффициентов: au = 1.142 и aT = 0.5865 [17],

1.3. Прямое статистическое моделирование

Все расчеты методом ПСМ проведены с помощью программного комплекса SMILE [18], основанного на схеме мажорантной частоты. Применяются две независимых сетки: одна для организации столкновений частиц, другая для выборки макроскопических величин. Обе сетки состоят из однородных прямоугольных ячеек, которые при необходимости разбиваются на более мелкие. На твердой стенке накладываются условия диффузного отражения,

2. Обсуждение результатов

Расчеты на основе уравнений Навье —Стокса и модельных кинетических уравнений проводились на пространственной сетке, содержащей 480 х 240 точек. Сетка в пространстве скоростей включала от 33 х 33 до 127 х 127 (в расчетах при = 10) точек. Таким образом, максимальный размер сетки в фазовом пространстве составлял примерно 1.86 х 109 узлов. Расчеты проводились па вычислительном кластере ИТПМ СО РАН, максимальное число использованных процессоров — 80,

(3)

В расчетах методом ПСМ использовано от 600 тыс, до 1 млн модельных частиц. Размер столкновительных ячеек и временной шаг выбирались заведомо меньшими локальной средней длины свободного пробега и локального среднего времени между столкновениями.

Поле числа Маха, полученное при Mœ = 2 с помощью ЭС-модели, представлено на рис, 1, Хорошо видны такие особенности течения как головная ударная волна, толстый пограничный слой, след за пластиной и веер воли разрежения за задней кромкой пластины. На рис, 2 приведены распределения плотности и давления в сечении X = 0.8, Модель Шахова (SBGK) и ЭС-модель (ESBGK) дают практически неразличимые профили, близкие к данным ПСМ (DSMC), Решения уравнений Навье —Стокса (XS) в целом хорошо согласуются с результатами ПСМ за исключением области вблизи поверхности пластины, в которой заметны значительные расхождения профилей. Результаты, полученные с помощью уравнения БГК (BGK), существенно отличаются от остальных, что объясняется неправильным значением числа Прандтля в данной модели.

Распределения поверхностных характеристик — коэффициента давления Cp — -, —, коэффициента трения С/ = -, — и коэффициента теплопередачи

Рис. 1. Течение около плоской пластины; поле числа Маха; Мо = 2

Рис. 2. Профили макровеличин в сечении X = 0.8; = 2: a, — плотность, б — давление

2

2

С\ = -у—-- вдоль пластины представлены па рис, 3, Как видно из рисунка, все три

~^Рооиоо

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

Расчеты были проведены также при числе Маха = 10, В случае гиперзвукового течения над пластиной эффекты разреженности проявляются значительно сильнее. Уравнение БГК в данных вычислениях пе использовалось в силу плохих результатов, полученных уже в предыдущих расчетах.

Что касается модели Шахова и ЭС-модели, то и в гиперзвуковом случае они демонстрируют достаточно хорошее согласие с данными ПСМ, в то время как уравнения Навье —Стокса приводят к результатам, даже качественно отличающимся от кинетических моделей.

Из рис. 4 можно видеть, что в расчетах, выполненных па основе континуального подхода, ударная волна находится дальше от пластины, имеет меньшую толщину, а максимальные значения плотности и давления, достигаемые за ударной волной, значительно выше этих показателей, полученных при кинетическом подходе. Кроме того, внутри пограничного слоя уравнения Навье —Стокса предсказывают совершенно иное

Рис. 3. Поверхностные характеристики на плоской пластине; = 2: а, — коэффициент давления Ор, б — коэффициент трения О/, в — коэффициент теплопередачи Оь

Рис. 4. Профили макровеличин в сечении X = 0.8; = 10: а — плотность, б — давление

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

Аналогичная ситуация наблюдается и с распределениями поверхностных коэффициентов (рис. 5). Решения модельных кинетических уравнений находятся в приемлемом согласии с результатами ПСМ вдоль всей длины пластины. Вместе с тем профили Навье—Стокса демонстрируют большое отклонение от данных, полученных при помощи кинетического подхода, хотя это отклонение постепенно и уменьшается с удалением от передней кромки.

Такое значительное расхождение континуального и кинетического подходов объясняется сильным нарушением равновесия в разреженном гинерзвуковом течении. На рис. 6 представлены поля отношения продольной температуры Тх к полной температуре Т, вычисленные с помощью ЭС-модели, Данная величина характеризует степень неравиовесиоети между различными поступательными степенями свободы молекул. Видно, что в случае умеренно сверхзвукового числа Маха течение неравновесно только в небольшой окрестности передней кромки, однако в гинерзвуковом потоке отсутствие равновесия наблюдается в обширной области над пластиной.

На рис. 7 показана скорость скольжения на поверхности пластины при Мо = 2 и Мо = 10, вычисленная с помощью всех трех подходов. При Мо = 2 ЭС-модель и метод ПСМ дают очень близкие результаты, в то время как скорость скольжения и в решении уравнений Навье — Стокса, вычисленная но формулам (3), существенно отличается. Следует, однако, учесть что скорость скольжения, фигурирующая в (3), — фиктивная

а б в

Рис. 5. Поверхностные характеристики на плоской пластине; Мо = 10 -.а — коэффициент давления Ср, б — коэффициент трения С/, в — коэффициент теплопередачи С^

0.2 0.15 0.1 0.05 0

0.2

0.4

X/L

0.6

0.8

Рис. 6. Отношение продольной температуры Tx к полной температуре T в расчетах на основе ЭС-модели; a — = 2, б — = 10

■I- NS, slip velocity

it- NS, gas velocity

ESBGK

0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1

I119(1(1*

0.2 0.4 0.6 0.8 л/L

Рис. 7. Скорость скольжения на пластине: а, — = 2, б — = 10

величина, определенная таким образом, чтобы решение уравнений Навье —Стокса совпадало с кинетическим решением вне кпудееповекого слоя, а не реальная скорость газа на стенке. Тем не менее, зная скорость скольжения и3, можно вычислить истинную скорость газа и9, используя приближенное решение уравнения Больцмана в кнудсеновском слое. Выкладки приводят к формуле |19|

Ug

V2Ä

1.1466

0.696 us.

(4)

Найденная с помощью этой формулы величина также приведена па рис. 7 и обозначена па нем как "XS, gas velocity". Видно, что она хорошо совпадает с данными ПСМ и ЭС-модели. Хорошее согласие между результатами ПСМ и решениями модельных кинетических уравнений сохраняется и при = 10. Однако здесь наблюдается значительное отличие этих результатов от данных, полученных из уравнений Навье — Стокса, и хотя

0

9М,€у)/Р(Х*)

0.3

0.3

0.2

0.2

0.1

0.1

о

о

а

б

Рис. 8. Редуцированная функция распределения д(х, £х, )/р(Х*) на поверхности пластины в точке X* ^ 0.2, полученная из расчета на основе ЭС-модели: а — М^ = 2, б — М^ = 10

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

Причина, по которой вблизи стенки отсутствует равновесие между поступательными степенями свободы, становится очевидной, если рассмотреть функцию распределения молекул по скоростям на поверхности пластины (рис. 8). Из рисунка можно видеть, что данная функция весьма далека от максвелловской: она имеет типичный "двугорбый" характер и разрывна при = 0. Это объясняется тем, что фактически данная функция является суммой двух отдельных функций распределения, соответствующих падающим (£у < 0) и отраженным (£у > 0) молекулам. Последние связаны между собой только условием нулевого потока молекул через стенку, причем функция распределения отраженных молекул является "половиной" максвеллиана с шириной распределения, определяемой температурой стенки. Следует отметить, что в случае М^ = 2 разница между максимумами функций распределения падающих и отраженных молекул составляет приблизительно 25%, в то время как при М^ = 10 они отличаются более чем в 10 раз. Кроме того, если сравнить ширину распределений вдоль оси £х, которая по сути определяет температуру Тх, то можно заметить, что в случае М^ = 2 ширина распределений падающих и отраженных частиц примерно одинакова, в то время как при М^ = 10 ширина распределения падающих частиц существенно больше. Все это указывает на гораздо более сильную неравновесность течения при М^ = 10.

Заключение

В работе выполнено численное моделирование разреженного течения одноатомного газа около конечной плоской пластины на основе трех различных подходов: уравнений

Навье —Стокса с граничными условиями скольжения и скачка температуры, прямого статистического моделирования и прямого конечно-разностного решения модельных кинетических уравнений релаксационного типа с помощью схем сквозного счета высокого порядка точности. Вычисления проведены для чисел Маха = 2 и = 10 и числа Кнудеена Kn = 0.01, Показано, что при = 2 все три подхода дают близкие результаты при условии, что модельное уравнение приводит к правильному значению числа Прандгля в континуальном пределе. Однако в случае = 10 в обширной области над пластиной отсутствует равновесие между поступательными степенями свободы. Это ведет к тому, что уравнения Навье —Стокса не позволяют корректно предсказать характеристики течения даже при использовании на поверхности пластины условий скольжения, В то же время решения модельных кинетических уравнений достаточно хорошо согласуются с данными ПСМ как при умеренно сверхзвуковых, так и при высоких гиперзвуковых скоростях. Таким образом, модельные кинетические уравнения релаксационного типа могут успешно применяться для численного моделирования неравновесных разреженных течений, в частности таких как гиперзвуковое обтекание тел с острой передней кромкой.

Авторы выражают благодарность проф. М.С, Иванову за полезные обсуждения полученных результатов.

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

fl] Bird G.A. Molecular Gas Dynamics and the Direct Simulation of Gas Flows. Oxford: Clarendon Press, 1994.

[2] Шахов E.M. Метод исследования движений разреженного газа. М.: Наука, 1974.

[3] Bhatnagar P.L., Gross Е.Р., Krook М. A model for collision processes in gases. I. Small amplitude processes in charged and neutral one-component systems // Phvs. Rev. 1954. Vol. 94. P. 511-525.

[4] Holway L.H. New statistical models for kinetic theory: Methods of construction // Phvs. Fluids. 1966. Vol. 9. P. 1658-1673.

[5] Шахов E.M. Об обобщении релаксационного кинетического уравнения Крука // Изв. АН СССР. МЖГ. 1968. № 5. С. 142-145.

[6] Huang А.В., Hartley D.L., Kinetic theory of the sharp leading edge problem in supersonic flow // Phvs. Fluids. 1969. Vol. 12. P. 96-108.

[7] Huang А.В., Hwangs P.F. Supersonic leading edge problem according to the ellipsoidal model // Ibid. 1970. Vol. 13. P. 309-317.

[8] Шахов E.M. Продольное обтекание пластины потоком разреженного газа. // Изв. АН СССР. МЖГ. 1973. № 2. С. 119-126.

[91 Vogenitz F.W., Broadwell J.е., Bird G.A. Leading edge flow by Monte Carlo direct simulation technique // AIAA J. 1970. Vol. 8, No. 3. P. 504-510.

[10] Elizarova T.G., Graur I.A., Chpoun A., Lengrand J.C. Comparition of continuum and molecular approaches for rarefied gas flows // Rarefied Gas Dynamics. 1995. Vol. 2. P. 780-786.

[11] Егоров И.В., Ерофеев А.И. Сопоставление моделирования гиперзвукового обтекания плоской пластины на основе метода Монте-Карло и уравнений Навье — Стокса // Изв. РАН. МЖГ. 1997. № 1. С. 135-144.

[12] Eropeev A.I., Friedlander O.G. Macroscopic models for non-equilibrium flows of monatomic gas and normal solutions // Proc. of 25th Int. Svmp. on RGD / Eds. M.S. Ivanov and A.K. Rebrov. Novosibirsk: Publ. House of the Siberian Branch of the Russian Academy Sciences, 2007. P. 117-124.

[13] Кузнецов M.M., Липатов И.И., Никольский B.C. Асимптотический анализ эффектов поступательной неравновесности в гиперзвуковом течении около плоской поверхности с острой передней кромкой // Письма в ЖТФ. 2008. Т. 34, вып. 8. С. 21-28.

[14] Aoki К., Kanba К., Таката S. Numerical analysis of a supersonic rarefied gas flow past a flat plate // Phvs. Fluids. 1997. Vol. 9. P. 1144-1161.

[15] Jiang G.-S., Shu C.-W. Efficient implementation of weighted ENO schemes //J. Comput. Phvs. 1996. Vol. 126. P. 202-228.

[16] Кудрявцев A.H., Поплавская T.B., Хотяновский Д.В. Применение схем высокого порядка точности при моделировании нестационарных сверхзвуковых течений // Математическое моделирование. 2007. Т. 19, № 7. С. 39-55.

[17] Коган М.Н. Динамика разреженного газа. М.: Наука, 1967. 440 с.

[18] Ivanov M.S., Markelov G.N., Gimelshein S.F. Statistical Simulation of Reactive Rare ed Flows: Numerical Approach and Applications. AIAA Paper 98-2669, 1998.

[19] Черчиньяни К. Теория и приложения уравнения Больцмана: Пер. с англ. М.: Мир, 1978. 496 с.

Поступила в редакцию 31 марта 2011 г.