Содержание к диссертации
Введение
1 Кинетически-соласованные схемы для одномерной газовой динамики 14
1.1 Кинетические схемы для уравнения Бюргерса 15
1.2 Кинетически согласованные схемы для уравнения Эйлера . 24
1.3 Кинетические согласованные схемы как схемы расщепления вектора потока 37
1.4 Кинетические схемы повышенного порядка 52
2 Кинетические схемы на неструктурированных сетках 62
2.1 Кинетически-согласованная аппроксимация на нерегулярной сетке 64
2.2 Кинетически-согласованные схемы повышенного порядка точности на нерегулярных сетках 72
2.3 Кинетические схемы для описания вязкого теплопроводного газа 77
3 Вывод и анализ квазигазодинамической системы уравнений 84
3.1 Модельное кинетическое уравнение 85
3.2 Квазигазодинамическая система 87
3.3 Асимптотика квазигазодинамической системы 91
3.4 Разностная аппроксимация и пример численного расчета. 94
Заключение 101
Литература 102
- Кинетически согласованные схемы для уравнения Эйлера
- Кинетические схемы повышенного порядка
- Кинетически-согласованные схемы повышенного порядка точности на нерегулярных сетках
- Квазигазодинамическая система
Введение к работе
Многие современные проблемы, стоящие перед наукой и техникой в той или иной мере связаны с решением уравнений газовой динамики. Полное решение большинства актуальных задач аэродинамики и, прежде всего, получение количественных характеристик предполагает использование новейших численных методов и их реализацию на быстродействующих вычислительных машинах. В последнее время математическое моделирование явлений и процессов является эффективным средством теоретического анализа задач, выдвигаемых наукой и техникой [1]. Вычислительный эксперимент является во многих случаях практически единственным способом решения сложных систем нелинейных дифференциальных уравнений, описывающих процессы движения жидкости и газа, уравнений Эйлера и Навье-Стокса. То есть, на основе математической модели при помощи непосредственного численного решения соответствующих уравнений количественно определяется поведение газодинамических течений в тех или иных условиях.
Общих методов решения газодинамических задач в настоящее время не существует. Особенности уравнений газовой динамики, такие как нелинейность и возможность появления разрывных решений, зачастую делают именно численные методы исследования наиболее предпочтительными и эффективными. В тоже время, именно нелинейность порождает многие явления, с которыми приходится считаться в практически важных случа-
ях — ударные волны, погранслойные, отрывные, осциллирующие течения, наконец, турбулентность. В последнее время интенсивно развивается новое направление в науке — вычислительная гидродинамика.
Одним из численных алгоритмов, используемых для моделирования течений вязкого теплопроводного газа, являются кинетически-согласованные разностные схемы (КСРС), полученные в 1983 году [2]. В отличии от традиционных подходов, в которых вычислительные алгоритмы для расчета газодинамических течений строятся на основе полной или определенным образом упрощенной систем уравнений Эйлера и Навье-Стокса, кинетически-согласованные разностные схемы учитывают тот факт, что сами уравнения гидродинамики могут быть получены как следствие более сложного уравнения переноса для одночастичной функции распределения [3]. В свою очередь КСРС получаются путем осреднения по скоростям молекул с сумматорными инвариантами дискретных моделей для одночастичной функции распределения. По сравнению с традиционными алгоритмами, в КСРС меняется порядок процедуры осреднения и дискретизации в последовательности, начинающейся с дифференциального уравнения Больцмана. То есть, в отличии от других алгоритмов решения уравнений газовой динамики, в данном случае вначале идет процедура разностной дискретизации, а потом уже осреднение разностной функции распределения [4, 5].
Связь уравнений газовой динамики с уравнением Больцмана хорошо известна, поэтому попытки построения вычислительных алгоритмов, обращаясь непосредственно к кинетическим моделям, осуществлялись и ранее. Например, в работе [6] была предложена феноменологическая модель,
описывающая перенос частиц в газе как последовательность процессов бес-столкновительного разлета молекул с последующей мгновенной максвел-лизацией и построена разностная схема решения задач газовой динамики в лагранжевой системе координат. Позднее появилась серия публикаций, в которых макроскопические газодинамические параметры определялись путем осреднения функции распределения после решения уравнения переноса [7, 8]. В работе [9] прослежена связь между одной из полученных таким образом схем и численным решением уравнения Больцмана.
Однако, хотя данные подходы показывали их принципиальную применимость для нахождения газодинамических параметров, но качество расчетов, с учетом вычислительных затрат, уступало уже существующим традиционным методам решения задач о течении идеального газа. В то же время, полученный в работе [2] алгоритм сразу же показал свою конкурентноспособность с традиционными подходами. Проведенные ранее на основе КСРС расчеты газодинамических течений показали перспективность их использования. Следует отметить, что аналогичные схемы были получены в 1986 году в работе [10, 11].
Уравнения Навье-Стокса являются теоретической основой для описания течений реального газа. Первоначально эти уравнения были получены феноменологическим путем, опираясь на известный экспериментальный закон о пропорциональности силы трения между движущимися слоями газа производной скорости. Впоследствии эти уравнения были получены, используя асимптотическое разложение Чемпена-Энскога решения кинетического уравнения Больцмана для одночастичной функции распределения [3]. Однако, известны неоднократно предпринимаемые попытки получить
.7
описание поведения вязкого газа уравнениями, отличными от традиционных уравнений Навье-Стокса [12, 13]. В работе [14, 15] впервые была получена обобщенная квазигазодинамическая система уравнений, объединяющая газодинамический и кинетический подход. В дальнейшем проводилось много работ в этом направлении, что отражено в монографии [16].
В работах [17, 18, 19] на основе кинетических моделей были получены квазигазодинамические дифференциальные уравнения и их дискретизация — кинетически-согласованные разностные схемы, использованные затем для численного моделирования течений вязкого теплопроводного газа. Данные уравнения отличаются от обычных уравнений Эйлера и Навье-Стокса наличием в правой части дополнительной вязкости специального вида, являющейся схемной искусственной вязкостью, влияющей на устойчивость разностной схемы.
Кинетически-согласованный подход обеспечивает естественное обобщение разностной схемы для линейного уравнения переноса на систему нелинейных газодинамических уравнений. С другой стороны, аппроксимацию уравнения переноса схемой с направленными разностями можно интерпретировать как физическую модель переноса кусочно-постоянной функции распределения. Получающиеся после осреднения разностные макроуравнения в дифференциальном приближении можно рассматривать как обобщение уравнений Эйлера и Навье-Стокса — квазигазодинамическая система уравнений [4, 5].
Кинетически-согласованные разностные схемы показали свою эффективность при решении задач газовой динамики. Целью диссертации является рассмотрение кинетически-согласованных разностных схем с точки
зрения классических методов и подходов, используемых при построении аппроксимаций уравнений газовой динамики. И на этом пути выявить возможные варианты модификации и улучшения существующих на сегодняшний день кинетических схем, а также возможность построения новых высокоэффективных алгоритмов на основе кинетического подхода.
Ниже приводится краткое содержание работы по главам. В первой главе рассматриваются базовые принципы построения кине-тически-согласованных разностных схем для невязкого уравнения Бюргерса в качестве модельной задачи и уравнений одномерной газовой динамики на основе использование модели расщепления по физическим процессам уравнения Больцмана для одночастичной функции распределения молекул по скоростям. Отмечено, что вид разностной схемы для макроскопических уравнений зависит не только от способа аппроксимации уравнения переноса для функции распределения, но и от вида функции распределения. Далее было показано, что КСРС можно интерпретировать как схему из класса схем вектора расщепления потока для аппроксимации уравнений Эйлера и предложена новая схема из этого класса. Затем на основе этого расщепления были построены кинетические схемы повышенного порядка точности, не прибегая к рассмотрению кинетического уравнения. Все построения иллюстрируются расчетами тестовых задач.
Во второй главе рассматриваются вопросы, связанные с обобщением полученных кинетических схем для аппроксимации конвективной части уравнений Навье-Стокса на неструктурированные треугольные сетки. Приведен вывод кинетически согласованной схемы на основе метода конечных объемов путем обобщения модели «бесстолкновительный разлет-
мгновенная максвеллизация» на произвольные треугольные сетки. Были сконструированы кинетические схемы повышенного порядка точности по пространственной переменной на треугольных сетках, аппроксимирующие невязкие уравнения Эйлера. Проведена аппроксимация диффузионных членов уравнения Навье-Стокса с использованием метода конечных объемов. В качестве примера использования кинетических схем для моделирования сложных задач был сделан расчет задачи обтекание двумерного аэродинамического профиля потоком вязкого газа с числом Маха набегающего потока равным 0.8 с использованием неструктурированных адаптивных сеток.
В третьей главе приводится вывод полудискретной квазигазод-намической системы уравнений, исходя из модели "бесстолкновительный разлет-максвеллизация". Далее проводится анализ полученной квазигазодинамической системы уравнений (КГУ) и показывается ее отличие от системы уравнений Навье-Стокса наличием членов второго порядка малости относительно времени бесстолкновительного разлета. Далее система КГУ аппроксимируется с использованием смешанного метода: метод конечных объемов и кинетическая схема для аппроксимации конвективного переноса, метод конечных элементов для вычисления вязких членов. Для сравнения поведения разностного решения уравнений Навье-Стокса и КГУ решается задача об обтекание сферы при числе Маха равном 0.1 и числе Рейнольдса равном 25.
В заключении приведены основные результаты диссертации. Цели и задачи диссертационной работы:
Основной целью данной диссертационной работы является обосно-
вание и улучшение существующих на сегодняшний день кинетически-согласованных разностных схем, а также проведение сравнения результатов, полученных при моделировании по уравнениям Навье-Стокса и системе КГУ, с целью выявле-ния классов течения, где их решения мало отличаются, а где система КГУ пре-восходит по точности моделирования уравнения Навье-Стокса. Для ее дости-жения решаются следующие задачи:
представление КСРС, как схемы принадлежащей классу схем расщепления вектора потока для аппроксимации системы невязких уравнений Эйлера;
изучение свойств КСРС относительно семейства схем расщепления вектора потока, а также построения схем повышенного порядка точности на основе базовых КСРС первого порядка;
построение кинетически согласованной аппроксимации на неструктурированных сетках (треугольных и тетраэдральных);
построение аппроксимации системы КГУ на тетраэдральной сетке и проведение расчетов обтекания сферы при различных параметрах течения и сравнение с аналогичными расчетами, полученными при решении уравнений Навье-Стокса.
Достоверность результатов
Достоверность результатов, полученных в работе, подтверждается тем, что для выполненных в работе численных расчетов наблюдается качественное совпаде-ние полученных результатов с точным решением со-
ответствующих задач, экспериментальными данными и результатами других авторов. Эффективность предложенных кинетических алгоритмов подтверждается результатами тестирования разработанных на их основе программ. Основные положения диссертации, выносимые на защиту
Кинетически-согласованных схем как новый класс схем расщепления вектора потока.
Кинетически-согласованных схемы повышенного порядка точности.
Аппроксимация. кинетически-согласованными схемами уравнений Эйлера и Навье-Стокса на неструктурированных сетках.
Кинетически-согласованных схемы повышенного порядка точности на неструктурированных сетках.
Аппроксимация трехмерной квазигазодинамической системы уравнений на тетраэдральной сетке.
Публикации автора по теме диссертации
Абалакин И.В., Четверушкин Б.Н. Кинетически-согласованные разностные схемы как модель для описания газодинамических течений, Математическое моделирование, 1996. Т.6, № 8.
Абалакин И.В., Жохова А.В., Четверушкин Б.Н. Кинетически-согласованная аппроксимация газодинамических уравнений на треугольных сетках, Тезисы докладов VII Всероссийской школы-семинара "Современные проблемы математического моделирования", Ростов-на-Дону, 1997, с. 1-4.
Абалакин И.В., Жохова А.В., Четверушкин Б.Н. Кинетически-согласованные разностные схемы на нерегулярных сетках, Математическое моделирование, 1997. Т.9, № 7.
Абалакин И.В., Жохова А.В., Четверушкин Б.Н. Кинетически согласо-ванный алгоритм для расчета газодинамических течений на треугольных сетках, Математическое моделирование, 1998. Т. 10, № 4, с.51-60.
Абалакин И.В., Жохова А.В. Кинетически-согласованные разностные схемы с коррекцией на треугольных сетках, Дифференциальные уравнения, 1998, Т.34, № 7, с.904-910.
Абалакин И.В., Жохова А.В., Четверушкин Б.Н. Разностные схемы на основе кинетического расщепления вектора потока, Математическое модели-рование, 2000. Т.12, № 4, с.73-82.
Абалакин И.В., Жохова А.В., Четверушкин Б.Н. Кинетически-согласованные схемы повышенного порядка точности, Математическое моделирование, 2001. Т.13, № 5, с.53-61.
I.V. Abalakin, A. Dervieux, Т.К. Kozubskaya A vertex centered high order MUSCL scheme applying to linearised Euler acoustics, INRIA report, № 4459, 2002.
Абалакин И.В., Суков С.А., Моделирование внешнего обтекания тел на многопроцессорных системах с использованием тетраэдрических сеток. В сб. "Фундаментальные физико-математические проблемы
и моделирование технико-технологических систем", вып. 7, под ред. Л.А. Уваровой. М., Изд-во "Janus-K", 2004, с. 52-57.
Автор выражает благодарность своему научному руководителю Борису Николаевичу Четверушкину за формирование своих научных взглядов и постоянное внимание в ходе работы над диссертацией; Татьяне Геннадиевне Елизаровой, чьи советы и полезные замечания способствовали написанию диссертации; Людвигу Вацлавовичу Дородницыну и Сергею Александровичу Сукову за помощь в постановке граничных условий и параллельной реализации кинетических схем на многопроцессорной вычислительных системах.
Кинетически согласованные схемы для уравнения Эйлера
В основе кинетического подхода к построению кинетически согласованных разностных схем (КСРС) лежит уравнение Правая часть уравнения (1.16), называемая интегралом столкновений, описывает взаимодействие молекул друг с другом. Так в случае уравнения Больцмана J (/, / ) есть интеграл, зависящий от относительной скорости молекул и их прицельного расстояния, а также от функций распределения до / и после / столкновения молекул.
Для интеграла столкновений J (/, / ) выполняется законы сохранения массы, импульса и энергии при столкновениях, то есть моменты интеграла столкновений для сумматорных инвариантов (1.21) равны нулю [3] (что очевидно для БГК интеграла столкновений (1.17)).
Моменты (1.19) функции распределения не зависят от ее вида и справедливы для любой одночастичной функции рапределения. А моменты (1.22) определены неоднозначно и для определения явного вида тензора напряжений и вектора теплового потока необходимо знать вид функции распределения.
Величина v () есть частота столкновений для молекул, имеющих скорость , а величина г = l/i/(), обратная частоте столкновений, равна времени между столкновениями молекул. Столкновения молекул в газе приводят к установлению максвелловского вида функции распределения. Для реальных газов установление происходит сразу после столкновения (на этом основан вывод БГК уравнения [3]). Определим следующую модель изменения функции распределения, основанную на методе расщепления по физическим процессам уравнения Больцмана или БГК уравнения [8].
Как было отмечено выше, столкновения молекул приводит к максвел-лизации функции распределения за момент времени т. Действительно, рассмотрим релаксационное уравнение с интегралом столкновений вида (1.17) и начальной на момент времени tk функцией распределения / (х, , t).
Следуя работе [25], в этой модели расщепления сделаем следующее допущение. Так как время свободного разлета молекул и время релаксации величины одного порядка, то детально процесс релаксации в этой модели рассматриваться не будет. Будем предполагать, что процесс установления происходит мгновенно (мгновенная максвеллизация).
Принимая это предположение, исходная модель расщепления сводится к решению начальной задачи для уравнения переноса (1.24) на отрезке времени Atk, где в качестве начальной функции распределения в момент времени tk берется локально-максвелловская (1.17).
Введем разностную сетку с узлами Xj и постоянным шагом Ах = Xj+i — Xj и определим расчетную ячейку интервалом [а _і/2, 2 +1/2] с граничными точками равными я_/±і/2 = {xj±i + xj)/%- Будем предполагать что функция распределения постоянна на ячейке на любой дискретный момент времени tk, а, следовательно, постоянны и макропараметры pj, Uj, Tj, определяемые этой функцией. Обозначим эту дискретную функцию распределения через / (XJ, , tk) = ff.
С одной стороны, из соотношения (1.33) после умножения на сумматорные инварианты и интегрирования по , следует система уравнений Эйлера для определения моментов р, и, Т, входящих в максвелловскую функцию распределения. С другой стороны, если предположить, что функция в правой части соотношения (1.33) известная величина, то (1.33) можно интерпретировать как уравнение пространственно-временной эволюции локально-максвелловской функции распределения.
После интегрирования разностной схемы с сумматорными инвариантами и использования свойства обращения в нуль соответствующих моментов для интеграла столкновения, получим разностные уравнения вида (1.28)-(1.31).
При таком выводе эти разностные уравнения представляют собой разностную схему, аппроксимирующую систему уравнений Эйлера с первым порядком точности по времени и пространству, где потоки вида Ffcc(Q) описывают аппроксимационную вязкость полученной схемы. Обозначим эту схему через КС1.
Отметим, что аппроксимация уравнения (1.33) любой подходящей схемой для уравнения переноса, после осреднения с сумматорными инвариантами будет определять соответствующую кинетическую разностную схему. Но представляется более естественным использовать аппроксимации, согласованные с кинетической моделью (1.23)-(1.24), а именно аппроксимации типа направленных разностей.
Кинетические схемы повышенного порядка
Опыт использования кинетических схем показал хорошее совпадение как с экспериментальными данными, так и с расчетами по другим методикам. Но для более точного отображения реальной картины течения при численном моделировании, требуется использовать разностные схемы повышенного порядка точности как по времени, так и по пространству. Данная параграф посвящен построению кинетически-согласованных разностные схем обеспечивающий повышенный порядок точности по пространственным переменным. Отметим, что в работе [36] была получена кинетическая схема для уравнений Эйлера второго порядка точности на основе схемы повышенного порядка для уравнения Больцмана. В данной работе повышение порядка аппроксимации пространственных производных обеспечивается за счет вычисления кинетически расщепленных векторов потока от линейно-интерполированных газодинамических переменных. Такого типа аппроксимация уравнений Эйлера относится к классу так называемых MUSCL (Monotonic Upstream Schemes for Conservation Laws) схем [37, 38]. Будем рассматривать кинетически-согласованную схему как одну из разновидностей метода расщепления вектора потока, как это было сделано в предыдущем пункте и аналогично работе [30].
Отметим, что данную аппроксимацию можно получить непосредственно из уравнений Эйлера, не прибегая к кинетическому уравнению. Основная идея повышения порядка аппроксимации заключается в том, что согласно работе [38], газодинамические параметры на расчетной ячейке задаются не как кусочно-постоянные функция, а как кусочно-линейные.
Одним из важных свойств разностных схем для уравнения переноса является свойство монотонности или TVD условия для нелинейного гиперболического уравнения [31]. Это свойство позволяет избежать нефизических осцилляции разностного решения в окрестности разрывных течений (ударных волн).
Для того чтобы разностная схема повышенного порядка аппроксимации была монотонной (удовлетворяла TVD свойству для нелинейных уравнений) должно выполнятся условие минимума значения градиента при задании линейного распределения функции на расчетной ячейке [39, 40, 41, ].
В конце параграфа рассмотрим два примера численного расчета Тестирование кинетических аппроксимаций уравнений Эйлера проведем на одномерном задачи о распаде произвольного разрыва и двумерной задачи о сверхзвуковом обтекании ступени в прямоугольном канале. Выбор этих тестов обусловлен наличием точного решения (распад разрыва) и большим количеством результатов, полученных по другим методикам для задачи об обтекании ступени. Это позволяет судить об адекватности получаемых результатов расчета.
Для расчета бралась схема KClm повышенного порядка точности по пространственным переменным с ограничителями «minmod » и «superbee». Интегрирование по времени проводилось по явной схеме первого порядка. Все расчеты проводились с числом Куранта равным 0.4. На рисунке 1.9 приведены профили плотности на момент времени t = 4, рассчитанные по схемам с разным выбором параметров. Так введенные процедуры ограничения наклонов («minmod » и «superbee») в кусочно-линейном распределении позволили получить разностное решение без паразитических осцилляции в окрестности течений с большими градиентами и с резкими профилями ударных волн. Такое поведение разностного решения характерно для расчетов с использованием TVD схем в одномерном случае. Как видно из рисунка 1.9 наименее диссипативная схема — это схема с ограничителем «superbee», обладающим наибольшей сжимаемостью [44]. Контактный разрыв разрешается с минимальным числом точек. Поведение разностного решения в окрестности ударной волны у всех представленных вариантов схемы приблизительно одинаково. Ударная волна размазывается на 3-4 точки. Вариант с ограничителем «superbee» и ограничителем «minmod» с параметром (3 = 1/3 дают более резкий профиль в окрестностях сочленения волн и зон постоянного течения.
Результат расчета данной тестовой задачи показывает, что кинетическая схема позволяет проводить расчеты сильных разрывных течений с хорошим разрешением фронта сильной ударной волны. Волна разрежения имеет плавный профиль без изломов в окрестности звуковой точки (правый фронт волны разрежения), что позволяет судить об удовлетворении кинетической схемой энтропийного условия.
Кинетически-согласованные схемы повышенного порядка точности на нерегулярных сетках
Для получения кинетически-согласованной разностной схемы повышенного порядка на неструктурированной сетке перепишем схему (2.9) в виде схемы расщепления вектора потока Q?+1 = Of - "FT S« (F+(Qy) + F-(Qi)) (2.13) Здесь F± = (F ± D)/ 2, а значения Q определены на грани расчетной ячейки. В случае кусочно-постоянного распределения вектор-функции Q равен ( = Q; и QJ- = (. Далее будем проводить прямую аналогию с одномерным случаем. Предположим, что искомая функция является линейной на каждом треугольнике с вершиной в точке I. Значения газодинамических переменных на грани 13 между ячейками I и 3 определим аналогично формуле (1.91) где градиенты (VQ ) и (VQ ) определяются как среднее взвешенное между градиентом (VQ)7 = (Qj — Qi) и градиентом (VQ)j определенным в узле расчетной сетки как среднее по по контрольному объему с центром в данном узле. Здесь 1у — вектор, идущий из точки і в точку j
Рассмотрим более подробно вычисление узлового градиента (VQ)j. Проведем осреднение по всем соседним узлам и проинтегрируем по всей расчетной ячейке (Рис.2.1). Такой подход к вычислению градиента объясняется тем, что в случае нерегулярной сетки непонятно какой из соседних узлов считать внешним (или «предыдущим»).
Для построения монотонной схемы на нерегулярной сетке проведем прямую аналогию с построением схемы для одномерного уравнения переноса [55] и рассмотрим случай положительных скоростей переноса (то есть вычисления значения Q -). Для этого необходимо определить значение градиента в «полуцелой» / точке, лежащей на прямой проходящей через точки і" и J. Причем точка / равноудалена от точек I к IJ — точки пересечения границы контрольного объема с отрезком, соединяющим точки / и J. Нам известно значение узлового градиента в точке / — VQ; (или значение VQ ) и значение градиента в точке IJ, соответствующего градиенту в точке i-\-1/2 одномерного случая (Рис. 2.2). Используя значения этих двух градиентов экстраполяцией вниз по потоку определим значение градиента в точке I (VQ/0" = 2 (VQ ) - (( - Q0 (2.22)
В заключении этого параграфа рассмотрим модельную задачу об отражении косого скачка уплотнения от жесткой стенки (угол падения скачка равен 29, число Маха входного потока — 2.9). Размер расчетной области — О г/ 1, 0 х 4.1. Число узлов расчетной сетки выбиралось равным 1200. Данный рачсчет показал на сколько точнее и аккуратнее моделируются задачи с разрывными решениями по схемам повышенного порядка по сравнению со схемами первого порядка на неструктурированных сетках (Рис. 2.2)
Распределение давления вдоль оси x/L (L = 4.1) при у = 0.5. (с функцией-ограничителем типа «minmod») на треугольных сетках, при разных значениях коэффициента р\ При сопоставлении полученных результатов с точным решением можно видеть, что для (3 = 1 схема является наиболее диссипативной. Для (3 — 0.1 влияние искусственной вязкости незначительно, так как при построении линейного распределения на ячейке становится преобладающим влияние интерполяционного (центрально-разностного) градиента.
Под кинетическими схемами для описания вязкого газа будем понимать кинетически-согласованную разностную схему, аппроксимирующую систему уравнений Эйлера или конвективную часть членов системы уравнений Навье-Стокса. Диссипативную часть членов уравнений Навье-Стокса будем аппроксимировать методом конечных объемов с использованием контрольного объема изображенного на Рис. 2.2.
В связи с вышесказанным кинетическая система, аппроксимирующая уравнения Навье-Стокса запишется в виде (2.9)-(2.12) с добавлением вязких потоков Qf+1 = Qk \Sij\ (Fij(Q, ШІ)) - Dy(Q, ny)) + H Q, Щ,-))), (2.25) где вязкие составляющие потоков Hjj = (О, Щ Щр Hf j выводятся следующим образом. Тензор вязких напряжений в двумерном случае записывается в локальной системе координат (п,т). Далее, применяя теорему Гаусса к объемному интегралу от дивергенции тензора, определим вязкие потоки, стоящие в разностных уравнениях движения. Для уравнения энергии проводим аналогичные действия, только вместо тензора вязких напряжений рассматриваем скалярную величину, представляющую собой сумму диссипативной функции и потока температуры.
В конце этого параграфа рассмотрим задачу об обтекании крылового профиля NACA0012 потоком вязкого теплопроводного газа с числом Маха М = 0.8 и числом Рейнольдса Re = 10000. Известно [56], что при таком числе Маха реализуется трансзвуковое течение. Движение газа вдоль профиля ускоряется и образуется местная сверхзвуковая зона, которая замыкается скачком уплотнения, за которым скорость вновь становится дозвуковой. Причем положение скачка уплотнения зависит от толщины профиля при фиксированном числе Маха внешнего потока.
Расчет этого течения проводился на адаптивной к решению сетке. Адаптация сетки к решению основана на построении метрики, связанной с матрицей Гесса от выбираемого газодинамического параметра. Такая методика была предложена в работе [57] и реализована в программном коде Ani2D любезно предоставленным автору Ю.В. Василевским.
Для течений, в которых присутствуют ударные волны и существенны вязкие эффекты, такие как пограничный слой, эффективным анализатором происходящих процессов является число Маха, которое зависит как от скорости (эффект пограничного слоя), так и от локальной скорости звука (эффекты ударных волн). Поэтому в качестве переменной, по которой строился гессиан разностного решения, выбиралось число Маха.
Процедура адаптации проводилась через 1000 временных шагов по времени до выхода разностного решения на устойчивый квазистационарный режим течения. Так на рисунках 2.7 и 2.9 показана финальная сетка после 4 этапов адаптации. Полученная сетка состоит из 41993 треугольников и имеет 21169 узел. Минимальный по площади треугольник финальной сетки с высотой равной h = 8.73 X 10 4 и находится в окрестности задней кромки профиля, где градиенты скорости достигают максимальных значений. На рисунке 2.7 можно видеть измельчение сетки в следе за профилем, где формируется дорожка Кармана, а также анизотропная деформация треугольников в направлении перпендикулярным скачку уплотнения.
На рисунке 2.10 показано распределение локального числа Маха на начальной сетке на момент формирования стационарного скачка уплотнения, на рисунке 2.11 изображено тоже поле локального числа Маха, посчитанного на сетке после финальной адаптации. Из сравнения этих рисунков видно, что адаптация сетки позволяет более аккуратно разрешать ударные волны и более точно рассчитывать течения в следе за обтекаемым телом.
Из результатов проведенного расчета можно сделать вывод, что аппроксимация системы уравнений Навье-Стокса кинетическими схемами на неструктурированных сетках, в том числе и на адаптивных, позволяет проводить расчеты сложных течений с отрывами пограничного слоя и областями взаимодействия ударных волн и пограничных слоев.
Квазигазодинамическая система
Вид вектора консервативных переменных Q, конвективных векторов потока Fj и диссипативных потоков Wy- приведен в Приложениях А и В. Нетрудно убедится, что диссипативные потоки W,j есть однородные функции степени один относительно вектора консервативных переменных (аналогично рассмотренному в Главе 1 свойству конвективных потоков и диссипативных потоков F; и F 2).
Исходя из вида системы (3.7) матрицы Якоби диссипативных потоков помноженные на величину г можно назвать матрицами вязкости. Несложно проверить (например, используя систему символьных вычислений MAPLE).
Сравнивания предыдущий вектор с диссипативными членами в одномерной системе уравнений Навье-Стокса получим вид диссипативных коэффициентов, таких как молекулярная ц и объемная вязкости и молекулярный коэффициент теплопроводности А. коэффициентов вязкости и теплопроводности получается при выводе уравнений Навье-Стокса методом Чепмена-Энскога из уравнения БГК [3]. При этом число Прандтля Pr = всегда бу„ое„.
Таким образом, определение диссипативных коэффициентов вида (3.11)-(3.13) и структура вектора типа (3.10) показывают, что матрицы гС действительно определяют те же самые диссипативные члены, что и в уравнениях Навье-Стокса при числе Pr = 1. Следовательно, квазигазодинамическая система уравнений включает в себя систему уравнений Навье-Стокса и отличается от нее дополнительными диссипативными членами.
В книге [16] показано, что дополнительные диссипативные члены (3.16) представляют собой самодиффузионный и термодиффузионный потоки в уравнении неразрывности, нелинейную добавку к тензору вязких напряжений в уравнениях импульса и дополнительного теплового потока в уравнении энергии. Как можно видеть из представления (3.16) и вида матриц приведенных в Приложении С, добавка к тензору вязких напряжений входит с коэффициентом вязкости зависящей, в частности, от квадратов скоростей.
Из представления системы КГУ в виде (3.14)-(3.16) дополнительные слагаемые имеют порядок малости 0(т) при соответствующей гладкости производных. Покажем, что в случае стационарных течений дополнительные добавки к системе имеют не первый, а второй порядок малости — 0(т2). Вывод этого результата будет отличаться от соответствующих рассмотрений в работах [16],[4], [5].
Отсюда следует, что стационарная квазигазодинамическая система отличается от стационарной системы уравнений Навье-Стокса диссипатив-ными добавками порядка 0(т2). Покажем, что дополнительная диссипация имеет порядок 0(т2) и для нестационарной системы квазигазодинамических уравнений. Для того разностную производную по времени с точностью до членов второго порядка малости по г можно представить в виде следующего разложения
Численная реализация системы КГУ проводится на неструктурированной тетраэдральной сетке. Для аппроксимации будем использовать систему КГУ, записанную в форме (3.7). Для аппроксимации конвективного потока использовалась кинетически-согласованная схема типа КС2 повышенного порядка точности на барицентрическом контрольном объеме. Аппроксимация вязких членов проводится методом конечных элементов, как и в случае аппроксимации вязких членов в уравнениях Навье-Стокса. Рассмотрим более подробно конечно-элементную аппроксимацию диссипа-тивных членов в квазигазодинамической системе уравнений.
Далее, дважды применяя к интегралу (3.25) теорему Грина и учитывая, что потоки Wnm определены в узлах тетраэдра, получаем следующую аппроксимацию квазигазодинамических потоков методом конечных элементов
В формуле (3.26) внутренне суммирование (индекс j) идет по вершинам тетраэдра, отличным от вершины щ, значения градиентов финитной функции определяется формулой (3.24), обозначение [/л/р] есть среднее арифметическое по вершинам тетраэдра.
Для численной проверки асимптотики системы КГУ, рассмотренной в параграфе 3.3 выбрана стационарная задача внешнего обтекания сферы с числом Рейнольдса Re = 25 и числом Маха равным М = 0.1. Малым параметром, определяющим асимптотическое поведение системы КГУ в безразмерном виде, можно считать число Кнудсена Кп = 0.1 равное отношению числа Маха к числу Рейнольдса Кп = M/Re [3]. В условиях этой задачи число Кп = 0.004, что соответствует плотному газу. Целями этого расчета являлись: тестирования численной аппроксимации системы КГУ на тетраэдральных сетках и проверка близости решений уравнений Навье-Стокса и системы КГУ, которая имеет место при малых числах Кнудсена.
Для корректности сравнения численных решений система уравнений Навье-Стокса аппроксимировалась по той же схеме, что система КГУ (конвективные потоки — кинетически-согласованная схема повышенного порядка, диссипативные члены — методом конечных элементов).
В заключении кратко суммированы основные результаты работы, выносимые на защиту:
1. Предложены новые методы кинетического расщепления. Проведен анализ кинетически-согласованных разностных как схем расщепления вектора потока. Выведены новые аппроксимации уравнений Эйлера на основе кинетиче-ского расщепления вектора потока.
2. Разработаны алгоритмы повышения порядка точности кинетически согласованных разностных схем с использованием методики MUSCL.
3. Построены кинетических схемы, аппроксимирующие уравнения Эйлера с повышенным порядком аппроксимации на неструктурированных треугольных и тетраэдральных сетках.
4. Проведен сравнительный анализ квазигазодинамической системы уравнений и системы уравнений Навье-Стокса как системы уравнений. Построена аппроксимация квазигазодинамической системы на неструктурированных тетраэдральных сетках.