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

кандидата физико-математических наук
Воронич, Иван Викторович
город
Москва
год
2005
специальность ВАК РФ
05.13.18
Диссертация по информатике, вычислительной технике и управлению на тему «Методы решения нестационарных задач газовой динамики»

Автореферат диссертации по теме "Методы решения нестационарных задач газовой динамики"

На правах рукописи

I оБЯЗАдаг-:-' 1

бесплатный _]

Воронич Иван Викторович

УДК 519.6:533.6

МЕТОДЫ РЕШЕНИЯ НЕСТАЦИОНАРНЫХ ЗАДАЧ ГАЗОВОЙ ДИНАМИКИ

(05.13.18 - математическое моделирование, численные методы и комплексы программ)

АВТОРЕФЕРАТ диссертации на соискание ученой степени кандидата физико-математических наук

Москва, 2005

Работа выполнена в Московском физико-техническом институте (государственном университете)

Научный руководитель: доктор физико-математических наук

профессор Хлопков Юрий Иванович

Официальные оппоненты: доктор физико-математических наук

профессор Толстых Андрей Игоревич

на заседании диссертационного совета К 212.156.07 при Московском физико-техническом институте (государственном университете) по адресу: 140180, г. Жуковский М.О., ул. Гагарина, д. 16.

С диссертацией можно ознакомиться в библиотеке ФАЛТ МФТИ.

кандидат физико-математических наук доцент Жаров Владимир Алексеевич

Ведущая организация: ГНЦ Центральный институт авиационного

моторостроения им П.И. Баранова

Защита диссертации состоится «_»

2005 года в_часов

Автореферат разослан «_»

2005 года.

Ученый секретарь диссертационного совета кандидат физико-математических наук, доцент

Киркинский А.И.

¿Ш-Ч

15914

ОБЩАЯ ХАРАКТЕРИСТИКА РАБОТЫ

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

Целью работы является построение вычислительно эффективного, надежного алгоритма решения нестационарных задач газовой динамики в широком диапазоне режимов, пригодного для раскрытия волнового поля течения, содержащего разномасштабные особенности, на основе систематизации свойств методов типа Годунова в сравнении с методами расщепления вектора потока и с привлечением адаптивных функций-ограничителей. К цели работы относится также применение развитых методов к решению содержательных задач, в том числе к задаче о взаимодействии вихря с профилем в дозвуковых и трансзвуковых потоках сжимаемого газа и генерации возмущений в ходе такого процесса.

Научная новизна выносимых на защиту результатов заключается в следующем:

1. получены новые данные о механизмах и характеристиках схемной диссипации методов типа Годунова и методов расщепления вектора потока;

2. создан новый комбинированный алгоритм вычисления аппроксимации вектора потока, позволяющий расширила' дилмзои прнмщщиия методов типа Годунова; *ос "ацн^нальнау

библиотека

С.Летер{ 09

3. создана новая адаптивная функция-ограничитель в рамках схем, сохраняющих монотонность, позволяющая повысить точность и вычислительную эффективность методов без повышения требований к расчетной сетке;

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

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

На защиту выносятся следующие положения:

1. результаты сравнительного анализа основных групп численных методов решения уравнений Эйлера, позволяющие определить область применения методов и уточнить механизмы схемной диссипации;

2. комбинированный алгоритм вычисления аппроксимации вектора потока, сочетающий применение двух аппроксимаций вектора потока, что позволяет расширить диапазон применения методов типа Годунова;

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

4. результаты численного решения задачи о взаимодействии вихря с профилем в дозвуковых и трансзвуковых потоках сжимаемого газа.

Апробация работы. Результаты работы докладывались и обсуждались на ряде научных конференций и семинаров: XLI-XLVI научные конференции МФТИ

XIV, XV Международные Школы по моделям механики сплошной среды

Конференции RRDPAE 1996 г., 1998 г. (Warsaw, Poland)

21st Int. Symposium on Rarefied Gas Dynamics - 1998 r. (Marseille, France)

Семинар кафедры физики ФАЛТ МФТИ - 2000 г.

Семинар отдела прикладной математической физики ВЦ РАН - 2000 г.

Школа-семинар «Актуальные проблемы аэрокосмической науки» - 2001 г.

Семинар кафедры прикладной математики МИФИ - 2004 г.

Семинар НИО 8 ЦАГИ - 2004 г.

Семинар ЦИАМ - 2005 г.

Публикации. Основные результаты диссертации опубликованы в 7 работах.

Структура и объем работы. Диссертация состоит из введения, четырех глав, заключения, списка цитированной литературы из 120 наименований. Материал изложен на 150 страницах, приведено 55 рисунков.

СОДЕРЖАНИЕ РАБОТЫ

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

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

В параграфе 1.1 рассматриваются различные формы уравнений Эйлера и их характеристические свойства. Рассмотрены свойства векторов потока и соответствующих матриц Якоби, включая свойства собственных векторов и собственных значений, необходимые для дальнейшего анализа. Приведены дифференциальные квазилинейные формы уравнений Эйлера, служащие базисом для получения характеристических форм, применяемых для построения численных методов и постановки граничных условий.

В параграфе 1.2 описана формулировка полудискретной задачи, получаемой при аппроксимации вектора потока через грань расчетной ячейки в рамках метода конечного объема.

В параграфе 1.3 рассмотрены задача Коши для гиперболической системы законов сохранения, понятие обобщенного решения, условие неубывания энтропии применительно к системе законов сохранения и к разностной схеме, ее аппроксимирующей. В связи с этим затрагиваются понятия аппроксимации, устойчивости и сходимости и рассматривается понятие схемы, согласованной с условием неубывания энтропии.

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

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

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

Решение описанных задач позволяет расширить диапазон применимости существующих методов и повысить их вычислительную эффективность.

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

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

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

Рассматриваемые варианты методов расщепления вектора потока {параграф 2.2) основаны на представлении вектора потока в виде суммы двух компонент, отвечающих характеристическим скоростям определенного знака. С одной стороны, это обеспечивает простоту таких методов, с другой приводит к повышенной схемной диссипации, в том числе при расчете контактных и тангенциальных разрывов. Конструкция таких методов не учитывает всех характеристических свойств законов

сохранения, вопрос об энтропийной согласованности также является открытым. Исключение составляет метод кинетического расщепления вектора потока, который построен на основе вычисления компонент вектора потока с помощью интегралов от плотности распределения. Для этого метода вопрос энтропийной согласованности решен положительно.

Анализ свойств собственных значений матрицы схемной диссипации методов дает много содержательной информации, ранее не приводившейся в полном объеме (рис. 1). Во-первых, указанные собственные значения для методов Стегера-Уорминга и ван Лира явно показывают существование проблемы звуковой точки, в которой одно из собственных значений обращается в нуль. Во-вторых, такие данные позволяют судить об относительном уровне схемной диссипации методов. В-третьих, наличие изломов и разрывов ветвей собственных значений говорит о некоторой нерегулярности механизма схемной диссипации методов Стегера-Уорминга и ван Лира. От них выгодно отличается метод кинетического расщепления вектора потока, который обладает регулярным механизмом схемной диссипации и энтропийной согласованностью. Традиционное ограничение на шаг по времени для рассматриваемых методов первого порядка точности в одномерном случае заключается в выполнении неравенства СП < 1, где СРЬ - число Куранта-Фридрихса-Леви. Более тонкие вопросы необходимости и достаточности этого ограничения для конкретного метода требуют дополнительного исследования.

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

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

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

Необходимое количество частиц в ячейке для моделирования определяется исходя из того, чтобы статистическая ошибка удовлетворяла требованию согласования ошибок: 5 ~ О(Ддг). С точки зрения вычислительных затрат точность результатов, получаемых с помощью методов ПСМ, определяется в первую очередь статистической погрешностью.

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

вычислительные затраты, а также демонстрирует взаимосвязь кинетических и сплошносредных моделей.

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

В параграфе 2.4 рассматриваются энтропийные свойства конечно-разностных методов, в первую очередь метода кинетического расщепления вектора потока. Для него можно получить явные выражения энтропийной функции и аппроксимации потока энтропии, а также доказать дискретное энтропийное неравенство. В кинетической теории рассмотрение ведется с помощью Я-функции, которая характеризует вероятность состояния системы, это понятие применимо к равновесным и неравновесным течениям. Я-теорема Больцмана является обобщением условия энтропии на неравновесные системы из большого числа частиц.

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

Параграф 2.5 посвящен развитию методов повышенной точности. В качестве базовой модели, используемой для построения метода повышенного порядка аппроксимации, традиционно рассматривается линейное уравнение переноса. Задача заключается в повышении порядка аппроксимации метода при соблюдении дополнительных ограничений, таких как условие сохранения монотонности. Эти ограничения призваны до некоторой степени заменить энтропийное условие. При этом разностная схема повышенной точности должна быть нелинейной в силу теоремы Годунова. Нелинейность разностной схемы связана с применением функций-ограничителей, используемых при реконструкции распределений переменных внутри ячейки. Подход, развитый для скалярного случая, обобщается на полную систему уравнений. Удовлетворить строго условию неубывания энтропии в этом случае не представляется возможным, тем не менее, можно надеяться, что сеточное решение будет высокого качества во

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

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

Аппарат функций-ограничителей допускает различные варианты таких функций. Под функцией-ограничителем понимается функция, зависящая от отношения приращений переменной на трехточечном шаблоне. Пространственный шаблон метода второго порядка точности является пятиточечным. Помимо набора обязательных требований к функции-ограничителю, таких как непрерывность, симметричность, равенство нулю при отрицательном значении аргумента, согласованность, ограниченный предел (<2) при стремлении аргумента к бесконечности, может быть наложен ряд дополнительных требований, регулирующих гладкость функции-ограничителя и разрешающие свойства метода. На практике необходим компромисс между «сжимающими» свойствами ограничителя (размывание разрывов) и его надежностью. Для расчета разрывных решений, как правило, используются ограничители, близкие к базовой функции пнптос!, предложенной В.П. Колганом:

где С/Д"1*2 - величины на грани расчетной ячейки после реконструкции, у, -переменная в ячейке /.

Практика показывает, что функции-ограничители с гладким переключением при Х- 1, лежащие близко к функции пиптос!:

(1)

= V, + , гу,

1+1/2 _ У/+1

тт(Ь,ЬХ, (1 + Х)/2), Х>0

Р.

ХйО

обладают повышенными разрешающими свойствами. Для расчета широкого класса течений можно использовать значения параметра «сжатия» Ъ в диапазоне [1,05, 1,3]. При b -> 1+0 эта функция переходит в функцию minmod. На рис. 2 изображены некоторые известные функции-ограничители F(X).

Свойства многомерных разностных схем, также как и применение такого рода способа повышения точности метода по пространственным переменным приводят к дополнительному ограничению на шаг по времени (в двумерном случае): CFL < ХЛ.

Аппарат адаптивных к локальной картине течения функций-ограничителей продолжает развиваться с целью более полного использования потенциала такого подхода. В работе М. Агога и P. Roe (1997) предложено учитывать локальное число Куранта в качестве дополнительного аргумента функции-ограничителя для компенсации потери точности при малых значениях локального числа Куранта. В работе А.Н. Гильманова (2000) в рамках подхода ограничения характеристических переменных предложено идентифицировать тип особенности (скачок, контактный разрыв, волна разрежения), возникающей на грани расчетной ячейки в каждом характеристическом поле, и в зависимости от этого использовать более или менее «сжимающие» ограничители. Идея, изложенная в работе G. Billet, О. Louedin (2001), согласуется с концепцией, согласно которой не всегда необходима «монотонизация» при вычислении реконструированных значений переменных (Ю.Б. Радвогин (1991), В.В.Остапенко (2000)). Вместо обычного способа реконструкции переменных предлагается применять схему центральных разностей (F= 1) в рассматриваемой ячейке, если в соседней ячейке в направлении сеточной линии имеется экстремум по ограничиваемой переменной {Х,,А < 0). Авторы работы аргументируют это тем, что таким образом появляются условия для развития каскада энергии, т. е. передачи энергии между модами с различными волновыми числами. Если изменение переменной на 5 узлах с центром в рассматриваемой ячейке является монотонным, применяется «сжимающая» функция-ограничитель. В случае X, < 0, как и ранее, F = 0. Таким образом, общий пространственный шаблон увеличивается до 7 узлов, но в отличие от методов 3-го порядка точности требования к качеству расчетной сетки не возрастают, так как дополнительные значения переменной нужны для качественного анализа. Опыт применения такого подхода показал, что необходимо учитывать и соотношение приращений значений переменной в соседней ячейке в случае наличия в ней экстремума. При \ХА\ \ = 1 (волновая картина) такой подход оправдан, при больших значениях тах(|Л",±||, 1/|ЛС±,|) (фактически это разрыв) нужна «монотонизация» сеточного решения. Построенная функция-ограничитель может быть выражена в виде:

^±|<0&гпах(|Х,±||,1/^,±||)<1,5

в оставшемся случае

1-1

Здесь Р\(Х) - функция-ограничитель, близкая к гшптос!, Р2(Х) -«сжимающая» функция-ограничитель. В качестве Р\(Х) может быть взята функция (2) с Ь=1,3, в качестве Р2(Х) - функция (2) с 6=1,3-^2 в зависимости от свойств рассматриваемого течения.

При интегрировании полудискретной задачи по времени в рамках явных схем второго порядка использованы двухшаговые методы Рунге-Кутга второго порядка точности:

где О! - вектор консервативных переменных на слое по времени г, Ьд -конечно-разностный оператор, р\, р2 - коэффициенты метода. В работе использованы два метода Рунге-Кутга, один из которых обладает свойством сохранения монотонности (р\ = \, р2 - 1/2), другой является низкодиссипативным методом (р1 = 1/2, р2 = 1). Практика показывает, что во многих случаях различия результатов, обусловленные применением этих методов, несущественны для коротких и средних времен моделирования. Нужно отметить также желательность согласования порядка точности метода интегрирования по времени с порядком аппроксимации пространственных производных.

Параграф 2.6 посвящен обсуждению вариантов построения комбинированного метода, сочетающего применение двух аппроксимаций вектора потока.

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

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

й'^ + дА'ьДш),

С?Г =а + Д/((1-А)ЬД({#}) +АЬД({Й'})), р,р2 = У2■ (4)

такими методами. Однако существует ряд ситуаций, которые не подпадают под такое объяснение. Одна из этих ситуаций рассматривается на примере задачи о течении с Маховским типом отражения ударной волны в главе 3. Кроме этого, общими для рассматриваемых методов являются проблема «перегрева» в области столкновения ударных волн и проблема расчета областей с низкой плотностью при разлете газа. Качественно первое явление можно наблюдать в вязком теплопроводном газе, второе явление попадает по существу в область кинетических моделей. Поэтому для построения адекватных численных методов в этом случае представляется необходимым применять модели, основанные на явном кинетическом диссипативном механизме.

Актуальна разработка комбинированного алгоритма, сочетающего разрешающие свойства для контактных и тангенциальных разрывов и надежность для интенсивных многомерных ударных волн. Поиск решения этой проблемы велся преимущественно в направлении понижения схемной диссипации методов расщепления вектора потока. J-M. Moschetta, D. Pullin (1997) для понижения схемной диссипации метода кинетического расщепления вектора потока применили метод Ошера приближенного решения задачи Римана на основе центрированных волн. Однако такой способ приводит к существенному росту вычислительных затрат. В работе Y. Wada, M.-S. Liou (1997) за основу принят метод расщепления AUSM (Advective Upwind Splitting Method), диссипация которого на контактных и тангенциальных разрывах понижается специальным выбором выражения для аппроксимации потока массы. Интересно, что в этом случае проявляется «карбункулярный» дефект. Методы этого класса основаны на расщеплении вектора потока на конвективную и звуковую части, что приводит к существенным дефектам сеточного решения на интенсивных скачках даже в одномерном случае. Поэтому в области интенсивных скачков было предложено использовать более диссипативный метод.

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

Ф*™. К + Г„,|/2 < 0,75 max (| ,\VnR илм

Ф=" (pL + pR)/2<0,15ma.x{pL,pR) (5)

Фв противном случае

Здесь ФКРУ5 - аппроксимация по методу кинетического расщепления вектора потока, - аппроксимация по методу расщепления разности вектора потока, У„ - нормальная к грани ячейки компонента скорости, р -давление. Таким образом, основным является метод типа Годунова, а при достаточно сильном перепаде нормальной скорости или давления применяется метод расщепления вектора потока. Такой способ основан на идентификации интенсивных разрывов, в том числе в начальный момент времени. Перепады касательной компоненты скорости и плотности не рассматриваются как условие переключения, так как контактные и тангенциальные разрывы не нуждаются в большом количестве схемной диссипации. Описанный подход нуждается в дальнейшем развитии.

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

Одномерные тестовые задачи (параграф 3.1) представляют варианты задачи Римана с различными начальными данными и включают решения со слабыми и интенсивными разрывами и волнами, с переходом через звуковую точку при различных условиях, с разлетом газа, со столкновением потоков, с медленно движущимся контактным разрывом. Подобные задачи широко используются для оценки свойств численных методов. Результаты решения одномерных задач позволяют проверить выводы, сделанные ранее о свойствах схемной диссипации методов, а также проверить пригодность комбинированного алгоритма и адаптивной функции-ограничителя на различных газодинамических конфигурациях. Список рассматриваемых задач следующий: акустическая задача, задача о столкновении потоков, задача о разлете газа, задача Лакса, задача Мах 3, задача об интенсивном разрыве, задача Сода, задача о контактном разрыве, задача Шу-Ошера. Расчетная область (0 < х < 1) разбивается на 80 ячеек (в последней задаче 0 < х й 10, 400 ячеек). Граничные условия на краях расчетной области ставятся в виде экстраполяции значений примитивных переменных в «мнимые» ячейки, применяемые в методе конечного объема. Момент времени, до которого проводится вычисление, выбран так, чтобы возмущения не доходили до границ расчетной области. Показатель адиабаты газа у = 1,4. Шаг по времени при расчете всех задач

определяется на основе числа Куранта, СРЬ = 0,45. Для каждой задачи рассматриваются распределения значений плотности и числа Маха в конечный момент времени. На рис. 3-5 приведены распределения для задачи Лакса, задачи об интенсивном разрыве и задачи Шу-Ошера.

Результаты расчета тестовых одномерных задач позволяют сделать следующие выводы: 1) в случае скачков и контактных разрывов разница между получаемыми по разным методам решениями в большинстве случаев оказывается незначительной; 2) метод типа Годунова и комбинированный метод дают существенно лучшее качество сеточного решения для слабого контактного разрыва и для медленно движущегося контактного разрыва, заданного в начальных условиях; 3) в случае волн разрежения все методы сравнимы по качеству сеточного решения, за исключением метода ван Лира, который дает дефект вблизи звуковой точки в интенсивных волнах разрежения, что подтверждает выводы главы 2; 4) применение комбинированного метода показало, что в рассмотренных задачах он наследует разрешающие свойства от метода типа Годунова и свойства надежности от метода кинетического расщепления вектора потока; 5) использование метода повышения точности на базе адаптивной функции-ограничителя позволяет повысить точность решений практически без увеличения вычислительных затрат.

В параграфе 3.2 рассматриваются двумерные задачи, имеющие как методическое, так и научное значение: 1) задача о переносе вихря, 2) задача о потенциальном обтекании цилиндра, 3) задача о стационарном течении со скачком Маха. Для расчетов использовались алгоритмы с применением функции-ограничителя вида (2) с Ь = 1,07 и адаптивного ограничителя (3) с 6=1,3 для /м(X) и 6= 1,35-И ,9 для Рг(Х) и метод интегрирования по времени (4) с рх = 1/2, р2 = 1. Шаг по времени при расчете всех задач определяется на основе числа Куранта, CFL = 0,45.

Первая из рассматриваемых задач представляет перенос изоэнтропического вихря, являющегося точным гладким (кроме центра вихря) решением уравнений Эйлера, что позволяет анализировать диссипативные и дисперсионные характеристики методов. Задача является нестационарной, так как рассматривается не в системе центра вихря, а в неподвижной системе координат. Расчетная область представляет собой единичный квадрат, на границах области используются периодические граничные условия. Вихрь задается в начальных данных, перенос осуществляется при скорости сносящего потока, соответствующей числу Маха М„ = 0,5, показатель адиабаты у = 1,4. Время эволюции выбрано из условия прохождения сносящим потоком пяти размеров расчетной области: Т= 5Ь/УХ. Для расчетов применялись как равномерные, так и хаотически возмущенные расчетные сетки. Относительная амплитуда сеточных возмущений во втором случае составляет 5%. Это позволяет оценить снижение точности решения для неравномерных сеток.

Размерность используемых расчетных сеток составляет 50x50, 75x75, 100x100, 200x200 ячеек. В качестве меры отклонения от точного решения используется ошибка в норме L2.

Изолинии полей плотности, полученных с помощью методов кинетического расщепления вектора потока (KFVS) и расщепления разности векора потока (FDS) на равномерной расчетной сетке 100x100, представлены на рис. 6. На сетках с недостаточным разрешением (50x50, 75x75) порядки сходимости численных решений близки к первому, что можно было ожидать исходя из способа повышения точности (параграф 2.5). Порядки сходимости решений на сетках с хорошим разрешением (100x100, 200x200) близки ко второму, в том числе на возмущенных сетках. Нерегулярность расчетной сетки увеличивает ошибку решения, рост ошибки зависит от свойств метода и более ощутим на подробных расчетных сетках. Метод KFVS дает во всех случаях большую ошибку решения, что подтверждает сделанный ранее вывод о диссипативных свойствах методов расщепления вектора потока. Применение адаптивной функции-ограничителя позволяет существенно повысить точность численного решения (снижение ошибки на 25-40%), что расширяет область применения методов второго порядка аппроксимации.

Вторая задача представляет плоское стационарное невязкое течение около кругового цилиндра при малой скорости набегающего потока, соответствующей практически несжимаемому режиму течения (р = const). Такое течение обладает свойствами симметрии и содержит детали, которые с трудом разрешаются при численном моделировании. Точность решения в этом случае зависит от свойств метода, качества расчетной сетки и реализации граничных условий. Расчет проводится при малом, но конечном значении числа Маха (М„ = 0,1, у= 1,4). Эффекты сжимаемости в этом случае незначительны, поэтому можно сопоставить решение с известным теоретическим решением для несжимаемого потенциального течения.

Реализация граничного условия непротекания на неподвижных твердых стенках в рамках метода конечного объема предполагает вычисление газодинамических параметров в «мнимых» ячейках, являющихся зеркальным отражением приграничных ячеек относительно их сторон, принадлежащих границе. При этом нужно учитывать кривизну твердых стенок, без чего невозможно аккуратно смоделировать течение вблизи стенки. На основе анализа вариантов реализации граничного условия непротекания в диссертации принят следующий подход (A. Dadone, В. Grossman (1994)):

Здесь и - нормальная к поверхности координата, Л - радиус кривизны стенки, я - энтропия, р - давление, - нормальная и касательная компоненты скорости, И - полная удельная энтальпия. Первое соотношение (6) предполагает использование уравнения импульса для определения перепада давления между приграничной и «мнимой» ячейками (с учетом знака кривизны стенки). Остальные соотношения представляют условия равенства нулю нормальных градиентов энтропии и полной удельной энтальпии на стенке, а также обращение в нуль нормальной компоненты скорости на стенке. Значения плотности и касательной компоненты скорости берутся из приграничной ячейки. Исходя из дискретизации этих соотношений на расчетной сетке, можно найти значения всех газодинамических переменных в «мнимой» ячейке. Для методов повышенной точности вводятся два и более слоев «мнимых» ячеек на границах, значения параметров в которых находятся аналогично. При наличии больших перепадов параметров вблизи стенки нужно следить за тем, чтобы давление в «мнимых» ячейках было положительным.

В качестве граничных условий на внешних границах расчетной области для задач внешнего обтекания целесообразно использовать те или иные волновые условия, обеспечивающие выход без отражения нестационарных возмущений из расчетной области при соблюдении условий набегающего потока. Одним из простых способов реализации таких условий на внешних дозвуковых границах является способ, основанный на использовании системы характеристических уравнений для касательной компоненты скорости, энтропии и инвариантов Римана, вдоль соответствующей сеточной линии, нормальной к границе. Зачастую нормальность сеточных линий к внешней границе не выдерживается строго, тем не менее, такой подход позволяет получить более корректное решение по сравнению со случаем «жестких» граничных условий.

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

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

Для задачи использовались сетки размерности 80x60, 120x90, 160x120, полученные с помощью функции растяжения Эйземана при Л = 0,5 (радиус цилиндра), Лтах = 10 (радиус расчетной области).

Изолинии полей давления, полученных с помощью методов кинетического расщепления вектора потока (КРУ8) и расщепления разности вектора потока (РОБ) на расчетной сетке 160x120, в сравнении с теоретическим решением приведены на рис. 7. Важную роль при решении задачии сыграл подход к постановке и реализации граничных условий на внешних границах и твердых стенках, что позволило с хорошей точностью раскрыть детали течения в согласии с теоретическим решением. Вывод, сделанный о повышенной диссипации, присущей методу кинетического расщепления вектора потока, подтверждается и в данном случае, уже применительно к стационарной задаче. Отметим повышение точности решения при использовании адаптивного ограничителя (снижение величины Сх, теоретически равной нулю, на 40-50%).

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

Течение возникает в результате втекания сверхзвукового потока (М»= 1,8, у= 1,4) в плоский канал (угол клина 0 = 15°) и является разрывным решением системы законов сохранения (рис. 8). Такую конфигурацию лучше рассматривать как половину конфигурации симметричного воздухозаборника (чтобы свести к минимуму влияние физической вязкости). При определенных условиях возникает конфигурация с тройной точкой (Т), состоящая из падающего скачка (18), отраженного скачка (ЯБ), скачка Маха (М8) и тангенциального разрыва (ББ). Кроме этого, веер разрежения (ЕР) взаимодействует с отраженным скачком (Яв) и тангенциальным разрывом (88), что вносит дополнительные особенности. Взаимное положение веера разрежения (ЕР) и падающего скачка (1Б) в рамках геометрических пропорций канала регулирует положение и высоту скачка Маха (М8). Под действием веера разрежения (ЕР) искривляется отраженный скачок (118) с формированием вниз по течению непрерывного энтропийного следа (на рисунке не показан). Тангенциальный разрыв (88) деформируется под действием прошедшего через отраженный скачок (ИЗ) веера разрежения (ЕР), в результате чего дозвуковое течение за скачком Маха (МБ) (который может быть слегка изогнут) переходит в сверхзвуковое течение на звуковой

линии (БЬ) в горле канала, образованного плоскостью симметрии и тангенциальным разрывом (ББ) аналогично течению в сопле Лаваля.

Формирование скачка Маха при различных геометрических конфигурациях наряду с вопросом перехода между регулярным и Маховским типами отражения (включая явление гистерезиса) является классической и в тоже время современной практически значимой проблемой. Течения с Маховским типом отражения скачков встречаются как в стационарных, так и в нестационарных ситуациях при различных условиях. Как правило, рассматривается вариант Маховского отражения при достаточно больших скоростях потока на входе (М* > 3-И). Это обусловлено запросами практики и изучением явления гистерезиса, которое возможно благодаря наличию области неединственности газодинамического решения при М» > 2,2 (2,66) при у = 1,4. Теоретический анализ в этом случае облегчается тем, что за отраженным скачком уплотнения (ЯБ) поток всюду сверхзвуковой, дозвуковая область сосредоточена только за скачком Маха (МБ). Интересно отметить, что при более умеренных числах Маха в зоне (4) за отраженным скачком (ИЗ) имеется дозвуковая область течения. Правая граница этой зоны при Мя < 2,2 совпадает с первой характеристикой веера разрежения (ЕР), прошедшей через отраженный скачок (ЯБ). При изменении М^ от 2,2 до 2,66 граница дозвуковой зоны сдвигается к тройной точке (Т) до исчезновения области дозвукового течения (4). Таким образом, при рассматриваемых условиях (М«, = 1,8) мы имеем дело с устойчивым решением, содержащим основные типы стационарных газодинамических особенностей.

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

Для решения задачи используется близкая к равномерной расчетная сетка, соответствующая геометрии канала, что позволяет более четко выявить свойства методов. Размер расчетной области составляет 4x2. Размерность расчетных сеток составляет 400x200, 800x400.

Результаты моделирования установившегося течения с использованием методов кинетического расщепления вектора потока (КБУБ), расщепления разности вектора потока (РОБ), комбинированного алгоритма (РОБш) и комбинированного алгоритма с использованием адаптивного ограничителя (Р08ш+) в виде изолиний поля плотности

представлены на рис. 9-12. На рис. 13 представлено решение на сетке 80рх400, полученное с помощью метода Р08т+.

Можно отметить хорошее раскрытие картины течения в согласии с теоретическими представлениями для всех методов. При этом положение особенностей и перепады газодинамических параметров на них достаточно хорошо согласуются с теоретическими значениями. Подтверждено, что в разрешении особенностей течения на неадаптированной расчетной сетке большую роль играют многомерные свойства методов. В этом смысле показательно разрешение искривленного тангенциального разрыва. Метод типа Годунова в рассматриваемом его варианте обеспечивает хорошее разрешение тангенциального разрыва, но дает явный дефект решения за скачком Маха, а также существенные ошибки в пристеночной области за косым скачком уплотнения и веером разрежения. При расчете других режимов течения с аналогичной конфигурацией первый дефект может не проявляться. Метод кинетического расщепления вектора потока обеспечивает надежный расчет интенсивных многомерных скачков, но страдает от излишней диссипации на тангенциальном разрыве. Качество сеточного решения с точки зрения раскрытия особенностей течения при расчете с помощью комбинированного алгоритма оказывается лучшим, чем для каждого из методов по отдельности. По-видимому, выбор критерия переключения (5) в комбинированном алгоритме можно считать разумным. Применение адаптивной функции-ограничителя (3) с Ь = 1,3 для Р\{Х) и с Ь= 1,35 для Рг{Х) дает существенное улучшение точности.

В выводах к главе 3 диссертации отмечены наиболее значимые результаты методического исследования: 1) способ конструирования аппроксимации вектора потока в численном методе для уравнений газовой динамики имеет принципиальное значение с точки зрения точности и «физичности» получаемого сеточного решения; 2) при расчете скачка и волны разрежения разница между получаемыми по разным методам решениями в большинстве случаев оказывается незначительной, при расчете контактного и тангенциального разрывов метод типа Годунова дает лучшее качество решения, тогда как методы расщепления вектора потока страдают от излишней диссипации; 3) диссипация методов расщепления вектора потока оказывается существенно большей, что согласуется с анализом главы 2; 4) задачи с интенсивными разрывами выдвигают экстремальные требования к конечно-разностным методам - в численных решениях появляются дефекты, общие для всех рассмотренных методов, для решения этих проблем требуется привлечение более общих моделей в областях столкновения ударных волн и глубокого разрежения; 5) разница в затратах компьютерных ресурсов на расчет задач с помощью рассматриваемых методов невелика - в пределах 10%, что связано с тем, что основная вычислительная нагрузка приходится на вычисление реконструированных значений переменных с учетом граничных условий,

поэтому вычислительная эффективность определяется в основном мерой ошибки решения; 6) с точки зрения расчета течений, содержащих тонкие детали на фоне неоднородного поля предпочтительно выглядят методы, учитывающие все характеристические свойства уравнений газовой динамики - методы типа Годунова, однако их недостатком является механизм многомерной неустойчивости при расчете интенсивных ударных волн, построенный комбинированный алгоритм позволяет удовлетворительно решить проблему методов типа Годунова; 7) использование адаптивной функции-ограничителя позволяет улучшить разрешающие свойства методов и точность решений при сохранении надежности и практически без увеличения вычислительных затрат и ужесточения требований к расчетной сетке, что повышает вычислительную эфективность методов.

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

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

Рассматривается обтекание простого симметричного 12% профиля, состоящего из двух дуг окружности. Параметры задачи следующие: невязкий сжимаемый совершенный газ с показателем адиабаты у = 1,4, значения числа Маха набегающего потока берутся равными Мж = 0,5; 0,65; 0,8, угол атаки а = 0°, вихревая структура является плоским изоэнтропическим вихрем с диаметром (по полю скорости), близким к хорде профиля, и давлением и плотностью в центре, составляющими около 0,9 от значений в невозмущенном потоке. Так как число Маха является фактором, который оказывает наибольшее влияние на взаимодействие и характеристики акустического поля, представляется необходимым исследовать это влияние в первую очередь. Максимальное число Маха течения, создаваемого вихрем, составляет 0,27.

Расчетная область взята в виде квадрата с размерами 10x10 хорд профиля, близкая к равномерной расчетная сетка типа «Н» состоит из

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

Первый и второй режимы течения (М*, = 0,5; 0,65) являются дозвуковыми с существенным влиянием сжимаемости. Третий режим (М„ = 0,8) трансзвуковой с замкнутыми изолированными сверхзвуковыми зонами на верхней и нижней поверхностях профиля. Для режима М» = 0,8 приведены четыре картины изолиний давления на последовательных этапах взаимодействия (рис. 14). Для всех режимов приведены зависимости от времени (в единицах Ыст, I - хорда профиля) аэродинамических коэффициентов Сх, Су, и М:, позволяющие получить представление о силовом воздействии на профиль (рис. 15).

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

различие структуры волновых возмущений, идущих от разных сторон профиля: разрежение-сжатие-разрежение (волны I типа) и сжатие-разрежение-сжатие (волны II типа), волны II типа имеют более отчетливую форму.

Анализ полученных данных позволяет определить базовые частотные и амплитудные характеристики. Значения амплитуд А\, Ац, длин волн Хь Хп и частот/ближнего поля представлены в таблице 1. Волны I и II типов отличаются амплитудами, но мало отличаются по длине. Длина волнового пакета варьируется вдоль его фронта, несколько возрастая в направлении от носка профиля к задней кромке, и уменьшается с ростом числа Маха. Амплитуда растет с ростом числа Маха, при этом амплитуда волн I типа больше, чем у волн II типа. Направление распространения имеет значение: максимальную амплитуду имеет часть волнового пакета, идущая в приблизительно поперечном к течению направлении.

Таблица 1

Параметры ближнего поля возмущений

Мое 0,5 0,65 0,8

Х\И = хн/ь 2-3 1.8-5-2.5 1.5+2

А,(%р„) =3% 24% =5%

А,, (%д,) =\% =2% £3%

/,Гц 550+850 650-950 8504-1150

Диапазон частот можно оценить исходя из скоростей распространения и масштабов длины. Для оценки принималось, что в системе координат, связанной с потоком, эти скорости близки к скорости звука 340 м/с. Хорда профиля для оценки взята равной 0,2 м. Как видно из таблицы, имеет место рост частоты волн при приближении числа Маха к единице, что обусловлено проявлением нелинейных эффектов.

Особенности силового воздействия на профиль заключаются в переменных знаках сил и моментов. При направлении вращения вихря против часовой стрелки Су и Mz на первом этапе взаимодействия испытывают существенные положительные приращения, тогда как при переходе ко второму этапу их знак резко меняется. Амплитуда Су на втором этапе меньше, амплитуда Mz, наоборот, увеличивается. Коэффициент Сх испытывает колебания противоположного знака с ощутимой относительной амплитудой. Релаксационный процесс возврата к невозмущенным значениям происходит за время ~L/(cK-Va). При увеличении числа Маха процесс взаимодействия протекает быстрее, при этом пиковые значения аэродинамических коэффициентов уменьшаются, хотя величины аэродинамических сил увеличиваются. Результат действия значительных сил и моментов до некоторой степени сглаживается кратковременностью явления (-0,001 с) и разными знаками воздействия.

Данные о значениях аэродинамических коэффициентов в ходе взаимодействия весьма ограничены, в том числе из-за быстроты протекания процесса. Качественно полученная картина (включая оценки амплитуд и частот поля возмущений) согласуется с данными исследований S.Lee, D.Bershader (1994), G.Djambazov, C.-H.Lai, A.Pericleous (1999), C. Loh, L. Hultgren, P. Jorgenson (1999), F.Schmitz, B.Sim (2001). Количественная оценка аэродинамических сил и моментов согласуется с данными работы Т. Wood, S. Grace (2001).

В выводах к главе 4 изложены основные результаты проведенного исследования: 1) взаимодействие носит нелинейный характер, начальный вихрь при нулевом прицельном расстоянии практически полностью разрушается; 2) число Маха является одной из определяющих величин взаимодействия; 3) выделены два типа волновых возмущений, образующихся в результате взаимодействия, которые отличаются по своей структуре и амплитуде; 4) полученные зависимости аэродинамических коэффициентов от времени позволяют подробно проследить динамику силового воздействия на профиль; 5) на основе результатов проводимых в настоящее время исследований можно предположить, что в качестве способа нейтрализации таких воздействий будут эффективны активные механические или тепловыделяющие устройства.

В заключении изложены главные результаты проведенного исследования.

ОСНОВНЫЕ РЕЗУЛЬТАТЫ РАБОТЫ

1. Систематизированы диссипативные свойства методов типа Годунова и методов расщепления вектора потока. Аналитическое и « численное исследования позволяют определить преимущества и недостатки методов и область их применения.

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

3. Построен вариант комбинированного алгоритма вычисления аппроксимации вектора потока, основанный на избирательном применении метода кинетического расщепления вектора потока и метода типа Годунова.

4. Проанализированы постановки численных граничных условий, обеспечивающие получение корректных решений.

5. Численно решен ряд одномерных и двумерных задач, что дает материал для сравнения методов и дальнейшего их совершенствования.

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

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

Основные результаты диссертации опубликованы в следующих работах:

1. В.В. Войков, И.В. Воронич Об одном подходе к построению кинетически-согласованных разностных схем газовой динамики. Численное моделирование в задачах аэродинамики и экологии: междуведомств, сб. науч. тр. М.: МФТИ, 1998, с. 3.

2. И.В. Воронич, В.В. Попов Статистическое моделирование невязких течений. Труды XIV Международной Школы по моделям механики сплошной среды, М: МФТИ, 1998, с. 39-45.

3.1. Voronich, V. Popov, Yu. Khlopkov Kinetic Calculation Model for Gas Dynamics. Proc. of 3rd Seminar on RRDPAE'98, Warsaw, Warsaw University of Technology Press, 1998, pp. 167-171.

4. И.В. Воронин Об эффективности схем с кинетическим расщеплением потока. Тезисы докладов Школы-семинар молодых ученых и специалистов «Актуальные проблемы аэрокосмической науки», Жуковский, ЦАГИ, 2001, с. 17-18.

5. И.В. Воронич Об особенностях моделирования динамики вихрей в сжимамом газе с помощью численных иметодов повышенной точности. Труды XLV научной конференции МФТИ (ГУ), ч. VI, Москва, 2002, с. 64.

6. И.В. Воронич Численное моделирование излучения акустических волн при взаимодействии вихревой структуры с аэродинамическим профилем в потоке сжимаемого газа. Электронный многопредметный научный журнал «ИССЛЕДОВАНО В РОССИИ», 2003, с. 1600-1615; URL: http://zhurnal.ape.reIarn.ru/articles/2003/136.pdf.

7. И.В. Воронич Сравнительный анализ группы численных методов решения уравнений газовой динамики. М.: МФТИ, 2005. (в печати)

а) б)

в) г)

Рис. 1. Собственные значения д\Ё\/д(2 для методов: а) Стегера-Уорминга, б) ван Лира, в) кинетического расщепления вектора потока, г) зависимости Ыс = 0,2) в методе Роу

Рис. 2. Функции-ограничители

От ВСТЬ Чж» Мм РЧ

> — Г1 1

11 08 - • Ю6 РОаг

* Р05т

о РОЗт - 1 1- О Ийп Ышт> 1

- Ими Мкш Г" я я*

/

/

• I -0? — /

ч 1 л.

«

0 " 0 ис. 3.1 1 " 5аспрс да * 0 пен ь " ие пл< ЗТН » ОС! ИИ 5 " 0 числа 2 " 0 Маха 4 ' Дл 1 0 138 1 9 ' Па» сса

для задачи об интенсивном разрыве

Пм) иг 4пЗ л

41 ^ * №У5 * РОЭ » РОЭп ■ V«! 1* Г о Ъша* м* * Нлт- V у Т1

Мш

• №\Г5

0 РСбт ««и Эм» 'чштт

гсб» мм*а

Рис. 5. Распределение плотности и числа Маха для задачи Шу-Ошера

в--

/

1 1

1 9

-

01 02 03

6« 0Г 08 09

06 07 0« 09

а)

б)

Рис. 6. Изолинии полей плотности в задаче о переносе вихря, полученных с помощью методов КТУ8 (а) и РББ (б)

а) б)

Рис. 7. Изолинии полей давления в задаче о потенциальном обтекании цилиндра, полученных с помощью методов КРУБ (а) и ГЭБ (б) в сравнении с теоретическим решением

Рис. 8. Схема течения в канале

04 1 1 5 2 29 3 35

Рис. 9. Изолинии поля плотности для метода КИУЗ

05 1 1 5 2 25 3 35

Рис. 12. Изолинии поля плотности для метода Р08т+

Рис. 13. Изолинии поля плотности для метода РЭ5т+ на сетке 800x400

Рис. 14. Изолинии поля давления в последовательные моменты времени в задаче о взаимодействии вихря с профилем (МЛ = 0,8)

Рис. 15. Зависимости Сх, Су, и Мот времени в задаче о взаимодействии вихря с профилем: a) Mw = 0,5, б) М„ = 0,65, в) Мж = 0,8

Подписано в печать 26.09.2005. Формат 60x90 1/16. Печать цифровая. Бумага «РегРэгтег». Печ. л. 2,0. Тираж 60. Заказ 6891.

Отпечатано в ФГУП «Производственно-издательский комбинат ВИНИТИ», 140010, г. Люберцы Московской обл., Октябрьский пр-т, 403. Тел. 554-21-86.

РНБ Русский фонд

2006-4 15914

Оглавление автор диссертации — кандидата физико-математических наук Воронич, Иван Викторович

Обозначения.

Введение.

1. Законы сохранения и основа численного метода.

1.1. Законы сохранения и их характеристические свойства.

1.2. Полудискретная задача.

1.3. Условие энтропии.

2. Численные методы решения уравнений Эйлера.

Введение.

2.1. Метод расщепления разности вектора потока.

2.2. Методы расщепления вектора потока.

2.2.1. Расщепление Стегера-Уорминга.

2.2.2. Расщепление ван Лира.

2.2.3. Кинетическое расщепление вектора потока.

2.3. Метод прямого статистического моделирования.

2.4. Энтропийные свойства методов.

2.5. Повышение точности базовых методов.

2.6. Комбинированный метод.

3. Расчет различных классов течений.

3.1. Одномерные тестовые задачи

3.2. Двумерные задачи.

Выводы.

4. Решение задачи о взаимодействии вихря с аэродинамическим профилем в потоке сжимаемого газа.

4.1. Постановка задачи.

4.2. Условия расчета и результаты.

4.3. Выводы.

Введение 2005 год, диссертация по информатике, вычислительной технике и управлению, Воронич, Иван Викторович

Численное моделирование движения газов и жидкостей является областью, активно развивающейся как в связи с растущим количеством практических приложений, так и в связи с увеличением производительности компьютерной техники и возрастающей ролью вычислительного эксперимента. Сложность представляет не только построение моделей движения жидкостей и газов (например, для турбулентных и многофазных течений), но и построение согласованных дискретных моделей, позволяющих получить корректное численное решение. Разработка точных, надежных и вычислительно эффективных численных алгоритмов, позволяющих моделировать течения в широком диапазоне режимов, является предметом обширных исследований [113].

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

При более чем 50-летней истории развития численных методов решения гиперболических систем уравнений [16-24], значительный прогресс в области построения численных методов решения уравнений Эйлера произошел в 1970-е - 1980-е годы [25-49]. Это можно объяснить востребованностью моделей такого уровня для решения практических задач, во многом благодаря росту производительности вычислительной техники.

Существенно, что с последующим развитием моделирования течений газов и жидкостей возросли требования к точности и надежности численных методов. В связи с этим развитие методов продолжилось, с новыми требованиями и преодолением новых проблем [14,15, 50-84]. Можно было бы думать, что тематика построения численных методов для уравнений Эйлера исчерпана, однако количество публикаций в этом направлении говорит о том, что это не так: ряд вопросов еще находится в стадии обсуждения. Некоторые из этих вопросов рассматриваются ниже.

Представляется важным выделить основные направления совершенствования численных методов решения уравнений газовой динамики. Во-первых, это способ построения аппроксимации вектора потока через грань расчетной ячейки в рамках метода конечного объема с учетом характеристических свойств законов сохранения и условия энтропии [11,42]. Такая аппроксимация должна обеспечивать хорошее разрешение разрывов (в частности, контактного и тангенциального, что важно для сдвиговых течений) и обладать достаточной надежностью при расчете интенсивных скачков и волн разрежения. Можно отметить широкое использование кинетических представлений для построения аппроксимации вектора потока применительно к моделированию течений сплошных сред и слаборазреженных газов [14, 33, 56, 57, 65, 75].

Во-вторых, большую роль играет способ повышения точности алгоритмов по пространственным переменным и времени. Эти два направления в настоящее время зачастую рассматриваются раздельно путем сведения к полудискретной задаче, представляющей систему обыкновенных дифференциальных уравнений по времени для переменных в ячейке с правой частью, полученной в результате дискретизации пространственных потоков [12]. Методы, использующие совместную аппроксимацию производных по пространственным переменным и времени, также продолжают развиваться [3,11,22]. Методам повышения пространственной точности посвящено большое количество работ [1-8,11], этой теме продолжает уделяться внимание. В настоящей работе рассматриваются методы второго (условно) порядка аппроксимации, потенциал которых представляется не исчерпанным. Такие методы широко используются в инженерной и научной практике как предсказуемые, терпимые к качеству расчетной сетки и допускающие более простые способы реализации граничных условий [7, 68]. Для интегрирования по времени разработаны группы как явных, так и неявных методов [12,47,48,53]. Для расчета нестационарных течений часто используются многошаговые методы Рунге-Кутта [83].

В-третьих, важным элементом вычислительного алгоритма является постановка и численная реализация граничных условий, так как точность и характер сеточного решения сильно зависят от реализации граничных условий [68,85-87]. С этой точки зрения можно выделить граничные условия на твердых стенках и внешних границах. Корректная постановка численных граничных условий позволяет избежать не только увеличения объема расчетной сетки и вычислительных затрат, но и получения ошибочного решения.

Расчетная сетка является необходимой компонентой процесса получения численного решения задачи. Желательно иметь адаптированную к геометрии задачи и к особенностям решения расчетную сетку при соблюдении ее качества. Адаптация расчетной сетки призвана способствовать повышению точности сеточного решения и ускорению сходимости расчета для стационарных задач [88-90] наряду с многосеточными методами [91] и методами предобуславливания [76].

Таким образом, вычислительный алгоритм включает в себя аппроксимацию вектора потока в системе законов сохранения, метод повышения точности по пространственным переменным, метод интегрирования по времени, численную реализацию граничных условий. Современные комплексы программ включают множество взаимосвязанных компонент, отвечающих за реализацию моделей различных процессов. Это осложняет верификацию таких комплексов и требует тестирования отдельных компонент [92]. В настоящей работе это учитывается путем анализа отдельных компонент алгоритма и их совместного тестирования. Этому способствует реализация алгоритмов с помощью объектно-ориентированных языков программирования. В случае, когда возникают сомнения в природе того или иного дефекта сеточного решения, необходимо также проводить расчеты с помощью методов первого порядка аппроксимации по пространству и времени.

Подход сквозного счета разрывных решений без явного выделения разрывов не свободен от недостатков. Прежде всего, это потеря точности сеточного решения в области за скачком, где традиционные методы сквозного счета обеспечивают не более чем первый порядок аппроксимации [34, 68,79]. Для преодоления этой проблемы разрабатывались методы выделения разрывов и адаптации расчетных сеток к разрывам [34]. Однако в случае движущихся и зарождающихся разрывов в многомерных течениях избежать снижения точности трудно. Альтернативой являются методы, построенные на базе понятия слабой аппроксимации, которые обеспечивают более чем первый порядок аппроксимации в области за скачком [79]. Недостатком таких методов является наличие нефизических осцилляций решения вблизи скачка. В настоящей работе используется подход сквозного счета разрывных решений при учете локальной картины течения для повышения разрешающих свойств и надежности методов. Подход сквозного счета позволяет достаточно полно выявить свойства аппроксимации вектора потока применительно к расчету различных типов течений на разнообразных расчетных сетках.

Ни один из вычислительно эффективных алгоритмов не лишен недостатков, как в силу конструктивных особенностей (связанных с аппроксимацией потоков), так и в силу трудностей, связанных с повышением точности по пространственным переменным и времени при ограничениях условия энтропии. Требования надежности и высокого разрешения являются в некотором смысле противоречивыми, так как требование повышения точности ведет к понижению схемной диссипации метода. Этот вопрос зачастую решается с помощью методов переменного порядка аппроксимации [21, 25,42]. Однако в этой связи можно отметить и востребованность комбинированных алгоритмов разного рода, так как разрешающие свойства зависят в первую очередь от базового метода первого порядка [14,55,58,61,74]. Под комбинированным алгоритмом здесь понимается такой алгоритм, который сочетает различные подходы к вычислению аппроксимации вектора потока в зависимости от локальной картины течения. Такие меры нужны в основном при расчете течений, содержащих интенсивные скачки и волны разрежения, так как методы, хорошо разрешающие контактные и тангенциальные разрывы (в первую очередь, методы типа Годунова), могут давать серьезные дефекты на интенсивных пространственных скачках [15, 55, 61, 74]. Традиционные методы расщепления вектора потока, наоборот, лучше пригодны для интенсивных скачков, но плохо разрешают контактные и тангенциальные разрывы, что важно для высокоскоростных сдвиговых течений [60].

Нужно отметить, что вопросы численного расчета многомерных течений с интенсивными скачками уплотнения с помощью методов сквозного счета являются дискуссионными до настоящего времени. Это связано с тем, что вопросы устойчивости ударных волн не исследованы в полной мере, поэтому не всякая неустойчивость численного расчета должна рассматриваться как дефект метода [15, 84, 93, 94].

Идеализированной целью является построение алгоритма, который имеет высокую точность в областях гладкого решения и хорошее разрешение разрывов (включая их положение и перепад параметров потока на разрыве) при сквозном расчете разрывных решений [50, 79]. По указанным выше причинам построение такого алгоритма невозможно без некоторого компромисса.

Представляется, что методы построения аппроксимации вектора потока для уравнений Эйлера можно классифицировать на базе понятия схемной вязкости, связанной со схемным диссипативным потоком - механизмом, регулирующим «физичность» получаемых сеточных решений [14,42, 50, 63]. Систематическое сравнительное исследование группы известных методов расщепления вектора потока и метода расщепления разности вектора потока (типа Годунова), использующих характеристические свойства законов сохранения, полезно для выделения взаимосвязи конструктивных особенностей этих широко используемых методов с качеством сеточных решений. Сравнительный анализ методов на основе свойств собственных значений схемного диссипативного потока и результатов решения ряда задач позволяет оценить различные подходы к построению методов и уточнить их область применимости [31, 51, 63].

Условие энтропии является существенным критерием при построении и анализе численных методов расчета обобщенных (разрывных) решений газодинамических задач [19,26,42]. Конструкция аппроксимации вектора потока определяет разрешающие и энтропийные свойства метода. Однако установить в общем случае строгую взаимосвязь между характеристиками схемной диссипации и энтропийными свойствами затруднительно. Действие схемной диссипации для энтропийно согласованного метода похоже на действие физической вязкости: процесс сходимости сеточного решения на последовательности измельчающихся расчетных сеток можно представить как предельный переход при стремлении физической вязкости к нулю. Настройка аппроксимации вектора потока позволяет в ряде случаев добиться улучшения качества сеточного решения [42, 58].

Для иллюстрации подходов к повышению точности по пространственным переменным рассматривается одномерный скалярный закон сохранения, для которого строится нелинейная разностная схема повышенной точности, удовлетворяющая некоторому условию, например условию сохранения монотонности или родственным условиям [25,44,49,67]. Нелинейность разностной схемы обусловлена применением функций-ограничителей приращений переменных для обеспечения корректной реконструкции распределений переменных внутри ячейки. Далее подход, развитый для скалярного случая, обобщается на полную систему уравнений. Удовлетворить строго условию энтропии в этом случае не представляется возможным, тем не менее, можно надеяться, что сеточное решение будет высокого качества во всей области течения, за исключением областей разрывов и снижения точности в областях их влияния. Для расчета многомерных течений существенен способ обобщения одномерной модели, как с точки зрения аппроксимации вектора потока [48,53,60,63], так и реконструкции распределений переменных в ячейке [28,29,45]. Для решения этих задач используются преобразования координат [48], методы факторизации конечно-разностных операторов [53], или обобщение на случай произвольной ориентации грани ячейки [45, 60, 63]. В настоящей работе рассматривается последний вариант как вполне естественный.

Роль линейной модели переноса показательна, так как она используется для описания движения среды на различных уровнях [10, 14,16,25]. Это относится к уравнениям движения сплошной среды и к кинетическим уравнениям как базису моделей сплошной среды. Метод прямого статистического моделирования является полноценной иллюстрацией этих взаимосвязей [95]. Существенен в этой связи также смысл условия энтропии с кинетической точки зрения.

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

Применение адаптивных функций-ограничителей рассматривается как способ повышения пространственной точности алгоритма [59,78,82].

Особенностью такой функции-ограничителя является то, что в зависимости от локальной волновой картины в распределениях различных переменных в некоторых случаях применяется центрально-разностный подход к вычислению значений переменных на грани ячейки вместо «монотонизации» решения по принципу минимальных приращений [25,42]. Нужно заметить, что такой метод отличается от применения схемы центральных разностей в областях гладкого решения [11]. Результаты показывают, что применение адаптивных функций-ограничителей повышает разрешающие свойства известных методов без ужесточения требований к качеству расчетной сетки, увеличивает их вычислительную эффективность с сохранением надежности [81].

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

Целью настоящей работы явилось построение вычислительно эффективного, надежного алгоритма решения нестационарных задач газовой динамики в широком диапазоне режимов, пригодного для раскрытия нестационарного волнового поля течения, содержащего разномасштабные особенности. Это требует систематизации свойств конечно-разностных методов (в том числе построенных с использованием кинетических моделей) расчета обобщенных (разрывных) решений уравнений газовой динамики с целью построения возможно более универсального алгоритма. Для решения этой задачи также необходимо исследование возможностей повышения разрешающих свойств методов с помощью аппарата функций-ограничителей без существенного ужесточения требований к расчетной сетке.

На защиту выносятся следующие положения: 1. результаты сравнительного изучения основных групп численных методов решения уравнений Эйлера, позволяющие определить область применения методов и уточнить механизмы схемной диссипации;

2. комбинированный алгоритм вычисления аппроксимации вектора потока, сочетающий применение двух аппроксимаций вектора потока, что позволяет расширить диапазон применения методов типа Годунова;

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

4. результаты численного решения задачи о взаимодействии вихря с профилем в дозвуковых и трансзвуковых потоках сжимаемого газа.

Научная новизна работы заключается в следующем:

1. получены новые данные о механизмах и характеристиках схемной диссипации методов типа Годунова и методов расщепления вектора потока;

2. создание нового комбинированного алгоритма вычисления аппроксимации вектора потока, позволяющего расширить диапазон применения методов типа Годунова;

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

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

Первая глава диссертации является вводной. В ней рассматриваются свойства законов сохранения газовой динамики, дискретизация по принципу конечного объема, полудискретная задача, понятия сходимости, аппроксимации, обобщенного решения, задача Коши для гиперболической системы законов сохранения, условие энтропии. Дается краткий обзор проблем построения универсальных методов повышенного порядка аппроксимации и формулируются задачи, решение которых дает возможность продвинуться в этом направлении.

Вторая глава посвящена детальному разбору методов. В ней рассматриваются следующие методы построения аппроксимации вектора потока (п. 2.1,2.2): метод типа Годунова - метод расщепления разности вектора потока, методы расщепления вектора потока - Стегера-Уорминга, ван Лира, кинетического расщепления вектора потока. Для этих методов рассмотрены свойства их схемной диссипации и способы энтропийной коррекции, сделаны выводы о характеристиках методов применительно к расчету различных классов течений. Для сопоставления континуального и кинетического подходов к построению численных методов расчета задач газовой динамики приведен анализ метода прямого статистического моделирования (п. 2.3). Обсуждению условия энтропии с кинетической точки зрения и установлению энтропийных свойств метода кинетического расщепления вектора потока посвящен п. 2.4. В п. 2.5 изложен подход к повышению точности методов по пространственным переменным и времени на базе условия сохранения монотонности. Обсуждается аппарат функций-ограничителей, на основе анализа существующих методов предложена адаптивная к локальному распределению переменных функция-ограничитель. Построению комбинированного алгоритма вычисления аппроксимации вектора потока с учетом существующих подходов посвящен п. 2.6.

В третьей главе изложены результаты применения описанных методов к решению ряда одномерных и двумерных по пространственным переменным задач. Одномерные задачи (п. 3.1) позволяют проверить выводы, сделанные ранее о свойствах схемной диссипации методов, а также проверить пригодность комбинированного алгоритма и адаптивной функции-ограничителя на различных газодинамических конфигурациях. В п. 3.2 рассматриваются 3 стационарные и нестационарные двумерные задачи. Разбираются вопросы построения расчетных сеток, реализации граничных условий и оценки точности и сходимости сеточных решений. Результаты решения двумерных задач позволяют подтвердить выводы о схемной диссипации методов, подтвердить эффективность комбинированного алгоритма и количественно проверить эффективность адаптивной функции-ограничителя.

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

В разделе заключение представлены основные результаты проведенных исследований.

Автор выражает признательность: научному руководителю - заведующему кафедрой компьютерного моделирования МФТИ профессору Ю.И. Хлопкову за руководство и поддержку исследований, доцентам В.А. Жарову и C.J1. Горелову за обсуждение постановок задач и результатов, профессору А.И. Толстых за обсуждение ряда методических вопросов, к. ф.-м. н. В.В. Власенко за ценные рекомендации и помощь в исправлении ряда неточностей,

В.В. Войкову за участие в совместных исследованиях, В.В. Ткаченко за помощь в подготовке работы, сотрудникам кафедры компьютерного моделирования МФТИ и лично М.В. Спиркиной за поддержку, своим родителям за под держку и терпение.

Работа выполнялась при содействии Государственной программы поддержки ведущих научных школ, грант НШ-1984.2003.1, руководитель -профессор М.Н. Коган.

Заключение диссертация на тему "Методы решения нестационарных задач газовой динамики"

Выводы

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

Выделены два типа акустических волн, образующихся в результате взаимодействия, которые отличаются по своей структуре и амплитуде. Оценки амплитуды и частоты ближнего поля возмущений для указанных режимов составляют ~1 КПа и ~1 КГц соответственно, амплитуда и частота растут с ростом числа Маха, рост частоты является нелинейным. В поперечном к потоку направлении амплитуда волн максимальна. Эти данные согласуются с имеющимися представлениями о спектре и амплитуде аэродинамических звуковых источников.

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

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

Заключение

Основными результатами исследования являются:

1. Систематизированы диссипативные свойства методов типа Годунова и методов расщепления вектора потока. Аналитическое и численное исследования позволяют определить преимущества и недостатки методов и область их применения.

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

3. Построен вариант комбинированного алгоритма вычисления аппроксимации вектора потока, основанный на применении диссипативного метода расщепления вектора потока вместо метода типа Годунова в окрестности интенсивных скачков.

4. Проанализированы постановки численных граничных условий, обеспечивающие получение корректных сеточных решений.

5. Численно решен ряд одномерных и двумерных задач, что дает материал для сравнения методов и дальнейшего их совершенствования.

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

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

Библиография Воронич, Иван Викторович, диссертация по теме Математическое моделирование, численные методы и комплексы программ

1. Б.Л. Рождественский, Н.Н. Яненко Системы квазилинейных уравнений. М.: Наука, 1978.

2. А.А. Самарский, Ю.П. Попов Разностные методы решения задач газовой динамики, (изд. 4-е) М.: УРСС, 2004.

3. К.М. Магомедов, А.С. Холодов Сеточно-характеристические численные методы. М.: Наука, 1988.

4. А.И. Толстых Компактные разностные схемы и их применение в задачах аэрогидродинамики. М.: Наука, 1990.

5. О.М. Белоцерковский Численное моделирование в механике сплошных сред (изд. 2-е). М.: Физматлит, 1994.

6. R. LeVeque Numerical Methods for Conservation Laws (2nd ed.). Basel: Birkhauser, 1996.

7. A. Jameson Essential Elements of Computational Algorithms for Aerodynamic Analysis and Design. ICASE Report N 97-68, 1997.

8. Upwind and High-Resolution Schemes. John Van Rosendale (ed.) N.Y.: Springer, 1997.

9. D. Wilcox Turbulence Modeling for CFD. DCW Industries, 1998.

10. О.М. Белоцерковский, A.M. Опарин Численный эксперимент в турбулентности: от порядка к хаосу (изд. 2-е). М.: Наука, 2000.

11. А.Г. Куликовский, Н.В. Погорелов, А.Ю.Семенов Математические вопросы численного решения гиперболических систем уравнений. М.: Физматлит, 2001.

12. Н. Lomax, Т. Pulliam, D. Zingg Fundamentals of Computational Fluid Dynamics. Berlin: Springer-Verlag, 2001.

13. T. Chung Computational Fluid Dynamics. Cambridge: Cambridge University Press, 2002.

14. К. Xu, L. Martinelli, A. Jameson Gas-Kinetic Finite Volume Methods, Flux-Vector Splitting, and Artificial Diffusion. J. Computational Physics, 1995, V. 120, pp. 48-65.

15. J.-Ch. Robinet, J. Gressier, G. Casalis and J.-M. Moschetta Shock Wave Instability and the Carbuncle Phenomenon: Same Intrinsic Origin? J. Fluid Mechanics, 2000, V. 417, pp. 237-263.

16. R. Courant, E. Isaacson, M. Rees On the Solution of Nonlinear Hyperbolic Differential Equations by Finite Differences. Comm. Pure Appl. Math., 1952, V. 5, pp. 243-255.

17. P. Lax Weak Solutions of Nonlinear Hyperbolic Equations and Their Numerical Computation. Comm. Pure. Appl. Math., 1954, V. 7, pp. 159-194.

18. О.А.Ладыженская О построении разрывных решений квазилинейных гиперболических уравнений, как пределов решений соответствующих параболических уравнений при стремлении «коэффициента вязкости» к нулю. Доклады АН СССР, 1956, Т. 111, N 2, с. 291-294.

19. О.А. Олейник Разрывные решения нелинейных дифференциальных уравнений. Успехи математических наук, 1957, Т. 12, № 3, с. 3-73.

20. С.К. Годунов Разностный метод численного расчета разрывных решений уравнений гидродинамики. Математический сборник, 1959, Т. 3, № 47, с. 271.

21. Р.П. Федоренко Применение разностных схем высокого порядка точности для гиперболических уравнений. ЖВМ и МФ, 1962, Т. 2, № 6.

22. P. Lax, В. Wendroff Difference Schemes for Hyperbolic Equations with High Order of Accuracy. Comm. Pure Appl. Math., 1964, V. 17, pp. 381-398.

23. B.B. Русанов Конечно-разностные схемы третьего порядка точности для сквозного счета разрывных решений. Доклады АН СССР, 1968, Т. 180, № 6, с. 1303-1305.

24. R. MacCormak The Effect of Viscosity in Hypervelocity Impact Cratering. AIAA Paper, 1969, N 69-354.

25. В.П. Колган Применение принципа минимальных значений производной к построению конечно-разностных схем для расчета разрывных решений газовой динамики. Ученые записки ЦАГИ, 1972, Т. III, № 6, с. 68-77.

26. P. Lax Hyperbolic Systems of Conservation Laws and the Mathematical Theory of Shock Waves. Society for Industrial and Applied Mathematics, 1972.

27. B. van Leer Towards the Ultimate Conservative Difference Scheme II. J. Computational Physics, 1974, V. 14, pp. 361-370.

28. В.П. Колган Конечно-разностная схема для расчета двумерных разрывных решений нестационарной газовой динамики. Ученые записки ЦАГИ, 1975, Т. VI, № 1, с. 9-14.

29. С.К. Годунов, А.В. Забродин, М.Я. Иванов, А.Н. Крайко, Г.П. Прокопов Численное решение многомерных задач газовой динамики. М.: Наука, 1976.

30. В.П. Колган Применение сглаживающих операторов в конечно-разностных схемах высокого порядка точности. ЖВМ и МФ, 1978, Т. 18, №5, с. 1340-1345.

31. G. Sod A Survey of Several Finite Difference Methods for Systems of Nonlinear Hyperbolic Conservation Laws. J. Computational Physics, 1978, V. 27, pp. 131.

32. B. van Leer Towards the Ultimate Conservative Difference Scheme V. J. Computational Physics, 1979, V. 32, pp. 101-136.

33. D. Pullin Direct Simulation Methods for Compressible Inviscid Ideal-Gas Flow. J. Computational Physics, 1980, V. 34, N 2, p. 231.

34. А.Н. Крайко, B.E. Макаров, Н.И. Тилляева Численное построение фронтов ударных волн. ЖВМ и МФ, 1980, Т. 20, № 3, с. 716-723.

35. P. Roe Approximate Riemann Solvers, Parameter Vectors, and Difference Schemes. J. Computational Physics, 1981, V. 43, N 10, pp. 357-372.

36. J. Steger, R. Warming Flux Vector Splitting of the Inviscid Gas Dynamics Equations with Application to Finite-Difference Methods. J. Computational Physics, 1981, V. 40, p. 263.

37. В. van Leer Flux Vector Splitting for the Euler Equations. Lecture Notes in Physics, 1982, V. 170, pp. 307-312.

38. S. Osher, F. Solomon Upwind Difference Schemes for Hyperbolic Systems of Conservation Laws. Math. Comput., 1982, V. 38, pp. 339-374.

39. G. van Albada, B. van Leer, W. Roberts A Comparative Study of Computational Methods in Cosmic Gas Dynamics. Astronomy and Astrophysics, 1982, V. 108, pp. 76-84.

40. В.И. Копченов, A.H. Крайко Монотонная разностная схема второго порядка для гиперболических систем с двумя независимыми переменными. ЖВМ и МФ, 1983, Т. 23, № 4, с. 848-858.

41. М.И. Волчинская, А.Н.Павлов, Б.Н. Четверушкин Об одной схеме расчета газодинамических уравнений. Препринт № 113 ИПМ им. М.В. Келдыша АН СССР, 1983.

42. A. Harten High Resolution Schemes for Hyperbolic Conservation Laws. J. Computational Physics, 1983, V. 49, N 3, pp. 357-393.

43. A. Harten, P. Lax, B. van Leer On Upstream Differencing and Godunov-Type Schemes for Hyperbolic Conservation Laws. SI AM Review, 1983, V. 25, N 1, pp. 35-61.

44. P. Sweby High Resolution Scheme Using Flux Limiters for Hyperbolic Conservation Laws. SIAM J. Numer. Anal., 1984, V. 21, p. 995.

45. Н.И. Тилляева Обобщение модифицированной схемы С.К. Годунова на произвольные нерегулярные сетки. Ученые записки ЦАГИ, 1986, Т. XVII, № 2, с. 18-26.

46. S. Deshpande Kinetic Theory Based New Upwind Methods for Inviscid Compressible Flows. AIAA Paper, 1986, 86-0275.

47. A.B. Родионов Повышение порядка аппроксимации схемы С.К. Годунова. ЖВМ и МФ, 1987, Т. 27, № 12, с. 1853-1860.

48. Н. Yee, A. Harten Implicit TVD Schemes for Hyperbolic Conservation Laws in Curvilinear Coordinates. AIAA Journal, 1987, V. 25, p. 266.

49. A. Harten, S. Osher Uniformly High-Order Accurate Nonoscillatory Schemes. SIAM J. Numer. Analys., 1987, V. 24, N. 2, pp. 279-309.

50. Ю.Б. Радвогин О квазимонотонных разностных схемах второго порядка. Препринт ИПМ им. М.В. Келдыша АН СССР, 1991.

51. В. Einfeldt, С. Munz, P. Roe, В. Sjogreen On Godunov Type Methods Near Low Densities. J. Computational Physics, 1991, V. 92, pp. 273-295.

52. S. Lele Compact Finite Difference Schemes with Spectral-Like Resolution. J. Computational Physics, 1992, V. 103, p. 16.

53. T. Pulliam Solution Methods in Computational Fluid Dynamics. Yon Karman Institute for Fluid Mechanics Lecture Series, 1994.

54. A. Harten Adaptive Multiresolution Schemes for Shock Computations. J. Computational Physics, 1994, V. 115, pp. 319-338.

55. J. Quirk A Contribution to the Great Riemann Solver Debate. Int. J. Numer. Meth. Fluids, 1994, V. 18, pp. 555-574.

56. T. Elizarova, I. Graur, J. Lengrand, A. Chpoun Rarefied Gas Flow Simulation Based on Quasi Gas-Dynamic Equations. AIAA Journal, 1995, V. 33, N 12, pp. 2316-2324.

57. S. Brown Approximate Riemann Solvers for Moment Models of Dilute Gases. PhD Dissertation. The University of Michigan, 1996.

58. R. Donat, A. Marquina Capturing Shock Reflections: an Improved Flux Formula. J. Computational Physics, 1996, V. 125, pp. 42-58.

59. M. Arora, P. Roe A Well-Behaved TVD Limiter for High-Resolution Calculations of Unsteady Flow. J. Computational Physics, 1997, V. 132, pp. 311.

60. J-M. Moschetta, D. Pullin A Robust Low Diffusive Kinetic Scheme for the Navier-Stokes/Euler Equations. J. Computational Physics, 1997, V. 133, pp. 193-204.

61. Y. Wada, M.-S. Liou An Accurate and Robust Flux Splitting Scheme for Shock and Contact Discontinuities. SIAM J. Sci. Comput., 1997, V. 18, N 3, pp. 633657.

62. R. LeVeque Wave Propagation Algorithms for Multi-Dimensional Hyperbolic Systems. J. Computational Physics, 1997, V. 131, pp. 327-353.

63. W. Kleb Comments Regarding Two Upwind Methods for Solving for Two-Dimensional External Flows Using Unstructured Grids. NASA TM N 109078, 1997.

64. X. Li Entropy Consistent, TVD Methods With High Accuracy For Conservation Laws. Electronic Journal of Differential Equations, Conference 01, 1997, pp. 171-191. ISSN: 1072-6691. URL: http://ejde.math.swt.edu.

65. T. Lou, D. Dahlby, D. Baganoff A Numerical Study Comparing Kinetic Flux-Vector Splitting for the Navier-Stokes Equations with a Particle Method. J. Computational Physics, 1998, V. 145, pp. 489-510.

66. R. MacCormack, T. Pulliam Assessment of a New Numerical Procedure for Fluid Dynamics. 29th ALAA Fluid Dynamics Conference, Albuquerque, NM, June 15-18, 1998.

67. В.Г. Крупа О построении разностных схем повышенного порядка точности для гиперболических уравнений. ЖВМ и МФ, 1998, Т. 38, № 1, с. 85.

68. А.Н. Минайлос Точность численных решений уравнений Навъе-Стокса. ЖВМ и МФ, 1998, Т. 38, № 7, с. 1220-1232.

69. С.А. Величко, Ю.Б. Лифшиц, И.А. Солнцев Расчет нестационарных течений с помощью схемы повышенной точности. ЖВМ и МФ, 1999, Т. 39, № 5, с. 817.

70. В.А. Гаранжа, В.Н. Конынин Численные алгоритмы для течений вязкой жидкости, основанные на консервативных компактных схемах высокого порядка аппроксимации. ЖВМ и МФ, 1999, Т. 39, № 8. с. 1378-1392.

71. М.В.Липавский, А.И.Толстых О сравнительной эффективности схем с нецентрированными компактными аппроксимациями. ЖВМ и МФ, 1999, Т. 39, № 10, с. 1705-1720.

72. S. Lui, К. Xu Entropy Analysis of Kinetic Flux Vector Splitting Schemes for the Compressible Euler Equations. ICASE Report N 99-5, Institute for Computer Applications in Science and Engineering, NASA Langley Research Center, 1999.

73. H. Reksoprodjo, R. Agarwal A Higher-Order Kinetic Wawe-Particle Flux-Splitting Algorithm for the Euler Equations. AIAA Paper 99-0920, 1999.

74. K. Xu Gas Evolution Dynamics in Godunov-Type Schemes and Analysis of Numerical Shock Instability. 1999, ICASE Report N 99-6.

75. Б.Н. Четверушкин Кинетически-согласованные разностные схемы в газовой динамике. М.: Изд-во МГУ, 1999.

76. Е. Turkel Preconditional Techniques in Computational Fluid Dynamics. Annual Rev. Fluid Mech., 1999, V. 31, pp. 385^116.

77. H. Yee, M. Vinokur, M. Djomehri Entropy Splitting and Numerical Dissipation. J. Computational Physics, 2000, V. 162, pp. 33-81.

78. А.Н. Гильманов Локально-характеристический подход в разностной схеме повышенного порядка аппроксимации. ЖВМ и МФ, 2000, Т. 40, № 4, с. 557-562.

79. В.В. Остапенко О построении разностных схем повышенной точности для сквозного расчета нестационарных ударных волн. ЖВМ и МФ, 2000, Т. 40, № 12, с. 1857-1875.

80. И.В. Абалакин, А.В. Жохова, Б.Н. Четверушкин Разностные схемы на основе кинетического расщепления вектора потока. Математическое моделирование, 2000, Т. 12, № 4, с. 73-82.

81. A. Rezgui, P. Cinnella, A. Lerat Third-Order Accurate Finite Volume Schemes for Euler Computations on Curvilinear Meshes. Computers&Fluids, 2001, V. 30, pp. 875-901.

82. G. Billet, О. Louedin Adaptive Limiters for Improving the Accuracy of the MUSCL Approach for Unsteady Flows. J. Computational Physics, 2001, V. 170, pp. 161-183.

83. S. Gottlieb, C.-W. Shu, E. Tadmor Strong Stability-Preserving High-Order Time Discretization Methods. SIAM Review, 2001, V. 43, N 1, pp. 89-112.

84. Y. Chauvat, J.-M. Moschetta and J. Gressier Shock Wave Numerical Structure and the Carbuncle Phenomenon. Int. J. Numer. Meth. Fluids, 2005, V. 47, pp. 903-909.

85. T. Poinsot, S. Lele Boundary Conditions for Direct Simulations of Compressible Viscous Flows. J. Computational Physics, 1992, V. 101, pp. 104^129.

86. A. Dadone, B. Grossman Surface Boundary Conditions for the Numerical Solution of the Euler Equations. А1АА Journal, 1994, V. 32, N 2, pp. 285-293.

87. М.А. Ильгамов, А.Н. Гильманов Неотражающие условия на границах расчетной области. М.: Физматлит, 2003.

88. P. Eiseman A Multi-Surface Method of Coordinate Generation. J. Computational Physics, 1979, V. 33, pp. 118-150.

89. J. Thompson, Z. Warsi, C. Mastin Numerical Grid Generation: Foundations and Applications. North-Holland, Elsevier, 1985.

90. А.Н. Гильманов Методы адаптивных сеток в задачах газовой динамики. М.: Физматлит, 2000.

91. P. Wesseling An Introduction to Multigrid Methods. Philadelpha: R.T. Edwards, Inc., 2004.

92. W. Kleb, W. Wood CFD: A Castle in the SancP. А1АА Paper, 2004, N2627

93. Л.Д. Ландау, E.M. Лифшиц Теоретическая физика. Т. VI. Гидродинамика. М.: Наука, 1986.

94. И.А. Знаменская О классификации типов неустойчивости ударных волн. Труды XIV Международной Школы по моделям механики сплошной среды. М.: МФТИ, 1998.

95. G. Bird Molecular Gas Dynamics and Direct Simulation of Gas Flows. Oxford: Oxford University Press, 1994.

96. Н.Я. Фабрикант Аэродинамика. M.: Наука, 1964.

97. JI.B. Овсянников Лекции по основам газовой динамики. М.: Наука, 1981.

98. В.П. Стулов Лекции по газовой динамике. М.: Наука, 2004.

99. М.Н. Коган Динамика разреженного газа. М.: Наука, 1967.

100. М.С. Иванов, С.В. Рогазинский Метод прямого статистического моделирования в динамике разреженного газа. Новосибирск: ВЦ СО АН СССР, 1988.

101. М.Н. Коган, А.С.Кравчук, Ю.И. Хлопков Метод «релаксация-перенос» для решения задач динамики газа в широком диапазоне разреженности среды. Ученые записки ЦАГИ, 1988, Т. XIX, № 2, с. 106.

102. I. Voronich, V. Popov, Yu. Khlopkov Kinetic Calculation Model for Gas Dynamics. Proc. of 3rd Seminar on RRDPAE'98, Warsaw, Warsaw University of Technology Press, 1999.

103. Чой Ен-Ин Применение метода прямого статистического моделирования к расчету течений невязкого газа. Магистерская диссертация. М.: МФТИ, 2002.

104. Ю.И. Хлопков, C.J1. Горелов Методы Монте-Карло и их приложение в механике и аэродинамике. М.: МФТИ, 1989.

105. В.В. Войков, И.В. Воронич Об одном подходе к построению кинетически-согласованных разностных схем газовой динамики. Численное моделирование в задачах аэродинамики и экологии: междуведомств, сб. науч. тр. М.: МФТИ, 1998.

106. J. von Neumann Refraction, Intersection and Reflection of Shock Waves. NAVORD Rep. N 203-45, Washington, 1945.

107. G. Ben-Dor Shock Wave Reflection Phenomena. N.Y.: Springer-Verlag, 1991.

108. H. Li, G. Ben-Dor A Parametric Study of Mach Reflection in Steady Flows. J. Fluid Mechanics, 1997, V. 341, pp. 101-125.

109. M. Ivanov, D. Vandromme, V. Fomin, A. Kudryavtsev, A. Hadjadj, D. Khotyanovsky Transition Between Regular and Mach Reflection of Shock Waves: New Numerical and Experimental Results. Shock Waves, 2001, V. 11, pp. 199-207.

110. A.B. Кузьмина Применение численных методов расчета обтекания трехмерных конфигураций для решения задач интерференции двигателя и планера транспортного самолета. Магистерская диссертация. М.: МФТИ, 2000.

111. М. Ван-Дайк Альбом течений жидкости и газа. М.: Мир, 1986.

112. N. Peake, D. Crighton Active Control of Sound. Ann. Rev. Fluid Mech., 2000, V. 32, pp. 137-164.

113. A.C. Гиневский, Е.В.Власов, P.K. Каравосов Акустическое управление турбулентными струями. М.: Физматлит, 2001.

114. S. Lee, D. Bershader Head-On Parallel Blade-Vortex Interaction. AIAA Journal, 1994, V. 32, N 1, pp. 16-22.

115. F. Schmitz, B. Sim Acoustic Phasing, Directionality and Amplifcation Effects of Helicopter Blade-Vortex Interactions. J. American Helicopter Society, 2001, V. 46, N 10, pp. 273-282.

116. G. Djambazov, C.-H. Lai, A. Pericleous Sound Generation by Vortex-Blade Interactions. Proceedings of 11th Int. Conf. on Domain Decomposition Methods, 1999, pp. 422-429.

117. C. Loh, L. Hultgren, P. Jorgenson Vortex Dynamics Simulation in Aeroacoustics by the Space-Time Conservation Element and Solution Element Method. AIAA Paper, 1999, N 99-0359.

118. T. Wood, S. Grace Wing Geometry Effect on Blade-Vortex Interaction Response Using BEM. Proc. of 2001 ASME International Mechanical Engineering Congress and Exposition, New York, 11-16 November 2001.

119. P. Chen, J. Baeder, R.Evans, J. Niemczuk Blade-Vortex Interaction Noise Reduction with Active Twist Smart Rotor Technology. Smart Materials and Structures, 2001, N 10, pp. 77-85.