«КВАЗИГАЗОДИНАМИЧЕСКИЕ УРАВНЕНИЯ И МЕТОДЫ РАСЧЕТА ВЯЗКИХ ТЕЧЕНИЙ Москва Научный Мир 2007 УДК 519.633:533.5 Т. Г. Елизарова. Квазигазодинамические уравнения и методы расчета вязких течений. Лекции по математическим ...»
Т. Г. Елизарова
КВАЗИГАЗОДИНАМИЧЕСКИЕ
УРАВНЕНИЯ И МЕТОДЫ РАСЧЕТА
ВЯЗКИХ ТЕЧЕНИЙ
Москва
Научный Мир
2007
УДК 519.633:533.5
Т. Г. Елизарова. Квазигазодинамические уравнения и методы расчета вязких течений. Лекции по математическим моделям и
численным методам в динамике газа и жидкости. М.: Научный Мир, 2007. – 350 с.
Монография посвящена современным математическим моделям и основанным на них численным методам решения задач динамики газа и жидкости. Приведены две взаимосвязанные математические модели, обобщающие систему уравнений Навье–Стокса и отличающиеся от нее дополнительными диссипативными слагаемыми с малым параметром в качестве коэффициента. Новые модели получили название квазигазодинамических и квазигидродинамических уравнений. На базе этих уравнений построены эффективные конечно-разностные алгоритмы расчета вязких нестационарных течений и приведены примеры численных расчетов. Универсальность, эффективность и точность построенных алгоритмов обеспечивается выполнением для них интегральных законов сохранения и теоремы о балансе энтропии.
Книга представляет собой курс лекций и предназначена для научных работников и инженеров, занимающимся построением численных алгоритмов и проведением практических расчетов течений газа и жидкости, а также студентам и аспирантам соответствующих вузов.
c Елизарова Т.Г., c Научный Мир, Оглавление Введение 1 Построение уравнений газовой динамики 1.1 Процедура осреднения.................. 1.2 Интегральные законы сохранения............ 1.3 Законы сохранения в дифференциальном виде.... 1.4 Уравнения Эйлера и Навье–Стокса........... 1.5 Квазигазодинамические уравнения........... 2 Элементы кинетической теории газов 2.1 Уравнение Больцмана................... 2.2 Равновесная функция распределения и система Эйлера 2.3 Уравнения Навье–Стокса................. 2.4 Уравнение Бхатнагара–Гросса–Крука......... 2.5 Средние характеристики движения частиц...... 2.6 Коэффициенты переноса в равновесном газе..... 2.7 Численное моделирование течений разреженного газа 2.8 Кинетически-согласованные разностные схемы.... 3 Квазигазодинамические уравнения 3.1 Регуляризованное кинетическое уравнение....... 3.2 Кинетический вывод КГД уравнений.......... 3.3 КГД уравнения в форме законов сохранения..... 3.4 Коэффициенты диссипации............... 3.5 Система Навье–Стокса как асимптотика КГД системы 3.6 Модель для течений с внешними источниками.... 3.7 Уравнение баланса энтропии............... 4 КГД уравнения и системы координат 4.1 КГД уравнения в произвольной системе координат.. 4.2 Декартова система координат.............. 4.3 Цилиндрическая система координат.......... 4 Оглавление 5 Алгоритмы решения задач газовой динамики 5.1 Система уравнений для плоских двумерных течений. 5.2 Система уравнений в цилиндрической геометрии... 5.3 Граничные условия.................... 5.4 Безразмерный вид уравнений.............. 5.5 Разностная аппроксимация................ 5.6 Введение искусственной диссипации.......... 5.7 Задача о распаде сильного разрыва........... 5.8 Задача о течении в окрестности цилиндра....... 5.9 Задача о течении в канале с уступом.......... 5.10 Численный алгоритм расчета дозвуковых течений.. 5.11 Устойчивость и точность КГД алгоритмов....... 6 Расчеты течений на неструктурированных сетках 6.1 Выбор сетки и построение контрольного объема... 6.2 Аппроксимация системы уравнений........... 6.3 Аппроксимация частных производных......... 6.4 Разностные схемы для двумерных течений...... 6.5 Аппроксимация граничных условий.......... 6.6 Расчет течения в окрестности цилиндра........ 7 Течения вязкой несжимаемой жидкости 7.1 Квазигидродинамическая система уравнений..... 7.3 Отрывное течение за обратным уступом........ 7.4 Тепловая конвекция в квадратной области...... 7.5 Тепловая конвекция при низких числах Прандтля.. 7.6 Конвекция Марангони в невесомости.......... 7.7 Течение в кубической каверне с подвижной крышкой 8.1 Молекулярная модель и функции распределения... 8.2 Системы координат и некоторые интегралы...... Оглавление 8.6 Примеры численных расчетов.............. 9.1 Исходная кинетическая модель............. 9.2 Построение моментных уравнений........... 9.3 Вычисление обменных членов.............. 9.4 Определение частот столкновений........... 9.5 Квазигазодинамические уравнения для смеси газов. 9.6 Одножидкостные приближения............. 9.7 КГДМ система для одномерного течения....... 9.8 Структура ударной волны в смеси гелия и ксенона.. B Течение вязкого сжимаемого газа в микроканалах Введение Предлагаемая читателю книга посвящена изложению нового эффективного подхода к численному моделированию широкого круга течений вязкого газа и жидкости. Этот подход основан на использовании двух математических моделей, обобщающих систему уравнений Навье–Стокса и отличающиеся от нее дополнительными диссипативными слагаемыми с малым параметром в качестве коэффициента. Новые модели получили название квазигазодинамических и квазигидродинамических (КГД) уравнений. Приведен способ построения и теоретическое обоснование КГД уравнений. На базе этих уравнений построены конечно-разностные алгоритмы расчета вязких нестационарных течений газа и жидкости. Универсальность, эффективность и точность построенных алгоритмов обеспечивается заложенными в них математическими свойствами и выполнением для них интегральных законов сохранения и теоремы о балансе энтропии. Простота и удобство численной реализации позволяют строить численные алгоритмы для расчета сложных течений и использовать параллельные технологии программирования для ускорения счета. Последнее особенно важно при расчетах нестационарных течений.
Описание течений газа и жидкости на основе уравнений Навье– Стокса имеет богатую историю. В настоящее время созданы и успешно применяются многочисленные коммерческие пакеты программ, реализующие численные алгоритмы решения этих уравнений. Однако используемые в них подходы нельзя считать совершенными. В разное время предпринимались попытки расширить возможности описания течений, заложенные в системе Навье– Стокса. Однако предлагаемые модели оказались существенно сложнее классической системы и не нашли применения в практических расчетах.
Система квазигазодинамических уравнений, расширяющая возможности модели Навье–Стокса, впервые появилась в ходе исследований, выполненных в восьмидесятых годах небольшой группой сотрудников Института прикладной математики АН СССР им. М.В. Келдыша под руководством профессора, а ныне членакорреспондента Российской Академии Наук, Б.Н. Четверушкина.
В самом начале этих исследований в 1982 г. автору посчастливилось включиться в эту работу, и с его непосредственным участием был выписан первый вариант квазигазодинамических уравнений.
Эти уравнения отличались от классических уравнений динамики газа дополнительными слагаемыми, имеющими вид вторых пространственных производных. Новые модели сразу позволили построить эффективные численные алгоритмы решения уравнений Эйлера, а впоследствии и уравнений Навье–Стокса.
Позднее в работах Ю.В. Шеретова квазигазодинамические уравнения были представлены в виде законов сохранения, детально исследованы и теоретически обоснованы. Кроме того, им была построена родственная этим уравнениям квазигидродинамическая система. Принципиальным и существенным отличием КГД подхода от теории Навье–Стокса явилось использование процедуры пространственно-временного осреднения для определения основных газодинамических величин плотности, скорости и температуры. Дополнительное сглаживание по времени явилось причиной возникновения в уравнениях дополнительных диссипативных слагаемых, которые формально отличают КГД системы от системы Навье–Стокса.
В стационарном случае обе КГД системы отличаются от уравнений Навье–Стокса и друг от друга дивергентными членами, имеющими формальные асимптотические порядки малости O(Kn2 ) при Kn 0, где Kn – число Кнудсена. Приближением пограничного слоя для новых уравнений, также как и для уравнений НавьеСтокса, служат уравнения Прандтля. Влияние добавочных членов незначительно для стационарных и квазистационарных газодинамических течений при малых числах Кнудсена. Однако для сильно нестационарных течений, а также при числах Kn, близких к единице, их вклад становится существенным. Именно в этом классе задач следует искать преимущества новых моделей. При численном моделировании дополнительные слагаемые проявляют себя как эффективные регуляризаторы.
Каждой из двух КГД систем соответствует свой способ решеВведение ния проблемы замыкания. Квазигазодинамические уравнения следует использовать при моделировании течений идеального политропного газа, а квазигидродинамические при исследовании движений газов и жидкостей с более общим уравнением состояния.
Монография состоит из девяти глав и четырех приложений.
В первой главе сформулированы общие идеи, позволяющие наряду с классической системой уравнений газовой динамики системой уравнений Навье–Стокса, построить две новые математические модели для описания вязких течений квазигазодинамическую и квазигидродинамическую системы уравнений.
Во второй главе приведены элементы кинетической теории, необходимые для дальнейшего изложения, и описан способ построения кинетически-согласованных разностных схем, которые послужили основой для первых вариантов КГД уравнений.
В третьей главе, которая имеет ключевой характер, приведены два возможных способа построения квазигазодинамических уравнений. Первый основан на использовании известной кинетической модели движения частиц, второй основан на интегральных законах сохранения, выписанных для малого, но конечного объема. Здесь же изложен способ записи полученных уравнений в виде законов сохранения, выписано уравнение баланса энтропии и прослежена связь КГД системы с уравнениями Навье-Стокса. Одним из результатов этой главы является получение приближенной формулы для вычисления коэффициента объемной вязкости.
В четвертой главе КГД уравнения выписаны для произвольной ортогональной системы координат.
Пятая, самая объемная глава, посвящена построению эффективных конечно-разностных алгоритмов решения КГД уравнений для численного моделирования газодинамических течений с использованием ортогональных сеток. Здесь же приведены примеры расчетов тестовых задач.
В шестой главе предложенные алгоритмы обобщаются на случай неструктурированных двумерных пространственных сеток. Седьмая глава посвящена исключительно второй КГД системе системе квазигидродинамических уравнений. На основе этих уравнений построены новые эффективные алгоритмы расчета нестационарных течений вязкой несжимаемой жидкости и приведены примеры моделирования двумерных и трехмерных течений.
В восьмой и девятой главах построены обобщения квазигазодинамических уравнений для моделирования течений газа с неравновесностью между поступательной и вращательной степенями свободы, и уравнения для описания бинарной смеси нереагирующих между собой газов.
В приложениях приведен детальный вывод квазигазодинамической системы уравнений для плоского одномерного течения, а также примеры расчета течения разреженного газа в тонком капилляре, течения в одномерной ударной волне и численное моделирование ламинарно–турбулентного перехода в потоке газа за обратным уступом.
Данная книга основана на курсе лекций, прочитанных автором для аспирантов кафедры математики Физического факультета МГУ в 2004 – 2008 гг. и изданных в 2005г. в качестве учебного пособия.
Всем коллегам, в разное время принимавшим участие в развитии КГД теории, я искренне благодарна за сотрудничество. Их вклад отражен в публикациях, послуживших основой этой монографии.
Я искренне благодарна Борису Николаевичу Четверушкину, идеи которого послужили основой КГД подхода. Я глубоко благодарна Юрию Владимировичу Шеретову за помощь в написании первой главы и ряда разделов в главах 2, 3, 5 и 7. Я признательна Александру Андреевичу Самарскому и Алексею Георгиевичу Свешникову за постоянную поддержку этого научного направления.
Глава Построение уравнений газовой динамики на основе законов сохранения В этой главе излагаются физические принципы, положенные в основу вывода уравнений классической газовой динамики и новых КГД уравнений квазигазодинамической и квазигидродинамической систем. КГД уравнения, также как и уравнения Навье–Стокса, являются следствием интегральных законов сохранения, имеют диссипативный характер и могут быть получены из общей системы уравнений сохранения. Принципиальным и существенным отличием КГД уравнений от уравнений Навье–Стокса представляется использование процедуры пространственно-временного осреднения для определения основных гидродинамических величин плотности, скорости и температуры. Использование пространственных средних приводит к системе уравнений Навье–Стокса. Для пространственновременных средних предложено два варианта замыкания общей системы уравнений, которые приводят к двум КГД системам. Выражения для векторов плотности потока массы, потока тепла и тензора вязких напряжений для КГД систем в этой главе приведены без вывода. Два способа построения этих величин для первой КГД системы представлены в главе 3. Обсуждается физический смысл вектора плотности потока массы. Изложение этой главы в основном основано на материалах [53, 111, 114, 120].
1.1 Процедура осреднения Рассмотрим одноатомный газ, состоящий из достаточно большого числа N атомов-шариков радиуса r0 и массы m0. В евклидовом пространстве R3 выберем декартову систему координат (x1, x2, x3 ) 1.1. Процедура осреднения и время t.
Движение каждого атома может быть описано уравнениями механики Ньютона. Однако такой подход к моделированию задач газовой динамики оказывается очень далеким от практики, поскольку N очень велико. Кроме того, возникают проблемы с определением начальных условий задачи и с последующим осреднением получающихся величин, которое необходимо для вычисления измеряемых величин плотности, скорости и температуры среды.
В классической гидродинамике используется другой подход, основанный на переходе от большого числа отдельных частиц к сплошной среде с помощью процедур осреднения. Эти процедуры могут быть выбраны по-разному.
1.1.1 Пространственные средние В теории Навье–Стокса используются так называемые мгновенные пространственные средние, которые определяются следующим образом.
Пусть V шар радиуса RV с центром в точке x (рис. 1.1).
Пусть атом находится в шаре, если его центр принадлежит этому шару. Пусть NV (t) число молекул в объеме V в момент времени t. Определим плотность, средний импульс и среднюю энергию единицы объема как где i (t) скорость i-той частицы в момент времени t, m0 ее масса. В приведенных выражениях = (x, t) удельная внутренняя энергия. Введем температуру T, которая определяется из выражения 12 Гл.1. Построение уравнений газовой динамики где cv = cp R удельная теплоемкость при постоянном объеме, cp удельная теплоемкость при постоянном давлении, R = kB /m0 = = R /M газовая постоянная, kB = 1.38 · 1023 Дж/K постоянная Больцмана, R – универсальная газовая постоянная, M – молекулярная масса газа, = cp /cv показатель адиабаты.
личин достаточно гладкими, то есть непрерывно дифференцируемыми столько раз, сколько это потребуется.
1.1.2 Пространственно-временные средние Определим пространственно-временные средние, которые отличаются от определений (1.1)–(1.3) дополнительным сглаживанием по времени. Для этого введем дополнительно некоторый интервал времени t и определим плотность, импульс и энергию единицы объема как 1.1. Процедура осреднения Предположим далее, что помимо двух пространственных масштабов Rmax и Rmin, существуют также два временных масштаба tmax tmin, таких, что при любом tmax > t > tmin значения указанных средних практически постоянны и не зависят от V и t. Характерные значения указанных масштабов связаны между собой как tmax Rmax /c, tmin Rmin /c, где c – скорость звука, определяющая скорость распространения возмущений в газе. Тогда соответствующие средние можно рассматривать как газодинамические величины плотность, импульс и полную энергию. Предположим также, что эти функции достаточно гладкие.
Иерархия характерных пространственных и временных масштабов для газодинамических течений детально обсуждается, в частности, в [63, 64] и [122].
Введение дополнительного сглаживания по времени при определении газодинамических величин представляется естественным по многим причинам. В экспериментах измерение всех газодинамических величин осуществляется за конечное время, что автоматически приводит к сглаживанию по некоторому временному интервалу. Число частиц в малом объеме V c прозрачными границами естественным образом меняется на временах порядка RV /c за счет частиц, хаотическим образом пересекающих его границу. При достаточно большом числе частиц в объеме осреднения V пространственные и пространственно-временные средние могут быть очень близки, что соответствует эргодической гипотезе об идентичности мгновенных пространственных и пространственно-временных средних.
В дальнейшем мы не будем отождествлять мгновенные пространственные и пространственно-временные средние, сохраняя 14 Гл.1. Построение уравнений газовой динамики при этом привычные обозначения для всех газодинамических величин. При необходимости будем обозначать мгновенные пространственные средние индексом s, а пространственно-временные средние – индексом st.
1.1.3 Преобразование Галилея Рассмотрим две инерциальные системы координат. Пусть K система координат, движущаяся относительно исходной системы координат K с постоянной скоростью U (рис. 1.2). Тогда координаты материальной точки x и время t в K связаны с координатами x и t в системе K соотношениями Рис. 1.2. Пространственно-временные средние и преобразование Галилея В момент t0 обе системы совпадают. Формулы (1.8)–(1.9) были названы Ф. Франком "преобразованиями Галилея" [80]. Преобразование Галилея сначала было выписано для материальной точки в классической механике Ньютона. При этом существенно использовалась гипотеза (1.9) об абсолютности времени. Согласно принципу Галилея инерциальные системы координат, движущиеся друг относительно друга с постоянными скоростями, равноправны с точки 1.1. Процедура осреднения зрения протекания в них механических явлений, то есть уравнения движения в этих системах координат инвариантны.
Однако для того, чтобы исследовать инвариантность уравнений гидродинамики, формул (1.8)–(1.9) недостаточно. Необходимо знать еще, как изменяются при переходе от K к K макроскопические параметры плотность, гидродинамическая скорость u и температура T. Ответ на последний вопрос зависит от используемой процедуры осреднения при определении этих макропараметров.
При использовании мгновенных пространственных средних имеют место равенства Инвариантность уравнений Навье–Стокса, которые строятся на основе пространственных средних, по отношению к преобразованиям (1.8)–(1.9) проверяется непосредственно.
Для пространственно-временных средних объемы, по которым проводится осреднение в неподвижной V и движущейся V системах координат, будут различаться – неподвижный для наблюдателя из системы координат K объем V представляется движущимся наблюдателю, связанному с системой координат K, и наоборот.
Поэтому равенства (1.10) выполняются не точно, а лишь приближенно:
Таким образом, плотность, скорость и температура оказываются относительными, и инвариантность по отношению к преобразованиям Галилея нарушается.
Сходная ситуация возникает в релятивистской гидродинамике, где несправедлива гипотеза (1.9) об абсолютности времени. Поэтому уравнения релятивистской гидродинамики также не инвариантны относительно преобразования Галилея.
Для существования инвариантности в уравнениях гидродинамики необходимо выполнение двух условий:
• определение гидродинамических величин как пространственных средних (то есть в неподвижной и подвижной системах 16 Гл.1. Построение уравнений газовой динамики координат используется один и тот же объем осреднения, что • абсолютность времени (время одинаково в неподвижной и движущейся системах координат, следовательно, t = t ).
Инвариантность нарушается, если не выполняется хотя бы одно из этих двух условий. Для уравнений Навье–Стокса выполнены оба условия. Для уравнений релятивистской механики не выполнено второе условие. Для уравнений, которые строятся на основе пространственно-временных средних, не выполняется первое условие. Детальному обсуждению справедливости преобразования Галилея для уравнений гидродинамики посвящена работа Ю.В. Шеретова [118].
1.1.4 Уравнение неразрывности В основу гидродинамики положен принцип сохранения массы, или уравнение баланса массы, которое записывается как где jm вектор плотности потока массы. Это уравнение носит название уравнения неразрывности. Интегральный вид уравнения неразрывности означает, что изменение массы в некотором замкнутом объеме V определяется потоком массы jm, протекающим через его границу. При этом при построении уравнения (1.12) определения газодинамических величин не используются.
Для пространственных средних полагается (см., например, [93, 114]), что плотность потока массы равна импульсу единицы объема В этом случае уравнение (1.12) имеет вид 1.1. Процедура осреднения и для уравнения неразрывности (1.15), также как и для самих пространственных средних, преобразование Галилея выполняется.
Покажем, что при использовании определений пространственновременных средних уравнение неразрывности вида (1.15) для этих средних не выполняется. Проинтегрируем уравнение (1.15) по малому интервалу времени t. Тогда, согласно определению st, первое слагаемое примет вид Второе слагаемое (1.15) примет вид поскольку Таким образом для пространственно-временных средних плотность потока массы jmst может не совпадать с импульсом единицы объема st ust. Это отражает тот факт, что даже за малое время t мгновенные значения плотности и импульса единицы объема успевают измениться.
В этом случае выражение для плотности потока массы можно записать в более общем виде. Введем малую добавку к скорости, которую обозначим как wst, и запишем плотность потока массы как Как было показано выше на основе качественных соображений, для уравнений газовой динамики, основанных на пространственновременных средних величинах, преобразование Галилея не выполняется.
1.2 Интегральные законы сохранения В евклидовом пространстве R3 выберем инерциальную декартову систему координат (x1, x2, x3 ). Пусть (e1, e2, e3 ) соответствующий ей ортонормированный базис единичных векторов, t время. Будем использовать следующие стандартные обозначения для величин, характеризующих течения сжимаемой вязкой теплопроводной среды:
= (x, t) давление, = (x, t) удельная внутренняя энергия, T = T (x, t) температура, s = s(x, t) удельная энтропия.
Предположим, что среда является двухпараметрической, то есть среди пяти термодинамических параметров, p,, T, s независимы лишь два. При этом заданы уравнения состояния Рис. 1.3. К выводу уравнений сохранения = (t), ориентированной полем внешних единичных нормалей n (рис. 1.3). Будем считать, что объем V (t) возникает из объема V0 = V (t0 ), где t0 начальный момент времени, путем непрерывной деформации, обусловленной его перемещением вдоль траекторий, определяемых некоторым векторным полем v. Запишем известное тождество Эйлера–Лиувилля [93]:
1.2. Интегральные законы сохранения в котором D = /t + v · дифференциальный оператор, = некоторое непрерывно дифференцируемое скалярное = (x, t) или векторное поле, dV Пусть в каждой точке x области течения в момент времени t определен вектор jm = jm (x, t), называемый плотностью потока массы. Пусть объем V0 перемещается вдоль траекторий, определяемых векторным полем v = jm /. Согласно гипотезе (1.12) это обеспечивает сохранение массы в объеме при его смещении. Тогда тождество Эйлера–Лиувилля примет вид Приведем постулаты, на основе которых будем строить уравнения газовой динамики.
В качестве первого постулата примем закон сохранения массы (1.12, 1.13), который для единообразия дальнейшего изложения запишем в эквивалентном интегральном виде как Вторым постулатом будет служить закон сохранения импульса где d элемент площади поверхности около единичного вектора n. Скорость изменения импульса в объеме V равна сумме приложенных к нему сил. Первый интеграл в правой части (1.20) есть объемная сила, действующая со стороны внешнего поля; второй определяет силы давления и внутреннего вязкого трения, приложенные к поверхности. Величина P = P (x, t) называется тензором внутренних напряжений. Символ (n · P ) обозначает свертку (скалярное произведение) вектора n и тензора второго ранга P, осуществляемую по первому индексу тензора P. Соответственно, запись (P · n) означает, что свертка P и n идет с участием второго индекса тензора P. В случае симметричного тензора P имеем (n · P ) = (P · n).
Третьим постулатом является закон сохранения полной энергии Первый интеграл в правой части (1.21) равен мощности внешних массовых сил, приложенных к объему V ; второй интерпретируется как мощность поверхностных сил давления и внутреннего вязкого трения. Последний член в (1.21) описывает приток энергии в единицу времени через поверхность за счет процессов теплопроводности. Конкретные выражения для векторных полей A = A(x, t) и q = q(x, t) будут приведены ниже.
Следующий постулат выражает закон сохранения момента импульса:
Он представлен в классической форме. Внутренние моменты, а также распределенные массовые и поверхностные пары, не учитываются. Символ используется для обозначения векторного произведения двух векторов.
Второй закон термодинамики, являющийся нашим пятым постулатом, имеет вид Поверхностный интеграл в правой части (1.23) определяет скорость изменения энтропии в объеме V за счет теплового потока. Он может быть как положительным, так и отрицательным. Последний интеграл всегда неотрицателен и дает увеличение энтропии за счет внутренних необратимых процессов. Величина X называется производством энтропии.
1.3. Законы сохранения в дифференциальном виде 1.3 Законы сохранения в дифференциальном Чтобы перейти от интегральных соотношений (1.19) – (1.23) к соответствующим дифференциальным уравнениям, воспользуемся формулой Эйлера–Лиувилля (1.18) о дифференцировании интеграла, взятого по движущемуся материальному объему. При этом будем считать, что все основные макроскопические параметры среды являются достаточно гладкими функциями пространственных координат и времени.
Полагая =, u, (u2 /2 + ), [x u] и s, и учитывая произвольность V, получим дифференциальные уравнения баланса массы импульса полной энергии момента импульса и энтропии Здесь (jm u) – тензор второго ранга, полученный в результате прямого произведения векторов jm и u. При вычислении дивергенции от тензора второго ранга свертка осуществляется по его первому индексу. В уравнении (1.27) символом Pij обозначен портрет тензора P в базисе (e1, e2, e3 ). По дважды повторяющимся индексам i и j идет суммирование.
Покажем, что полученная система уравнений (1.24) – (1.26) является диссипативной. Пусть все величины, входящие в эту систему, определены. Предположим, что течение газа происходит в замкнутом сосуде объема V0 с непроводящей тепло твердой стенкой 0.
Добавим к системе (1.24) – (1.26) начальные условия а также граничные условия Здесь 0 = 0 (x) > 0, u0 = u0 (x), T0 = T0 (x) > 0 заданные значения плотности, скорости и температуры в момент времени t = 0.
Первое из условий (1.30) означает, что газ прилипает к стенкам сосуда; второе обеспечивает отсутствие потока массы через границу;
третье влечет равенство нулю на 0 нормальной составляющей теплового потока. Интегрируя (1.28) по объему V0 и принимая во внимание (1.29), (1.30), приходим к неравенству для полной термодинамической энтропии S(t) = V0 s dx. Из (1.31) следует, что величина S(t) является неубывающей функцией времени.
Система (1.24) – (1.28) не является замкнутой. Необходимо представить величины jm, P, q, A, X как функции макроскопических параметров среды и их производных. Проблема замыкания может быть решена различными способами.
1.4 Уравнения Эйлера и Навье–Стокса Изложим сначала классический подход, в котором для определения гидродинамических величин используются мгновенные пространственные средние [72, 79, 93]. В этом случае вектор плотности потока массы jm в любой точке (x, t) совпадает со средним импульсом 1.4. Уравнения Эйлера и Навье–Стокса единицы объема u, и первое замыкающее соотношение имеет вид Далее вводится представление о силах давления и внутреннего вязкого трения, мгновенно действующих на поверхность материального объема. Закон движения последнего выбирают таким же, как и в механике твердого тела. Это допущение называют принципом отвердевания. Уравнение баланса момента импульса (1.27) является следствием закона сохранения импульса (1.25) при условии симметричности тензора напряжений P. В теории ньютоновских сред для P = PN S используется выражение где тензор второго ранга, называемый навье–стоксовским тензором вязких напряжений, I – единичный тензор-инвариант второго ранга.
Верхним индексом Т обозначена операция транспонирования.
Тепловой поток q = qN S задается в соответствии с законом Фурье Для идеальных одноатомных газов при малых числах Кнудсена гипотезы (1.34), (1.35) подтверждаются кинетическим расчетом.
Работу в единицу времени поверхностных сил давления и внутреннего вязкого трения вычисляют по той же формуле, что и в механике твердого тела, а именно Считают, что удельная термодинамическая энтропия подчиняется дифференциальному тождеству Гиббса Уравнение баланса энтропии (1.28) может быть получено на основе тождества Гиббса с использованием законов сохранения массы, импульса и энергии (1.25) – (1.26). При этом производство энтропии X = XN S имеет вид где (N S : N S ) = 3 (N S )ij (N S )ij – двойное скалярное произi,j= ведение двух одинаковых тензоров. Величина называется диссипативной функцией, величина которой определяет диссипацию энергии за счет сил вязкого трения. Заметим, что правая часть равенства (1.38) неотрицательна. Подстановка выражений (1.32) – (1.36) в уравнения (1.24) – (1.26) дает классическую систему Навье–Стокса для вязкой сжимаемой теплопроводной среды импульса и полной энергии Первое соотношение (1.39) называется уравнением баланса массы или уравнением неразрывности. Равенства (1.40) и (1.41) выражают законы сохранения импульса и полной энергии, соответственно.
Система становится замкнутой, если ее дополнить граничными и начальными условиями и уравнениями состояния а также выражениями для вычисления положительных коэффициентов динамической вязкости µ, второй вязкости и теплопроводности.
1.4. Уравнения Эйлера и Навье–Стокса Для случая идеального политропного газа, состоящего из упругих шариков, зависимости (1.42) выбираются в виде Первое соотношение (1.43) называется уравнением Менделеева– Клапейрона, или уравнением состояния идеального газа. Второе соотношение характеризует газ как политропный. В этом случае удельная термодинамическая энтропия выражается формулой Здесь cv = R/( 1), cp = R/( 1).
Зависимости µ = µ(, T ) и = (, T ) могут быть найдены либо экспериментально, либо методами кинетической теории газов. Для идеального политропного газа вязкость и теплопроводность зависят только от температуры и могут быть аппроксимированы функциями в которых µ1 известное значение коэффициента динамической вязкости при температуре T1, заданный показатель температурной зависимости из промежутка [0.5, 1], P r число Прандтля.
Приближенно коэффициент второй (объемной) вязкости можно аппроксимировать формулой Этот коэффициент всегда положителен и связан с наличием внутренних степеней свободы молекулы. Для одноатомного газа = 5/ и = 0. В противном случае 1 < < 5/3, и > 0. Указанная формула была получена на основе кинетической теории в [193] для газа с вращательными степенями свободы. Эта же формула получается на основе анализа КГД уравнений для произвольного (см. с. 66). Влияние коэффициента второй вязкости на форму профиля плотности в ударной волне обсуждается в приложении B.
26 Гл.1. Построение уравнений газовой динамики Для других сред (например, газа Ван-дер-Ваальса) зависимости (1.43) и (1.45) могут видоизменяться.
Уравнения (1.39) – (1.41) являются инвариантными относительно преобразований Галилея. Это соответствует принципу относительности Галилея об одинаковом виде законов движения в различных инерциальных системах отсчета.
Выписанная система уравнений удовлетворяет закону сохранения момента импульса и уравнению баланса энтропии в виде где диссипативная функция Вытекающий отсюда закон неубывания полной энтропии в замкнутом адиабатически изолированном объеме указывает на необратимый (диссипативный) характер системы Навье–Стокса.
Если в уравнениях (1.39) – (1.41) пренебречь эффектами вязкости и теплопроводности, то придем к классической системе уравнений Эйлера.
1.5 Квазигазодинамические и квазигидродинамические уравнения Для мгновенных пространственных средних справедливо равенство (1.32). Для пространственно-временных средних это равенство в общем случае не выполняется (см. раздел 1.1.4). Возможный выбор величин jm, P, A, q, X в предположении, что jm, вообще говоря, не равен u, приведен далее.
По аналогии с выражениями (1.33), (1.36) будем полагать, что тензор напряжений и работа сил давления и вязкого трения связаны с тензором вязких напряжений соотношениями вида:
1.5. Квазигазодинамические уравнения Для пространственно-временных средних было построено два варианта замыкания общей системы уравнений (1.24) – (1.28). Получающиеся в результате этого системы были названы квазигазодинамической и квазигидродинамической системами уравнений. Сокращенно обе системы именуются одинаково КГД системы. Присутствующие в КГД уравнениях добавки, пропорциональные малому параметру, связаны с дополнительным осреднением (сглаживанием) по времени при определении газодинамических параметров.
Первая система описывает поведение идеального политропного газа. Первый вариант этой системы был получен на основе кинетической модели в 80-е годы [50, 51, 105]. Позднее эта система была представлена в виде законов сохранения [111, 114]. Вторая система была получена позднее Ю.В. Шеретовым на основе анализа уравнений сохранения в дифференциальном виде [110, 111]. Эта система описывает течение газа с более общим уравнением состояния, и в приближении = const может использоваться для моделирования течений вязкой несжимаемой жидкости.
Там, где это не приводит к путанице, будем использовать аббревиатуру КГД. В специальных случаях будем называть эти системы полностью.
1.5.1 Квазигазодинамическая система Пусть за некоторое физически бесконечно малое время мгновенные значения средней плотности, среднего импульса и энергии единицы объема успевают измениться. Для идеального политропного газа, то есть для газа с уравнением состояния был найден способ замыкания общей системы (1.24) – (1.28), который приводит к квазигазодинамической системе уравнений. Вариант построения этой системы на основе кинетической модели будет изложен в главе 3. Здесь приведем сразу замыкающие соотношения, которые имеют вид где – некоторый малый коэффициент размерности времени, который в дальнейшем будем называть параметром релаксации, или сглаживания. При = 0 приведенные выше выражения для jm, P и q вырождаются в соответствующие величины для уравнений Навье– Стокса. Способы нахождения будут обсуждаться далее.
Вектор A и неотрицательная величина производства энтропии Х записываются в виде Заметим, что производство энтропии для КГД системы представляет собой производство энтропии для уравнений Навье–Стокса с дополнительными членами, которые являются квадратами левых частей классических уравнений Эйлера в стационарном случае с положительными коэффициентами. Для КГД уравнений производство энтропии неотрицательно.
1.5. Квазигазодинамические уравнения Подставляя выписанные выше значения векторов и тензора вязких напряжений в общую систему уравнений (1.24) – (1.28) получим квазигазодинамическую систему в виде Здесь = div(u) приближенное значение плотности в точке (x, t + ). Величина определяется выбором значения плотности в первом слагаемом правой части (1.20) в сдвинутой по времени точке с помощью соотношения = + /t, где /t + div(u) = 0.
Уравнение баланса энтропии было получено на основе недивергентного вида КГД системы в форме (1.28) в [115] Вариант его построения для течений с внешними источниками энергии приведен в последнем параграфе главы 3.
1.5.2 Квазигидродинамическая система Второй способ решение проблемы замыкания системы (1.24) – (1.28) был предложен Ю.В. Шеретовым в работах [111, 114].
Пусть за любое физически бесконечно малое время успевает измениться только мгновенное значение среднего импульса единицы объема, а изменениями мгновенных значений плотности и температуры можно пренебречь. В этом случае для газа с уравнениями состояния (1.16) и удовлетворяющему тождеству Гиббса (1.37) величины jm, P, q, A и X были построены в виде:
причем Подставив выражения (1.58), (1.59) и (1.61) вместо величин jm, P и A в (1.24) – (1.26), получим квазигидродинамическую систему уравнений Шеретова КГД система (1.64) – (1.66) становится замкнутой, если дополнить ее уравнениями состояния (1.16), а коэффициенты µ, и представить как функции макроскопических параметров среды.
Подстановка выражений (1.58), (1.60) и (1.62) в (1.28) дает уравнение баланса энтропии в котором неотрицательная диссипативная функция.
1.5. Квазигазодинамические уравнения 1.5.3 Вектор плотности потока массы и параметр Формально во всех уравнениях КГД системы присутствуют дополнительные по сравнению с системой Навье–Стокса диссипативные слагаемые, представляющие собой вторые пространственные производные от плотности, скорости и давления, перед которыми стоит численный коэффициент. При 0 КГД уравнения переходят в уравнения Навье–Стокса.
Коэффициент = (, T ), названный параметром релаксации, или сглаживания, связан с включением в определение газодинамических величин осреднения по времени. Это дополнительное осреднение позволяет учесть влияние малых флуктуаций числа частиц в объеме V, которым в классической гидродинамике пренебрегается. Величина параметра сглаживания может изменяться в широких пределах в зависимости от типа рассматриваемого течения.
Для определения величины параметра рассмотрим структуру вектора плотности потока массы. Для определенности остановимся на квазигазодинамической системе уравнений для стационарного течения идеального политропного газа. В этом случае плотность потока массы вычисляется согласно (1.48), (1.51) Преобразуем это выражение к виду Первый член в правой части описывает плотность потока массы, связанного с конвективным движением газа. Второе слагаемое поток, определяемый движением частиц во внешнем поле. Третье слагаемое поток массы за счет самодиффузии. Четвертое слагаемое так называемый термодиффузионный поток. Последнее слагаемое вклад в поток массы за счет градиента скорости. В реальных течениях все эти потоки тесно связаны между собой и не могут быть разделены.
Остановимся более подробно на третьем слагаемом. Плотность потока массы за счет самодиффузии имеет вид 32 Гл.1. Построение уравнений газовой динамики Коэффициент самодиффузии D для многих сред хорошо известен из экспериментов с изотопами. Согласно, например, [134] и [182], в политропном газе коэффициент самодиффузии связан с коэффициентом вязкости как где Sc число Шмидта. Согласно [134], число Шмидта в газе близко к единице и приближенно может быть получено как Сравнивая коэффициент самодиффузии (1.70) и выражение для этого коэффициента в (1.69), получаем, что для газа с уравнением состояния p = RT параметр релаксации равен Таким образом, с точностью до числа Шмидта величина совпадает с так называемым максвелловским временем релаксации m = = µ/p, то есть близка к среднему времени свободного пробега частиц в газе.
Рассмотрим теперь второе слагаемое уравнения (1.69). Если рассматривать каждую молекулу как броуновскую частицу, то плотность потока массы этих частиц можно связать с массовой плотностью внешних сил соотношением где коэффициент b называется подвижностью молекулы. Подвижность молекулы связана с коэффициентом самодиффузии соотношением Эйнштейна Подставляя в выражение для плотности потока массы во внешнем поле подвижность молекулы через соотношение Эйнштейна, вновь придем в формуле (1.72) для коэффициента релаксации.
1.5. Квазигазодинамические уравнения Термодиффузионный поток, который описывается четвертым слагаемым, согласно [72, 78] представляется в виде где kT безразмерная величина, называемая термодиффузионным отношением, которое определяет связь коэффициентов термодиффузии DT и самодиффузии D в виде DT = DkT. Сопоставляя четвертое слагаемое в (1.69) и выражение для jT, опять приходим к уже полученному нами ранее выражению для параметра релаксации (1.72) с точностью до коэффициента kT в виде = kT µ/(pSc).
В [114] для вычисления параметра релаксации была предложена более общая формула При учете формулы Лапласа для скорости звука в газе c2 = p/, выражение (1.73) преобразуется к виду (1.72).
Для плотных газов и жидкостей величина параметра сглаживания, выбранная в соответствии с (1.73), оказывается весьма малой, и влиянием членов в КГД уравнениях содержащих можно пренебречь. Например, для воздуха при температуре T = 200 = 1.4, Sc = 0.74, c = 3.4 · 104 см/сек, = µ/ = 0.15 см2 /cек, и = = 2.45·1010 сек. Для воды при аналогичных условиях = 1, Sc = 1, c = 1.45 · 105 см/сек, = 0.01см2 /cек, и = 4.75 · 1013 сек.
Для течений разреженного газа параметр сглаживания может быть достаточно большим. Вариант его вычисления для разреженных течений в каналах рассматривается в Приложении C. При описании быстропеременных или турбулентных течений вклад дополнительных вязких членов также может оказаться существенным.
Численное моделирование турбулентного течения за уступом в плоском канале и способ выбора параметра сглаживания в этой задаче рассматривается в Приложении D.
При проведении расчетов слагаемые с коэффициентом могут использоваться как удобные и эффективные регуляризаторы для нахождении численного решения. При этом величина параметра 34 Гл.1. Построение уравнений газовой динамики уже не связывается с молекулярными свойствами газа, а определяется шагом пространственной сетки и выбираться из условий сходимости и точности разностного решения задачи.
Выражение для плотности потока массы (1.68) включает в себя пространственную производную от давления, что делает уравнение неразрывности для обеих КГД систем уравнением второго порядка.
Поэтому при постановке начально-краевой задачи для КГД систем требуется дополнительное по сравнению с системой Навье–Стокса граничное условие. Это дополнительное условие может быть получено из условий для вектора плотности потока массы jm (1.48), (1.58) на границе.
Предположим, что граница представляет собой непроницаемую твердую стенку, и что внешняя сила равна нулю. Тогда условие непротекания для потока массы (jm · n) = 0 и условия непротекания для скорости (u · n) = 0 приводят к условию для давления на границе в виде p/n = 0.
Из формул для вычисления тензора вязких напряжений и теплового потока для обеих КГД систем (1.49) – (1.51) и (1.59) – (1.60), (1.63) следует, что условие непротекания для скорости (u · n) = приводит к тому, что пропорциональные дополнительные слагаемые в тепловом потоке и тензоре вязких напряжений обращаются в нуль, и для КГД систем на стенке выполняется условие Тем самым формулы для вычисления теплового потока и силы трения на твердой стенке для КГД уравнений совпадают с традиционными выражениями, полученными в рамках уравнений Навье– Стокса.
КГД системы отличаются от других систем, которые в разное время предлагались в работах [4, 15, 63, 64, 96, 122]. КГД уравнения принципиально отличаются от уравнений Барнетта [73, 104], добавочные 1.5. Квазигазодинамические уравнения слагаемые в которых имеют вид третьих пространственных производных и носят не диссипативный, а дисперсионный характер.
Для целого ряда задач сравнение численных результатов, полученных на основе модели Навье–Стокса и КГД уравнений, приведено в последующих главах. Для расчета разреженных течений такое сопоставление проведено и с расчетами по кинетическим моделям.
Было выполнено детальное сравнение результатов, полученных в расчетах по квазигазодинамическим уравнениям, с данными, полученными на основе уравнений Навье–Стокса и кинетических подходов для задач о стационарных течениях в окрестности пластины [146], диска [140] и полого цилиндра [161]. В этих расчетах получено, что все три модели дают очень близкие результаты для течений достаточно плотных газов. При этом численный алгоритм, основанный на КГД уравнениях, оказывается существенно проще в численной реализации и обладает большим запасом устойчивости.
С увеличением разреженности результаты расчетов начинают различаться. В этом случае данные КГД модели как правило, располагаются между данными, полученными с использованием кинетического подхода, и результатами расчета с использованием уравнений Навье–Стокса.
В качестве наглядного примера сопоставления КГД уравнений с моделью Навье–Стокса рассмотрим классическую задачу гидростатики о распределении давления в идеальном политропном газе, находящимся в однородном поле тяжести [72, 79]. В состоянии равновесия макроскопическая скорость u равна нулю, и распределение параметров газа не зависит от времени. В поле тяжести Земли F = = g, где g = 9.8 · 102 см/сек2.
В этом случае обе КГД системы существенно упрощаются и принимают одинаковый вид 36 Гл.1. Построение уравнений газовой динамики Отсюда непосредственно следуют условия механического равновесия, вытекающие из системы уравнений Навье–Стокса Пусть температура газа постоянна. Тогда приходим к классической формуле, определяющей распределение давления в газе где p0 заданное значение давления в точке x = 0. Формула (1.74) называется барометрической формулой, или формулой Лапласа. Таким образом формула Лапласа является точным решением уравнений Навье–Стокса и обеих КГД систем.
Другим точным решением, общим для КГД систем и системы Навье–Стокса, является решение задачи о течении Куэтта, приведенное в монографии [114]. В этой же монографии рассмотрен ряд других точных решений задач гидродинамики и прослежена связь таких решений, построенных в рамках классической модели и КГД систем.
Исторические сведения Отдельные экспериментальные факты течения жидкости и газа были установлены Архимедом (287–212 до н. э.), Паскалем (Blaise Pascal, 1588–1651), Торичелли (Evangelista Torricelli, 1608–1647) и Ньютоном (Isaac Newton, 1643–1727).
Основоположниками теоретической гидродинамики можно считать двух выходцев из Швейцарии, работавших в том числе и в России Леонарда Эйлера и Даниила, или Даниеля Бернулли.
Термин "гидродинамика"был введен Бернулли (Daniel Bernoulli, 1700–1783). Даламбер (Jean le Rond D’Alembert, 1717–1783) ввел закон сохранения массы для жидкости в виде уравнения неразрывности.
1.5. Квазигазодинамические уравнения Леонард Эйлер (Leonard Euler, 1707–1783) в 1755 г. выписал уравнения движения идеальной жидкости и развил их математическую теорию. Он вывел уравнение неразрывности, выражающее свойство сохранения массы в движущемся вместе с жидкостью материальном объеме. Он же получил уравнение баланса импульса в локальной форме без учета влияния вязкости. В те времена жидкость и газ рассматривали как сплошную среду в буквальном смысле слова. Молекулярный состав вещества в рассмотрение не принимался. Плотность определялась как формальный математический предел отношения массы жидкости в момент времени t в объеме к величине этого объема при его стремлении к нулю.
Работы Л. Эйлера продолжил Лагранж (Joseph Louis Lagrange, 1736–1813).
Клод Луи Навье (Claude Louis Navier, 1785–1836) вывел уравнения движения вязкой жидкости, пользуясь гипотезой взаимодействия молекул. История уравнений вязкой жидкости отсчитывается с того момента, когда Навье в 1822 г. сделал доклад об их простейшем варианте в несжимаемом случае. Соответствующая статья была опубликована через пять лет.
Джордж Стокс (George Gabriel Stokes, 1819–1903) получил уравнения движения вязкой жидкости на аксиоматической основе. По современным представлениям в качестве постулатов использовались интегральные законы сохранения массы, импульса и полной энергии в материальном объеме, движущемся вдоль интегральных кривых поля скорости. Его можно считать основателем современной гидродинамики.
Осборн Рейнольдс (Osborne Reynolds, 1842–1912), изучая движение вязкой жидкости, ввел понятия ламинарного и турбулентного течений и указал возможность резкого перехода от одного вида течения к другому.
Кинетическое обоснование уравнений гидродинамики было построено на основе уравнения Больцмана. Это уравнение для описания поведения функции распределения частиц моноатомного газа с бинарными столкновениями было выписано австрийским физиком Людвигом Больцманом (Ludwig Boltzmann, 1844–1906) в 1872 г.
[97, 109].
Глава Элементы кинетической теории газов В этой главе излагаются некоторые аспекты кинетической теории.
Эти сведения используются при выводе квазигазодинамических уравнений в главе 3, построении их обобщений (главы 8 и 9) и рассмотрении задач о структуре ударной волны и течении в микроканалах (приложения B и C). Приведено схематическое описание кинетического алгоритма DSMC, который в настоящее время широко применяется в численном моделировании течений разреженного газа. Расчеты в рамках этого подхода использовались для верификации КГД алгоритма при моделировании умеренно-разреженных течений. В последнем разделе изложен способ построения кинетически-согласованных разностных схем, дифференциальные аналоги которых послужили основой первых вариантов КГД уравнений. Изложение в этой главе опирается на работы [30, 31, 50–52] и [65, 68, 114, 120, 134, 182].
2.1 Уравнение Больцмана В 1872 г. Л. Больцман предложил интегро-дифференциальное кинетическое уравнение, которому суждено было стать классической моделью в теории разреженных одноатомных газов [12,62,65,73,134].
Это уравнение имеет вид и описывает эволюцию одночастичной функции распределения f = = f (x,, t). Здесь скорость отдельной частицы, которую будем рассматривать как атом-шарик массой m0, F действующая на частицы внешняя сила, отнесенная к единице массы, оператор Гамильтона. Функция f нормирована так, чтобы равенство 2.1. Уравнение Больцмана определяло вероятное. или ожидаемое, число dN частиц в элементе объема dx d около точки (x, ) фазового пространства координат и скоростей в фиксированный момент времени t.
Интеграл столкновений I(f, f ) есть нелинейный функционал, определяющий изменение функции распределения в результате парных столкновений. Конкретный вид этого интеграла можно найти в книгах [12, 65, 73, 134].
Важным и нужным нам в дальнейшем свойством интеграла столкновений является его ортогональность любому из так называемых столкновительных, или сумматорных, инвариантов То есть можно записать Это соотношение выражает законы сохранения массы, импульса и энергии частиц при их парном столкновении. Здесь и далее интегрирование выполняется по всему пространству скоростей частиц.
Зная функцию распределения f, можно определить гидродинамические величины плотность, скорость u, давление p, температуру T, удельную внутреннюю энергию, тензор вязких напряжений и тепловой поток q с помощью выражений Здесь c = u скорость хаотического движения частицы газа, или тепловая скорость, cv = 3R/2 удельная теплоемкость при постоянном объеме для одноатомного газа.
Интегрируя (2.1) c весами 1,, 2 /2 и пользуясь свойством (2.2), получим систему уравнений для макроскопических параметров которая, однако, не является замкнутой.
Для положительных решений f = f (x,, t) уравнения (2.1), в предположениях, что таковые существуют, обладают необходимыми свойствами гладкости и стремятся к нулю при | |, Л. Больцман доказал свою знаменитую H-теорему.
Предположим, что одноатомный газ находится в ограниченном объеме V0 с зеркально отражающей внутренней стенкой. Пусть заданы соответствующие начальные и краевые условия для функции распределения частиц в этом объеме. Тогда для функции Больцмана при всех t 0 справедливо неравенство Таким образом, рассматриваемое движение газа в сосуде сопровождается невозрастанием с течением времени величины H(t), что указывает на его необратимый характер.
2.2 Равновесная функция распределения и система уравнений Эйлера Точным решением уравнения Больцмана является функция распределения, называемая локально-максвелловской равновесной функцией, которая в размерных величинах имеет вид 2.3. Уравнения Навье–Стокса Для функции f (0) справедливо соотношение I(f (0), f (0) ) = 0, и она связана с f соотношениями Непосредственной подстановкой можно убедиться, что для локально-максвелловской функции распределения Функция f (0) называется также локально-равновесной функцией распределения.
Интегрирование уравнения Больцмана с весами 1,, 2 /2 в нулевом приближении, то есть когда f считается равной f (0), позволяет замкнуть систему (2.4)–(2.6) и получить систему уравнений Эйлера.
2.3 Уравнения Навье–Стокса В 1916–1917 годах С. Чепмен и Д. Энског предложили асимптотический метод решения уравнения Больцмана, позволяющий замкнуть систему (2.4)–(2.6) и получить систему уравнений первого приближения для описания течений вязкого теплопроводного газа систему уравнений Навье–Стокса [65, 73, 134].
Суть метода заключается в том, что решение приведенного к безразмерному виду уравнения (2.1) ищется в виде формального асимптотического ряда по степеням малого положительного параметра числа Кнудсена Kn, в виде где Здесь средняя длина свободного пробега частиц в невозмущенном потоке, L характерный размер области течения. В качестве нулевого приближения используется локально-максвелловская функция (2.8).
В первом приближении по числу Kn вычисления с помощью метода Чепмена–Энскога приводят к так называемой локально-навье– стоксовской функции распределения Величины N S и qN S были выписаны ранее (см. (1.34) и (1.35)).
Последовательно интегрируя уравнение Больцмана с сумматорными инвариантами 1,, 2 /2 в предположении, что f совпадает с fN S, получим систему уравнений Навье–Стокса, выписанную ранее в разделе 1.4. Процедура Чепмена–Энскога позволяет провести приближенный расчет коэффициентов вязкости и теплопроводности.
Для газа твердых сфер приближенный расчет этих коэффициентов приводит к выражениям где r0 – радиус частицы. Число P r в (2.12) оказывается равным 2/3, а сами коэффициенты зависят только от температуры, что согласуется с известными экспериментальными данными.
Используя следующее приближение в разложении функции распределения по числу Кнудсена, можно получить уравнения Барнетта. Эти уравнения включают в себя третьи пространственные производные, что вызывает существенные трудности при их численном решении и постановке граничных условий.
2.4. Уравнение Бхатнагара–Гросса–Крука 2.4 Уравнение Бхатнагара–Гросса–Крука В работе П. Бхатнагара, Е. Гросса и М. Крука в 1954 году [132] было предложено приближенное кинетическое уравнение вида то есть уравнение (2.1), в котором столкновительный интеграл I(f, f ) аппроксимировался с помощью выражения В настоящее время уравнение (2.13) называют модельным кинетическим уравнением Бхатнагара–Гросса–Крука (БГК). Примерно в то же самое время уравнения (2.13), (2.14) независимо были опубликованы П. Веландером [73]. Положительный параметр в правой части равенства (2.14) интерпретируется как характерное время релаксации функции f к локально-максвелловскому распределению f (0), определяемому формулой (2.8), и считается заданной функцией плотности и температуры. Величина совпадает по порядку величины со средним временем свободного пробега молекул в газе.
Макроскопические величины, входящие в формулу для вычисления, также являются квадратурами от f. Для модели БГК справедлив аналог H-теоремы Больцмана.
Применение метода Чепмена–Энскога к уравнению БГК также приводит к системе Навье–Стокса [65, 73, 134]. При этом коэффициент динамической вязкости µ и коэффициент теплопроводности вычисляются по формулам Из приведенных формул следует, что в БГК приближении число Прандтля равно единице.
В настоящее время разработаны усовершенствованные модели типа БГК приближения. В частности, предложена S-модель Шахова, которая позволяет ввести в рассмотрение реальное значение числа Прандтля. При этом вместо равновесной функции распределения в интеграле столкновений (2.14) выбирается функция распределения вида где - некоторая функция [108].
Имеются обобщения БГК приближения на случай, когда характерное время релаксации зависит от скорости частиц = (). Вариант релаксационного уравнения, учитывающий неравновесность по внутренним степеням свободы частиц [134], используется в главе 8 для построения газодинамических уравнений. Вариант БГК приближения для смеси газов и его использование для построения моментных уравнений обсуждается в главе 9.
Приведем определения основных величин, характеризующих хаотическое движение частиц в газе с функцией распределения f. Полученные выражения будут использоваться в дальнейшем.
Средняя тепловая скорость частиц вычисляется как Средняя относительная скорость частиц определяется в виде где cr = 1 2 величина относительной скорости двух сталкивающихся молекул, cr = [(cx1 cx2 )2 + (cy1 cy2 )2 + (cz1 cz2 )2 ]1/2, f и f2 соответствующие функции распределения.
В случае равновесной функции распределения f = f (0), определенной согласно (2.8), интегралы (2.16) и (2.17) вычисляются аналитически и значения соответствующих средних определяются как 2.5. Средние характеристики движения частиц Наиболее вероятная скорость частиц определяется "шириной"функции распределения и в этом случае составляет Средняя частота между столкновениями определяется как где - сечение столкновения молекул, или сечение рассеяния Для газа с максвелловской функцией распределения f1 = f1, f2 = = f2 в VHS приближении [134] где ref – величина сечения столкновения при температуре Tref.
Для газа твердых сфер ( = 0, = 0.5), (2.20) принимает вид Среднее время между столкновениями вычисляется через обратную частоту столкновений как Для газа твердых сфер В VHS приближении максвелловское время релаксации связано со средним временем между столкновениями [134] Для газа твердых сфер ( = 1/2) = 5/4.
Cредняя длина свободного пробега частиц определяется средним временем свободного пробега c и средней тепловой скоростью (2.16) в виде Для газа твердых сфер с максвелловской функцией распределения подстановка соотношений (2.18) и (2.23) приводит к выражению Приведем некоторые оценки характерных параметров для воздуха: при атмосферном давлении на уровне моря число частиц составляет n = 2.4 · 1025 1/м3, среднее расстояние между частицами r = n1/3 = 3 · 107 м, = 107 м, c = 2.5 · 1010 сек, характерная скорость < c >= 300 м/сек, = 1018 м2.
На высоте 300 км от поверхности Земли концентрация частиц n = /m0 1015 м3, среднее расстояние между молекулами r 3 · 105 м, средняя длина свободного пробега 103 м, c сек.
Из приведенных оценок наглядно видно, что значения средних величин в газе сильно меняются с изменением плотности частиц.
В частности, это касается соотношения средней длины свободного пробега и среднего расстояния между молекулами. Это поясняет тот факт, что масштабы пространственного V и временнного t осреднения, введенные в рассмотрение в п. 1.1 первой главы, могут существенно изменяться в зависимости от конкретной рассматриваемой задачи.
2.6 Коэффициенты переноса в равновесном Хаотическое перемещение частиц, рассматриваемое на микроскопическом уровне, сопровождается перераспределением их числа, а также переносом каждой частицей своего импульса и энергии. Тем самым на макроскопическом уровне описания движение молекул поКоэффициенты переноса в равновесном газе рождает три взаимосвязанных процесса переноса – это диффузия или самодиффузия, вязкость и теплопроводность.
Эти три транспортных процесса тесно связаны между собой. Во всех трех случаях потоки пропорциональны градиенту соответствующей величины. Согласно [182], эти процессы могут быть единообразно описаны с использованием основных понятий кинетической теории а именно, в терминах средней скорости частиц и средней длины свободного пробега.
Действительно, пусть длина свободного пробега много меньше характерного размера задачи, связанного с градиентами макроскопических величин плотности, скорости и давления. Рассмотрим перенос некоторой скалярной величины A через площадку единичной площади, расположенную перпендикулярно оси z, за единицу времени. Тогда нормальная составляющая потока A, протекающего за счет хаотического движения частиц через эту площадку в единицу времени, определяется как Здесь n – плотность частиц, v – средняя скорость хаотического движения частиц, – средняя длина свободного пробега. Коэффициент 1/3 выбирается из предположения, что все три координатных направления при хаотическом движении равноправны.
Связывая величину A с концентрацией A = n1 /n, найдем диффузионный поток через площадку Отсюда получаем выражение для коэффициента диффузии в виде Для вычисления переноса импульса через единичную площадку представим величину A как A = m0 u, где u макроскопическая скорость течения газа вдоль площадки. Тогда поток величины A, переносимый при случайном блуждании частиц, будет соответствовать компоненте тензора вязких напряжений Таким образом, получаем кинетическую оценку для коэффициента вязкости – формулу Максвелла Связывая значение A с тепловой энергией частицы A = m0 v 2 /2 и заменяя v 2 на наиболее вероятную скорость v 2 = c2 = 2RT, получим выражение для потока тепла через площадку То есть коэффициент теплопроводности имеет вид где cv = 3kB /m0 теплоемкость единицы массы при постоянном объеме для одноатомного газа. Полученные на основе простых кинетических оценок значения коэффициентов диффузии, вязкости и теплопроводности оказываются связанными между собой и пропорциональны длине свободного пробега и средней скорости хаотического движения частиц v. Определяя v как среднюю тепловую скорость (2.18), получим известные выражения для всех трех коэффициентов переноса при Sc = 1 и P r = 1, а также связь длины свободного пробега с коэффициентом вязкости, которая отличается от формулы Чепмена (см. раздел 3.4) численным коэффициентом порядка единицы.
Таким образом, упрощенное рассмотрение процессов переноса в равновесном газе позволяет получить качественно правильные выражения для коэффициентов диффузии, вязкости и теплопроводности. Все три процесса переноса представляются равноправными и все три явно присутствуют в КГД уравнениях. В систему Навье– Стокса входят только два из них перенос импульса и тепловой энергии, связанные с коэффициентами вязкости и теплопроводности.
2.7. Численное моделирование течений разреженного газа разреженного газа 2.7.1 Общие замечания Характеристикой степени разреженности газодинамического течения является число Кнудсена Kn = /L, представляющее собой отношение средней длины свободного пробега молекул к характерному линейному размеру задачи L. Обычно газ рассматривается как плотный, если Kn 0 (на практике Kn < 0.01). Условия, при которых Kn (на практике Kn > 10), характерны для свободномолекулярных течений, когда столкновения между частицами практически отсутствуют. При промежуточных числах Kn газ считается разреженным.
Методы расчета свободномолекулярных режимов к настоящему времени достаточно хорошо разработаны. Для этих задач столкновениями частиц между собой можно пренебречь и учитывать только взаимодействие частиц со стенками. Распределение частиц по скоростям с большой точностью можно считать известным, например, равновесным с функцией распределения f (0). Основной проблемой при этом является описание взаимодействия частиц со стенками.
Этот процесс можно приближенно описать с помощью коэффициента аккомодации, который обозначается через. Простейшими моделями здесь являются: полная аккомодация частиц на стенке, или так называемое диффузное отражение, которое соответствует значениям = 1, и модель зеркального отражения, в которой полагается Под умеренно-разреженным течением газа понимают такие течения, когда число Кнудсена лежит в диапазоне порядка от 0.01 до 0.1, в зависимости от рассматриваемой задачи. Течения умеренноразреженного газа представляют собой область, находящуюся на границе применимости кинетического подхода и подхода, связанного с решением моментных уравнений. Расчет таких течений методами кинетической теории требует неоправданно больших вычислительных ресурсов, что обусловлено высокой плотностью газа. В то же время уравнения Навье–Стокса, полученные в приближении Kn 0, теряют свою точность при анализе указанных режимов.
Для расчета течений умеренно-разреженного газа в рамках моментных уравнений возникает необходимость учета отклонения от режима сплошной среды, в первую очередь вблизи обтекаемой поверхности. Для этого используются специальные граничные условия.
Для всех чисел Кнудсена, как бы малы они ни были, вблизи стенки существует слой газа, толщина которого имеет порядок средней длины свободного пробега молекул так называемый слой Кнудсена. Для того, чтобы в рамках макроскопических уравнений учесть влияние этого слоя на поле течения, вводятся граничные условия, представляющие собой условия скольжения для скорости и скачка для температуры.
Первый вариант таких условий был выписан Максвеллом в предположении диффузного характера отражения молекул от стенки1. В настоящее время в литературе имеется много вариантов таких условий, которые приведены, например, в книгах [65, 78, 134]. Все они имеют одинаковую структуру и отличаются между собой лишь численными коэффициентами порядка единицы. Приведем здесь условия в форме Смолуховского, где по сравнению с классическими условиями Максвелла учтено значение коэффициентов аккомодации для скорости u и энергии e, которые могут быть различными, а также учтено влияние градиента температуры вдоль стенки:
где un и us нормальная составляющая скорости газа вблизи стенки и скорость скольжения вдоль стенки, n и s координаты вдоль внешней нормали к стенке и вдоль стенки, Ts температура газа вблизи стенки и Tw температура стенки.
Для большинства материалов в условиях сверхзвукового обтекания коэффициенты аккомодации для скорости и энергии можно Maxwell J.C. Philos. Trans. R. Soc., London, 1879, v.170, p. 2.7. Численное моделирование течений разреженного газа полагать одинаковыми и близкими к единице. В формуле для скорости скольжения второе слагаемое начинает играть заметную роль лишь при числах Кнудсена, приближающихся к единице [19].
Течения газа в диапазоне чисел Кнудсена 0.1 – 10 представляют собой существенную сложность для аналитического исследования и численного моделирования, так как в этом диапазоне не удается выделить малый параметр по числу Кнудсена типа Kn или 1/Kn.
В этом диапазоне применяются методы кинетической теории. Численный анализ течений проводится либо путем непосредственного решения уравнения Больцмана или его упрощенных вариантов, либо на основе методов прямого численного моделирования – методов Монте–Карло, или DSMC методов [12], [134].
Сложности использования этих подходов связаны с большими затратами машинного времени при моделировании этапа столкновений частиц и большой размерностью задачи в целом, которая рассматривается в пространстве 7 измерений (x,, t). Дополнительная сложность вызвана необходимостью вычисления осредненных характеристик для получения измеряемых газодинамических величин - скорости, плотности и давления. В рамках DSMC подхода дополнительные трудности представляют расчеты нестационарных течений и течений с малыми, то есть дозвуковыми скоростями, поскольку такие расчеты сопровождаются значительными флуктуациями вычисляемой плотности частиц.
2.7.2 Метод Монте-Карло В широком смысле слова методами Монте–Карло называют основанные на моделировании случайных величин методы решения различных задач из таких областей, как статистическая физика, вычислительная математика, теория игр, математическая экономика и многие другие [56].
Метод прямого численного моделирования для расчета течений разреженного газа (DSMC Direct Simulation Monte–Carlo, или ПММК прямое моделирование Монте–Карло) был разработан в 1960-х годах австралийским ученым Г.А. Бердом (G.A. Bird) и усовершенствован впоследствии [12, 134].
В газовой динамике нашел применение вариант метода Монте– Карло, основанный на моделировании реального течения газа посредством относительно небольшого числа молекул. То есть проводится численный эксперимент, в котором прослеживается история ограниченного числа частиц, каждая из которых является представителем большого числа W реальных молекул. Величина W называется "весовым множителем"(weighting factor).
Для каждой из молекул запоминаются ее координаты, скорость и энергия. По этим величинам путем осреднения по всем частицам определяются газодинамические параметры течения.
Для стационарных задач расчет начинается с задания некоторого достаточно произвольного распределения частиц в расчетной области, которое с течением времени эволюционирует к своему равновесному состоянию. Далее перечислены основные этапы DSMC метода.
Дискретизация и моделирование движения частиц Область течения разбивается на пространственные ячейки, причем такие, чтобы изменение газодинамических параметров течения в каждой ячейке было малым. Размер ячейки имеет порядок средней длины свободного пробега. Для эффективности счета число частиц в каждой пространственной ячейке не должно сильно различаться и составлять порядка нескольких десятков.
Моделирование физического движения молекул проводится посредством дискретных шагов по времени t, малых по сравнению со средним временем между столкновениями молекул, t < c. Движение молекул и межмолекулярные столкновения на временном интервале моделируются последовательно. На каждом шаге по времени t движение частиц разбивается на два этапа и описывается в рамках кинетической модели, которая представляет собой циклическиповторяющийся процесс бесстолкновительного разлета частиц и последующих столкновений, которые рассматриваются как мгновенные. Эта модель соответствует двум этапам расчета.
• Перемещение На первом этапе все молекулы перемещаются на расстояние, 2.8. Кинетически-согласованные разностные схемы определяемое их скоростями t. Учитываются пересечения молекулами поверхностей твердых тел, линий и плоскостей симметрии и границ течения. При наличии потока внутрь области на соответствующих границах генерируются новые молекулы. Если молекула покидает область расчета, то она исчезает.
• Столкновения На втором этапе моделируются столкновения между молекулами с последующей коррекцией молекулярных скоростей.
Выбор очередной сталкивающейся пары частиц проводится в пределах одной ячейки и производится на основе данных генератора случайных чисел. Предполагается, что сталкиваются только те частицы, которые находятся в одной пространственной ячейке.
Важной частью метода прямого моделирования является вычисление числа столкновений. Частота столкновений определяется свойствами реального газа, для которого решается задача, и именно эта величина определяет диссипативные характеристики течения коэффициенты вязкости µ и теплопроводности моделируемого газа.
Расчет макроскопических характеристик Для вычисления макроскопических параметров газа, u, p, T запоминаются и аккумулируются данные для всех молекул. Затем происходит дополнительное осреднение по последовательности расчетов, чтобы сгладить статистические флуктуации, возникающие в процессе вычислений.
2.8 Разностная аппроксимация уравнения Больцмана и кинетически-согласованные разностные схемы Численное решение уравнения Больцмана (2.1) представляет собой существенные сложности, связанные с большой размерностью задаГл.2. Элементы кинетической теории газов чи и проблемами аппроксимации и вычисления интеграла столкновений. Развитые для этой задачи численные алгоритмы приведены, в частности, в [9, 124].
Левая часть уравнения (2.1) имеет вид уравнения переноса функции f со скоростью, и может быть аппроксимирована с помощью разностной схемы первого порядка точности с направленными разностями схемы Куранта–Исааксона–Рисса. При согласованном выборе шагов по времени и по пространству h = ||t эта разностная схема описывает перенос заданного возмущения функции f без искажений [87, 88].
Приведем вид соответствующей разностной схемы для одномерного пространственного течения в отсутствие внешних сил на равномерной пространственной сетке с шагом h = xi+1 xi и шагом по времени t = tj+1 tj Здесь k I(f, f )h,i разностный аналог интеграла столкновений.
Разностная схема (2.34) тождественным образом преобразуется к эквивалентному виду Пусть для разностного аналога интеграла столкновений I(f, f )h,i выполняются условия ортогональности сумматорным инвариантам h(). Разностную схему (2.35) можно осреднить в пространстве скоростей, предполагая, что f = f (0). При этом все интегралы удается вычислить аналитически и сразу построить разностную схему для 2.8. Кинетически-согласованные разностные схемы газодинамических величин – плотности, скорости u и энергии E [30, 51, 105, 106, 114].
Полученная схема оказывается довольно громоздкой, поскольку она включает в себя интегралы ошибок, возникающие при осреднении модулей скоростей k. Тем не менее эта схема оказалась весьма эффективной для решения уравнений Эйлера и была использована при проведении расчетов как одномерных, так и некоторых двумерных сверхзвуковых газодинамических течений. Позднее близкие по структуре разностные уравнения были получены в работе [144].
Заменим в разностной схеме (2.35) коэффициент перед последним слагаемым, полагая || c, где c – скорость звука. Вводя в рассмотрение время = h/2c, характеризующее скорость пересечения частицей разностной ячейки, получим Схема (2.35) примет вид Выражение (2.37 ) в случае f = f (0) также допускает аналитическое осреднение с сумматорными инвариантами h() и приводит к построению более изящной по сравнению с предыдущим случаем системы разностных уравнений для описания макроскопических величин в газе. Полученные таким образом схемы были названы кинетически-согласованными разностными схемами КСРС.
С использованием принятых в [92] обозначений для центральноразностных производных первого и второго порядка по пространству, перепишем разностную схему (2.37) в виде Здесь f = (fi+1 fi1 )/h – центральная разностная производная первого порядка, fx = (fi+1 fi )/h и fx = (fi fi1 )/h – левая и правая разностные производные, fxx = (fx fx )/h – центральная разностная производная второго порядка.
Применяя введенные обозначения, запишем КСРС для случая плоского одномерного течения в виде Модифицированные позднее КСРС оказались очень эффективными при численном моделировании широкого круга газодинамических течений [6–8, 24, 25, 30, 31, 50–52, 105, 106].
Дифференциальные аналоги этих схем легли в основу первых вариантов КГД уравнений.
Глава Квазигазодинамические уравнения В этой главе приведено два варианта построения квазигазодинамической системы уравнений, которые позволяют получить конкретный вид выписанных ранее без вывода выражений для векторов плотности потока массы jm, тензора вязких напряжений и вектора теплового потока q.
Первый вариант основан на использовании достаточно простой кинетической модели движения частиц. Эта модель представляет собой циклически-повторяющийся процесс бесстолкновительного разлета и последующих столкновений частиц с установлением максвелловского равновесия. Изложенный вывод КГД уравнений устанавливает связь параметра сглаживания с максвелловским временем релаксации. Построенные таким образом КГД уравнения представляются затем в виде законов сохранения. При этом выражения для jm, и q записываются как соответствующие выражения для уравнений Навье–Стокса с добавками, пропорциональными.
Показано, что для стационарных течений эти дополнительные слагаемые имеют второй порядок малости по. На основе этого рассмотрения получен вид диссипативных коэффициентов и выписаны их обобщения.
Второй вариант базируется на рассмотрении законов сохранения массы, импульса и энергии для малого, но конечного неподвижного объема, где при выполнении осреднения учитывается нестационарность газодинамических величин во времени. Этот подход позволяет построить КГД систему для течений в присутствии внешних сил и источников тепла. Здесь же приведен вывод уравнения баланса энтропии. Это уравнение демонстрирует диссипативный характер возникающих дополнительных слагаемых.
Изложение этих результатов приведено в соответствии с работами [43, 49–53, 114, 115, 146, 158].
3.1 Регуляризованное кинетическое уравнение Решение интегро-дифференциального уравнения Больцмана является сложной задачей, поэтому к настоящему времени развиты упрощенные кинетических модели, которые позволяют находить приближенные решения отдельных задач. Многие модели основаны на идее расщепления задачи по физическим процессам, которая иначе называется принцип суммарной аппроксимации. Именно такой подход лежит в основе DSMC метода, описанного в п. 2.7.2.
Рассмотрим кинетическую модель, которая использовалась в [9] для численного моделирования разреженных течений в рамках уравнения Больцмана, и на основе которой впервые была получена система квазигазодинамических уравнений.
этап мгновенной максвеллизации. Схематическое изобраРис. 3.1. Схема, поясняющая кинети- жение этой модели приведено ческую модель в [51, 52, 105, 106] и здесь не обсуждается.
Пусть в некоторый момент времени t = t1 функция распределения имеет локально-максвелловский вид:
Затем в течение времени t происходит бесстолкновительный разРегуляризованное кинетическое уравнение лет молекул, который описывается уравнением Больцмана для свободномолекулярного течения Это уравнение представляет собой линейное уравнение переноса и имеет точное решение где f0 = f0 (x, ) известная функция распределения в момент времени t = 0.
Далее, в момент времени t = t2 = t1 + t функция распределения f мгновенно вновь становится локально-максвелловской (3.1), но уже с новыми значениями макроскопических параметров, u и T.
Мгновенная максвеллизация имитирует этап столкновения молекул, который в уравнении Больцмана описывается интегралом столкновений I(f, f ). Затем оба этапа циклически повторяются.
Считая время бесстолкновительного разлета t достаточно малым, разложим функцию распределения (3.2) в ряд Тейлора, полагая функцию f0 максвелловской f0 = f (0). Ограничиваясь членами второго порядка малости по t, получим:
Параметр разложения t · при больших значениях скоростей нельзя считать малым, и отбрасываемые при таком разложении члены могут оказаться существенными. Однако в данном случае их вкладом при || 2RT можно пренебречь, так как сама функция распределения f (0) и все ее производные экспоненциально убывают с ростом.
В равенстве (3.3) перенесем все слагаемые в левую часть и разделим обе части уравнения на t, заменяя разностную производную по времени в первом слагаемом на дифференциальную. В результате получим уравнение для функции распределения в виде Здесь в правую часть полученного уравнения добавлен интеграл столкновений в БГК приближении, обеспечивающий релаксацию функции распределения на новом временном слое к максвелловскому распределению.
Отождествляя временной интервал t/2 с максвелловским временем релаксации, получим окончательный вид регуляризованного кинетического уравнения Заметим, что это уравнение представляет собой дифференциальный аналог разностной схемы Куранта–Иссаксона–Рисса для БГК уравнения, выписанной в виде (2.38) для случая = const.
При осреднении уравнения (3.4) с сумматорными инвариантами, как будет показано далее, получается система КГД уравнений для Уравнение (3.4) формально можно выписать исходя из уравнения Больцмана в БГК приближении:
заменяя функцию распределения в конвективном слагаемом уравнения (3.5) на ее приближенное значение вида Заметим, что формальная замена функции распределения f в конвективном слагаемом (3.5) на функцию распределения Навье– Стокса (2.11) позволяет, после осреднения с сумматорными инвариантами, получить систему уравнений Навье–Стокса.
В [158] было показано, что соотношения (2.11) и (3.6) могут быть представлены в идентичном виде как где P 3 (i ) полином третьей степени. При этом коэффициенты полиномов для обеих функции распределения оказываются близкими.
3.2. Кинетический вывод КГД уравнений Регуляризованное кинетическое уравнение (3.4) построено здесь с использованием большого числа допущений, однако оно имеет целый ряд свойств, сходных со свойствами БГК уравнения. В частности, для стационарных течений показано, что если функция распределения f удовлетворяет БГК уравнению (3.5), то она удовлетворяет и уравнению (3.4) с точностью O( 2 ). И наоборот, если функция распределения f удовлетворяет регуляризованному кинетическому уравнению, то с точностью O( 2 ) она удовлетворяет и стационарному уравнению БГК. Таким образом, решения этих уравнений должны быть близки. Для уравнения (3.4) доказан аналог H-теоремы Больцмана [114].
На основе регуляризованного кинетического уравнения (3.4) удается построить КГД систему уравнений и родственные ей моментные уравнения для описания слабо-неравновесных газодинамических течений (главы 8 и 9).
3.2 Кинетический вывод КГД уравнений Приведем способ построения квазигазодинамической системы уравнений на основе регуляризованного кинетического уравнения, выписанного в предыдущем параграфе.
Последовательно интегрируя кинетическое уравнение (3.4) с весами 1,, 2 /2 и используя выражения получим систему КГД уравнений в следующем виде:
Замыкание вида (3.7), основанное на локально-максвелловской функции распределения, автоматически приводит к построению системы моментных уравнений для газа твердых сфер, то есть для = 5/3. Обобщение на случай произвольных молекул проводится путем выделения слагаемых с удельной внутренней энергией в предположении = 5/3, и затем обобщения этого выражения на случай произвольного значения.
Система (3.8)–(3.10), выведенная для случая произвольных криволинейных эйлеровых координат, замыкается уравнениями состояния идеального политропного газа При записи уравнения энергии в последнем слагаемом имеются члены, не содержащие операций дивергенции от скорости. Из этого слагаемого можно сразу выделить слагаемые в форме теплового 3.3. КГД уравнения в форме законов сохранения потока Навье–Стокса qN S и ввести число Прандтля P r = Пример построения КГД уравнений для плоского одномерного течения газа детально выписан в Приложении A, где наглядно прослежены все использованные преобразования.
На основе осреднения регуляризованного кинетического уравнения с учетом внешней силы F Ю.В. Шеретов построил КГД уравнения для описания течений газа в поле внешних сил [113, 114].
Левые части системы (3.8) – (3.10) представляют собой левые части уравнений Эйлера. Правые части имеют дивергентный вид и пропорциональны малому параметру. При этом построенные уравнения (3.8) – (3.10) не имеют вид законов сохранения, то есть в них в явном виде не выделены тензор вязких напряжений и векторы плотности потока массы и теплового потока. Диссипативные слагаемые входящие в уравнение Навье–Стокса также явно не выделены.
Преобразование построенной выше КГД системы к форме законов сохранения проведена в следующем разделе.
Представим КГД систему в виде законов сохранения (1.54)–(1.56), которые для случая F = 0 имеют следующий вид:
Уравнение неразрывности и вектор плотности потока массы Сравнивая уравнения переноса массы (3.8) и (3.13), найдем вектор плотности потока массы Или, обозначив добавку к скорости через представим jm в виде, предложенном ранее в первой главе, и удобном для дальнейшего анализа уравнений Уравнение импульса и тензор вязких напряжений Сопоставляя уравнения импульса (3.9) и (3.14), найдем вид тензора вязких напряжений. Для этого представим div(jm u) = div(u u) div( div(u u) u) div( p u).
Тогда (3.9) можно переписать в виде где использовано тождество 3.3. КГД уравнения в форме законов сохранения Из (3.17) и (3.14) следует, что тензор вязких напряжений имеет вид Представим этот тензор в виде суммы тензора вязких напряжений Навье–Стокса и некоторой добавки. Для определения вида этой добавки воспользуемся тождествами Тогда примет вид Преобразуем полученное выражение, прибавляя и вычитая величины Ip div u и (2/3) pI div u Группируя слагаемые, получаем:
Сравнивая (3.21) с видом тензора Навье–Стокса (3.19), увидим, что Из первого соотношения сразу следует, что имеет смысл максвелловского времени релаксации [134].
Из полученной формулы для коэффициента объемной (второй) вязкости следует, что коэффициент неотрицателен и связан с наличием внутренних степеней свободы молекулы, что соответствует теоретическим представлениям, изложенным в [18, 57, 72, 78]. Действительно, для одноатомного газа = 5/3 и = 0, в противном случае, при наличии колебательных или вращательных степеней свободы молекулы, < 5/3 и > 0.
В результате этих преобразований получаем, что тензор вязких напряжений в КГД системе уравнений представляется в виде тензора вязких напряжений Навье–Стокса N S, (3.19) и добавки, которая пропорциональна и обращается в нуль при u = Уравнение полной энергии и вектор теплового потока В заключение рассмотрим третью пару уравнений (3.10), (3.15), из которых определим вид векторов A и q. Для этого запишем (3.10) в 3.3. КГД уравнения в форме законов сохранения виде а (3.15) перепишем в виде Сравнивая последние два выражения, найдем здесь для w использовано представление (3.16). Приводя подобные члены, получим где мы воспользовались тождеством Приведем подобные слагаемые Применяя тождество получим Сгруппируем слагаемые 3.3. КГД уравнения в форме законов сохранения Применив тождества получим Из полученного выражения выделим тензор в виде (3.21) В теории Навье–Стокса [72] представляет собой работу поверхностных сил давления и внутреннего вязкого трения в единицу времени. По аналогии положим Остальные слагаемые в (3.25) представим как тепловой поток и обозначим через Выделим в тепловом потоке q часть, связанную с потоком Навье– Стокса Примем в качестве уравнений состояния уравнения состояния идеального политропного газа представим формулу для теплового потока в виде Сопоставляя формулу для q и qN S, получим, что коэффициент теплопроводности в КГД модели определяется как и выражение для теплового потока q в КГД уравнениях имеет вид Таким образом, полученная на основе кинетического уравнения (3.4) КГД система вида (3.8)–(3.10) представлена в виде законов сохранения (3.13)–(3.15). При этом вектор плотности потока массы, тензор вязких напряжений и вектор теплового потока представляются в виде суммы соответствующих величин в форме Навье– Стокса и малых добавок существенно нелинейного вида, пропорциональных коэффициенту (3.16, 3.23, 3.26).
3.4. Коэффициенты диссипации 3.4 Коэффициенты диссипации 3.4.1 Формулы для диссипативных коэффициентов и На основании кинетического вывода КГД уравнений и представления их в виде законов сохранения удается получить вид диссипативных коэффициентов µ, и. Эти коэффициенты связаны между собой через параметр и выражаются как Именно в таком виде коэффициенты вязкости µ и теплопроводности получаются при выводе системы уравнений Навье–Стокса методами Чепмена–Энскога из уравнения БГК. Коэффициент объемной вязкости, совпадающий с формулой (3.27), был получен в [193] на основе БГК приближения для немоноатомного газа, обладающего вращательными степенями свободы. Действительно, приведенная в [193] формула для этого коэффициента имеет вид где i число вращательных степеней свободы. Выражая общее число степеней свободы n = i + 3 через как n = i + 3 = 2/( 1), получим коэффициент объемной вязкости вида (3.27).
Вводя число Прандтля P r = 1, запишем коэффициент теплопроводности в виде В п. 1.5.3 с точностью до численного коэффициента была найдена связь релаксационного параметра и коэффициента динамической вязкости µ в виде где Sc число Шмидта. Для газа его величина близка к единице.
Уточнение коэффициента объемной вязкости для газа, обладающего вращательными степенями свободы, приведено в следующем разделе.
3.4.2 Коэффициент объемной вязкости Диссипативные эффекты, связанные с релаксацией внутренних степеней свободы в газе, могут оказывать заметное влияние на течения с ударными волнами и на быстропеременные по времени процессы.
В гидродинамических моделях эти диссипативные эффекты наиболее просто описываются в терминах так называемой второй, или объемной, вязкости. Это описание является достаточно приближенным, и его использование ограничено течениями, в которых характерные времена релаксации внутренних степеней свободы малы по сравнению с характерными гидродинамическими временами задачи [57, 72].
Как уже отмечалось ранее, коэффициент объемной вязкости существенно связан с наличием внутренних степеней свободы молекул и обращается в нуль для одноатомных газов. Между тем величина этого коэффициента в ряде случаев может быть сопоставима с величиной коэффициента динамической вязкости µ.
Способы вычисления коэффициентов второй вязкости и получающиеся при этом выражения в общем случае достаточно сложны (см., например, [18, 78, 102]).
Относительно простые выражения для коэффициента второй вязкости имеются только в приближении, когда внутренние степени свободы в газе являются вращательными степенями свободы.
Пример такого выражения для приведен в предыдущем разделе.
Более точное выражение для коэффициента второй вязкости, зависящей от вращательных степеней свободы, предложено в [57] в виде:
здесь rot доля внутренней энергии, содержащейся во вращательных степенях свободы.
Покажем, что это выражение отличается от полученной на основе КГД уравнений формулы для (3.27) только численным коэффициентом. Для этого выразим rot через :
3.4. Коэффициенты диссипации где i число вращательных степеней свободы, n общее число степеней свободы. Если принять во внимание, что n = i + 3, где 3 число поступательных степеней свободы, то формулу для rot можно переписать в виде Учитывая связь числа степеней свободы с показателем адиабаты получим, что Время релаксации вращательных степеней свободы rot выражается следующем образом:
где c среднее время свободного пробега молекул; Zrot коэффициент обмена энергий между вращательными и поступательными степенями свободы, который равен среднему числу межмолекулярных столкновений, необходимых для релаксации к равновесному состоянию вращательной моды. Выражения для Zrot приведены, например, в [57, 134]. Эти выражения для Zrot использовались для расчетов неравновесных течений в струях CO2 [172] и N2 [159], а также для численного моделирования структуры ударной волны в азоте [55]. Например, для азота где Z = 23, T = 91.5K, если температура газа в невозмущенном потоке составляет 273K [57]. В этом случае значение Zrot изменяется в пределах 5 16 в области температур от 300 до 6000 K (рис. 3.2).
Выразим c через динамическую вязкость газа µ. Cредняя длина свободного пробега частиц в газе определяется средним временем свободного пробега частиц и средней тепловой скоростью < c > как Величина средней тепловой скорости частиц в равновесном газе составляет (п. 2.5) < c >= 8RT /. Таким образом, Средняя длина свободного пробега может быть связана с вязкостью газа как где формула Чепмена [3] или формула Берда [134]. Здесь показатель степени в законе вязкости µ T.
Соотношения (3.31) и (3.32) позволяют выразить c через динамическую вязкость газа µ в виде Подставив выражения для cv, rot (3.29) и c (3.33) в (3.28), запишем исходную формулу для второй вязкости (3.28) в виде Последняя формула (3.34) с точностью до численного множителя совпадает с формулой для коэффициента второй вязкости (3.27), полученной ранее при анализе КГД уравнений. На рис. 3.2 приведены значения B (пунктир) и Zrot (сплошная линия) в зависимости 3.5. Система Навье–Стокса как асимптотика КГД системы Рис. 3.2. Коэффициенты B (пунктир) и Zrot (сплошная линия) для азота от температуры для азота. Видно, что коэффициент B изменяется в пределах от 1 до 4.
Определяющую роль коэффициента объемной вязкости демонстрируют приведенные в приложении B расчеты структуры ударной волны в азоте.
3.5 Система Навье–Стокса как асимптотика Диссипативные слагаемые N S и qN S, входящие в КГД уравнения, при условии достаточной гладкости соответствующих производных, имеют порядок малости O(µ) O( ), или, в безразмерном виде, O(Kn). Добавки к тензору вязких напряжений и векторам плотности потока массы и теплового потока, выписанные в п. 3.3, имеют такой же порядок малости. Покажем, что для стационарных течений эти КГД добавки имеют уже не первый, а второй порядок малости по, то есть O( 2 ).