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

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

Бахвалов Павел Алексеевич

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

Специальность 05.13.18 — математическое моделирование, численные методы и комплексы программ

Автореферат диссертации на соискание учёной степени кандидата физико-математических наук

5 ДЕК 2013

005541463

Москва - 2013

005541463

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

Научный руководитель:

Официальные оппоненты:

доктор физико-математических наук, старший научный сотрудник Козубская Татьяна Константиновна

кандидат физико-математических наук Титарев Владимир Александрович Вычислительный центр им. А. А. Дородницына РА отдел механики, сектор кинетической теории газов, ведущий научный сотрудник кандидат физико-математических наук, старший научный сотрудник Миронов Михаил Арсеньевич, Акустический институт им. акад. Н. Н. Андреева, теоретическая лаборатория, начальник теоретической лаборатории

Ведущая организация: Институт вычислительной математики РАН

Защита состоится 2013г. в часов на за-

седании диссертационного совета Д 212.156.05 при Московском физико-техническом институте (государственном университете) по адресу 141700, Московская обл., г. Долгопрудный, Институтский пер., д. 9, ауд. 903 КПМ.

С диссертацией можно ознакомиться в библиотеке Московского физико-технического института (государственного университета).

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

_____2013 года.

Ученый секретарь

диссертационного совета Д 212.156.05

7

Федько О.С.

Общая характеристика работы

Актуальность темы.

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

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

Среди схем высокого порядка на неструктурированных сетках, пригодных для решения аэроакустическнх задач, активно развиваются два главных направления: конечно-элементные схемы, главным образом — метод Галёр-кина с разрывными базисными функциями (DG), и конечно-объёмные схемы, основанные на полиномиальной реконструкции переменных. Эти методы близки по своим свойствам. Они характеризуются возможностью построения схем сколь угодно высокого порядка точности, низкой чувствительностью к качеству сетки. Недостатком обоих классов схем является их ресурсоёмкость, особенно при необходимости воспроизведения разрывов.

В 1998 году A. Dervieux и С. Debiez предложили группу разностных схем под названием «Mixed-element-volume MUSCL methods with weak viscosity» (смешанные конечно-элементные конечно-объёмные низкодиссипа-тивные MUSCL-методы). Их идея заключалась в построении консервативной схемы 2-го порядка на произвольной неструктурированной сетке, которая в. случае декартовой сетки вырождалась бы в схему 5-го порядка. В той же работе были предложены TVD-моднфикации данной схемы. Методика была применена к двумерному уравнению переноса.

Обобщение данной методики на нелинейные задачи было предложено в работе A. Dervieux, И. А. Абалакпна и Т. К. Козубской в 2006 году в ра-

боте «High Accuracy Finite Volume Method for Solving Nonlinear Aeroacoustics Problems on Unstructured Meshes» (высокоточный конечно-объёмный метод для решения нелинейных аэроакустических задач на неструктурированных сетках). На декартовой сетке предложенная схема вырождалась в конечно-разностную схему 5-го порядка. Можно сказать, что эта работа открыла новый класс разностных схем, а именно, консервативных конечно-разностных схем повышенной точности на неструктурированных сетках. В последующих работах эти схемы были названы схемами с реконструкцией переменных вдоль направления ребра (Edge-Based Reconstruction, EBR), а некоторые из них — схемами с квазиодномерной реконструкцией.

Позднее в работах N. Gourvitch, G. Roge, В. Koobus, F. Alauzet, И. Абалакина, Т. Козубской полученный класс схем был обобщён на 3-мерный случай, также были продолжены исследования по построению лимитеров и сохранению неотрицательности плотности и давления газа. В работах В. Koobus, M.-V. Salvetti, S.Camarri, И. Абалакина, А. Горобца, А. Дубеня, Т. Козубской и др. этот класс схем был применён к промышленным задачам (моделирование звукопоглощающих конструкций, турбулентных течений в щелях и кавернах и др.).

В работах вышеперечисленных авторов было замечено, что рассматриваемые численные схемы обладают лучшими характеристиками, чем предсказывают известные аналитические оценки.

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

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

1) уточнение аналитических оценок точности существующей EBR-схемы;

2) получение экспериментальных оценок точности существующей EBR-схемы на новом наборе тестов;

3) построение более эффективных схем в рамках квазиодномерного подхода;

4) построение объёмно-центрированных схем с квазиодномерной реконструкцией переменных;

5) исследование характеристик акустического рупора с помощью разработанных разностных схем.

Научная новизна.

1. В настоящей работе аналитические оценки, выписанные А. Беппеих и С. БеЫег на декартовой сетке в случае орто- и барицентрических ячеек, обобщены на случай произвольных ячеек и произвольной сетки, полученной однородным разбиением параллелограммов (параллелепипедов).

2. Предложены новые формулировки схемы в рамках квазиодномерного подхода.

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

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

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

• научный семинар факультета управления и прикладной математики МФТИ (2012);

• семинар сектора вычислительной аэроакустики института прикладной математики им. М. В. Келдыша РАН (ИПМ им. М. В. Келдыша РАН) (2011, 2012);

• семинар «Математическое моделирование» ИПМ им. М. В. Келдыша РАН (2013),

а также на следующих конференциях:

• Международная научно-техническая конференция "Мехатроника, автоматизация, управление - 2009" (МАУ-2009) (Дивноморское (Геленджик), 2009);

• Всероссийская открытая конференции по авиационной акустике (Звенигород, 2009);

• Третья открытая всероссийская конференция "Вычислительный эксперимент в аэроакустике"(Светлогорск, 2010);

• 53-я научная конференция МФТИ "Современные проблемы фундаментальных и прикладных наук"(Долгопрудный, 2010);

• Вторая всероссийская открытая конференция по авиационной акустике (Звенигород, 2011);

• XIII международный семинар "Супервычнсления и математическое моделирование" (Саров, 2011);

• Parallel CFD 2011 (Barcelona, Spain, 2011);

• Четвёртая открытая всероссийская конференция "Вычислительный эксперимент в аэроакустике"(Светлогорск, 2012);

• Третья всероссийская открытая конференция по авиационной акустике (Звенигород, 2013);

• European Workshop on High Order Nonlinear Numerical Methods for Evolutionary PDEs: Theory and Applications (HONOM 2013), March 2013, Talence, Bordeaux, France.

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

Публикации. Основные результаты по теме диссертации изложены в 13 печатных работах, три из которых [2][3][4] изданы в журналах, рекомендованных ВАК.

На защиту выносятся основные результаты диссертационной работы.

Объем и структура работы. Диссертация состоит из введения, пяти глав, заключения и списка использованных источников. Полный объем диссертации 103 страницы текста с 37 рисунками и 4 таблицами. Список использованных источников содержит 104 наименования.

Диссертационная работа была поддержана Российским фондом фундаментальных исследований (проекты 12-01-00486-а, 12-01-33022).

Содержание работы

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

6

диссертации. Формулируются цель и задачи работы, её научная новизна и практическая значимость.

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

Параграф 1.1 является вспомогательным. В нём определяется базовое семейство конечно-разностных схем для решения одномерного уравнения переноса

ди ди

на равномерной сетке. Временная дискретизация проводится методом Рунге-Кутты. Аппроксимация пространственной производной записывается в консервативном виде, например, для схемы 5-го порядка,

ди дх

ч 1 13 47 9 1

= ^+1/2 = -и_2--и_1 + -гю+25«1--и2 (2)

Приводятся графики амплитудной и фазовой ошибки в зависимости от длины волны.

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

В параграфе 1.3 изучаются консервативные конечно-разностные схемы на неравномерной сетке для уравнения переноса. Описывается используемый в ЕВЛ схемах метод перехода к разделённым разностям, позволяющий обеспечить точность на линейной функции и при этом вырождение в схему (2) на равномерной сетке. Проиллюстрируем применение одного из возможных вариантов этого метода для построения схемы на 5-точечном шаблоне. Для этого перепишем реконструкцию (2) в следующем виде:

^•+1/2 = Щ + ~

-¿(«-1 - и-2) + - и-х) + - и0) - 1(«2 - щ)

(3)

после чего заменим первые разности разделёнными разностями, а коэффициент 1/2 — расстоянием от узла г до центра сегмента г + 1/2:

р- _ ^г+1 ~ Щ (__1 1 - Цг-2 11 Щ —

¿+1/2 - 2 ^ 15а.._1 _ х._2 + 30а;. _ х,_ +

| 4 им - щ _ 1 щ+2 - щ+1 \ 5 — 10 £¿+2 — £¿+1)

Проводится экспериментальное сравнение точности получаемых таким образом схем с определением на 3-, 5- и 7-точечном шаблонах, а также различные модификации схем в рамках метода разделённых разностей. Также проводится сравнение с неконсервативной схемой, в которой производная определяется путём дифференцирования интерполяционного полинома на скошенном шаблоне.

—■— 3-point reconstruction

-5-point reconstruction

—о— 3lh order polynomial

Average mesh step

—■— 3-point reconstruction

-5-point reconstruction

---- 7-point reconstruction

—о— 3th order polynomial

Average mesh step

Рис. 1: Сходимость схем к точному решению на неравномерной сетке

На рисунке 1 приведены графики сходимости численного решения к точному на примере одномерного уравнения переноса (1) с начальными данными в виде гауссиана с шириной 6 = 3. Слева приведены результаты расчётов на время Т = 2, справа — на время Т = 2000. Схемы с 3-, 5- и 7-точечной реконструкцией формально показывают 2-й порядок точности. Однако из рисунка 1 видно, что на относительно грубой сетке при длительном расчёте схема с 3-точечной реконструкцией практически неотличима от полиномиальной схемы 3-го порядка аппроксимации, а схема на 5-точечном шаблоне позволяет добиться ещё большей точности. Дальнейшее же увеличение шаблона на неравномерной сетке в рамках метода разделённых разностей, как правило, бесполезно.

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

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

щ = и(х¿) + сци"(х^

где коэффициенты а* находятся из условия точности вычисления пространственной производной на полиномах 2-го порядка. Это условие представляет собой систему уравнений с 4-диагональной матрицей, сумма элементов в

8

каждом столбце которой (кроме одного) равна 0. Домножение этой системы уравнений на матрицу Т с элементами Т.у = {1,г < j;0,i > j} приводит к 3-диагональной матрице (с дополнительным заполненным первым столбцом) с элементами порядка 0(1), которая при отношении соседних шагов сетки не выше ка-и = 1-148 имеет диагональное преобладание. В негативной норме правая часть исходной системы имеет порядок О^^щ,), а, следовательно, правая часть домноженной системы имеет тот же порядок в норме С. Доказывается, что число обусловленности матрицы домноженной системы ограничено, и поэтому щ = 0(/г^ож). Исходя из этого легко показывается, что в смысле нового оператора проецирования схема имеет 2-й порядок аппроксимации, а, следовательно, в предположении устойчивости схемы, и 2-й порядок сходимости. Если же мерить ошибку в смысле точечного значения, то возникает дополнительная ошибка величиной 0(/г^аа;), что не ухудшает общей оценки.

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

Будем называть сегментом дк общую границу ячеек д и к, если она имеет ненулевую (й — 1)-мерную внутренность. Обозначим множество всех отличных от д ячеек, имеющих общий сегмент с ячейкой д, за N(g).

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

где Уд — объём ячейки д, — ориентированная площадь сегмента дк, внешняя к ячейке д.

Чтобы определить потоки проведём прямую через центры ячеек д и к. В силу равномерности сетки, эта прямая пройдёт через центры других ячеек, и центр каждой следующей ячейки отмерит на этой прямой одинаковое расстояние, см. рисунок 2. Это позволяет определить Р^ как линейную комбинацию значений функции £ в 2М + 1 точках (центрах ячеек) на одномерном шаблоне но формуле (2). Такую схему будем называть схемой с одномерной реконструкцией переменных.

Рис. 2: Схема 5-точечной реконструкции для равномерной двумерной сетки

В параграфе 2.2 доказывается следующее утверждение. Объёмно-центрированные схемы с одномерной реконструкцией переменных на равномерной сетке имеют порядок аппроксимации 0(к2М+1), где М — полуширина ! шаблона реконструкции.

Доказательство проводится в 2 шага. На первом шаге показывается, что для схем с одномерной реконструкцией на равномерной сетке условие произвольного высокого порядка аппроксимации равносильно точности на линейной функции формулы Гаусса-Грина

(4)

9 кеМ(д)

Это следует из того, что граница ячеки равномерной сетки состоит из пар противолежащих сегментов дк и дк, и разделённая разность потоков (Едк — ¥у~к)/(гд — гчерез каждую пару противолежащих сегментов с порядком 2М + 1 аппроксимирует производную по соответствующему направлению (А- — гэ)/|г£ — г9|.

На втором шаге показывается, что на равномерной сетке формула Гаусса-Грина (4) точна на линейной функции. Это очевидно в случае, если вектор (г5 -1- г^/2 совпадает с центром масс сегмента дк и требует дополнительных пояснений в общем случае.

Также исследуется старший член ошибки аппроксимации.

Параграф 2.3 носит вспомогательный характер и содержит описание метода контрольных объёмов для построения консервативных вершинно-центрированных схем на неструктурированной тетраэдральной сетке.

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

объёмов для таких сеток приведены на рисунке 3.

10

Риг. 3: Невыпуклые ячейки, замощающие трёхмерное пространство

Рассматривая в качестве ячеек контрольные объёмы, а в качестве центров ячеек — узлы сетки, вершинно-центрированные схемы можно представить в форме объёмно-центрированных схем. При этом, если исходная сетка — трансляционно-симметрическая и методика построения контрольных объёмов однородна по пространству, то все контрольные объёмы получаются одинаковыми. Отсюда вытекает следующее утверждение. Вершинно-центрированные схемы с одномерной реконструкцией переменных на трансляционно-симметрической сетке имеют порядок аппроксимации 0(/г2Л/+1), где М — полуширина шаблона реконструкции.

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

Рассматривается класс схем, определяемый следующим образом.

1. Метод линий, т. е. счёт по схеме заключается в решении систем ОДУ по времени, источник в которых формируется за счёт пространственной аппроксимации.

2. Консервативность, т. е. схема записывается в потоковом виде:

(все обозначения введены выше).

3. Поток через сегмент на единицу площади сегмента Ру вычисляется по распадной схеме относительно предраспадных значений Р^ и и/или

-¿Г + V • Е(ф) = О

т

4. Предраспадные значения вычисляются как комбинации значений сеточных функций по некоторому набору ячеек:

Fij = atj,к^к

Набор ячеек Т^ будем называть шаблоном реконструкции.

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

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

Параграф 3.2 носит вспомогательный характер. В нём рассматривается известное семество EBR схем, являющихся вершинно-центрированными схемами с квазиодномерной реконструкцией. Будем обозначать их как EBRn, где п — используемое количество узлов реконструкции (3, 5 или 7). Приводится доказательство их точности на линейной функции при использовании барицентрических контрольных объёмов (другое доказательство было приведено в работе Barth, 1991). Использование такого подхода для других, например, ортоцентрических ячеек приводит к потере точности, что сказывается в первую очередь при решении трёхмерных задач. Соответствующие эксперименты будут приведены в главе 4.

В параграфе 3.3 приводится модификация процедуры реконструкции для схем EBR5 и EBR7 на более узком шаблоне, построенная по аналогии с приведённой выше объёмно-центрированной схемой. Для 3-мерной задачи шаблон новой «5-точечной» реконструкции содержит 14 узлов, тогда как шаблон исходной реконструкции формально содержит 26 узлов (некоторые из которых могут совпадать). Это позволяет сэкономить около 50% от времени счёта, затрачиваемое на вычисление предраспадных значений, и таким образом до 25% от общего процессорного времени.

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

1. Определим прямую, проведённую через центр масс ячейки г и центр

масс сегмента ij. Будем говорить, что центр масс сегмента ij находится

12

Рис. 4: Схема двух 5-точечных реконструкций для определения потока на сегменте Ц

«правее» центра масс ячейки г, и обозначать единчный, направленный «вправо» направляющий вектор этой прямой за e¿j.

2. На прямой отметим дополнительно 4 точки, две — справа от центра масс сегмента гj и две — слева от центра масс ячейки г, см. рисунок 4. Будем называть их, вместе с центром ячейки г, узлами реконструкции вдоль прямой 1]. Выбор точек осуществляется таким образом, чтобы значение в них можно было определить линейной интерполяцией по некоторым двум центрам ячеек (отличных от г), имеющих друг с другом общую точку, а на равномерной сетке эта линейная интерполяция становилась бы тривиальной.

3. Предраспадное значение вычисляется как линейная комбинация значений в 5 узлах реконструкции. Учёт неравномерности проводится методом разделённых разностей, описанной в главе 1.

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

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

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

ные. Производные по заданному направлению можно выразить через производные по направлениям реконструкции едк по формуле

при надлежащем определении &дк (см. выше). Производные по направлению едк определяются методом разделённых разностей по 5 или 7 точкам на соответствующей прямой, значения в которых вычисляется линейной интерполяцией. Кратко описываются проблемы, возникающие пути применения данного подхода в практических приложениях.

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

Целью работы, представленной в данной главе, было построение эталонной высокоточной схемы на неструктурированной сетке для оценки точности ЕВ11 схемы. Среди таких схем находят широкое применение конечно-объёмные схемы с полиномиальной реконструкцией и конечно-элементные схемы типа 01зсоп1;тюиз Са1егкш. В сравнении с конечно-разностной ЕВЯ схемой, конечно-объёмные схемы являются значительно более ресурсоёмкими, так как требуют для определения потока через каждый сегмент вычисления значений в каждой точке гауссовой квадратуры на этом сегменте. Упрощение возможно только в линейном случае. К недостаткам конечно-объёмной схемы высокого порядка также можно отнести необходимость проводить осреднение начальных данных и точного решения по ячейке. Примерно те же трудности возникают при использовании метода Галёркина с разрывными базисными функциями. Всё это послужило поводом для разработки сравнительно быстрой конечно-разностной полиномиальной схемы для задач с гладкими решениями, построение которой становится возможным за счёт отказа от требования консервативности.

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

дет приведён результат демонстрационного расчёта, но в линейной постановке.

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

Шаг 1. Дивергенция функции потока выражается через производные по направлениям едь по формуле 5.

Шаг 2. На лучах по направлениям определяются узлы реконструкции так же, как и в схеме с квазиодномерной 7-точечной реконструкцией. Предположим, что значения в них определены и равны ит, т = 0..6. Тогда производная по направлению е^ определяется методом как

где коэффициенты /3 определяются из условия минимизации суммы квадратов при условии точности на многочленах до 5-го порядка:

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

Шаг 3. Будем называть ячейки соседями 1-го порядка, если они имеют общую точку, соседями 2-го порядка, если они имеют общего соседа 1-го порядка и т. д. Значения в узлах реконструкции определяются 3-мерной интерполяцией 5-го порядка. Шаблоном для определения интерполяционного

б

полинома является совокупность центров ячеек — соседей до 3-го порядка включительно от ячейки д. На этот полином также накладывается условие точности на многочленах до 5-го порядка. Из получившегося многопараметрического семейства искомый полином выбирается из условия, чтобы коэффициенты интерполяции отличались от коэффициентов простейшей линейной интерполяции (равных 0 во всех узлах кроме <1, где (I — размерность пространства) минимально в среднеквадратичном смысле. Такой способ определения значения обеспечивает точность на полиноме 5-го порядка и при этом квазиодпомериость. Таким образом, завершено построение обобщения центрально-разностной схемы на неструктурированную сетку.

Шаг 4. Добавление диссипации. На каждом сегменте дк с двух сторон определяются интерполированные значения по полиномам 4-го порядка, определённым по соседям 2-го порядка от ячеек дик соответственно. Полиномы строятся методом наименьших квадратов. Чтобы сохранить квазиодномерность, в ячейках, входящих в шаблон реконструкции ЕВИ5 для соответствующего предраспадного значения, условие накладывается строго. Диссипативный поток через сегмент на единицу площади определяется по формуле

где 5 — коэффициент диссипации. Хотя формально диссипативные потоки вносят ошибку порядка 0(Л4), благодаря консервативному определению дис-сипативных членов схема показывает 5-й порядок точности.

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

В параграфе 4.2 описываются детали реализации и изучается ресурсо-ёмкость предложенной схемы. В 3-мерном случае время счёта по полиномиальной схеме при явном интегрировании по времени примерно на 75% больше, чем по известной схеме ЕВИ и примерно в 2.5 раза больше, чем по предлагаемой в параграфе 3.3 модификации схемы ЕВГ1. Основной вклад (около 75%) в затраты процессорного времени вносит вычисление диссипативных потоков. Расход памяти составляет на хранение коэффициентов составляет примерно 15 КБ/узел, из них около 6 КБ/узел — коэффициенты центральной части и около 9 КБ/узел — коэффициенты диссипативной части.

В параграфе 4.3 предложенная полиномиальная схема сравнивается со схемами с квазиодномерной реконструкцией для 3-мерных линеаризованных уравнений Эйлера на неструктурированной сетке. Данный параграф дополняет экспериментальное исследование точности ЕВИ. схем, приведённое в параграфе 3.4.

Рис. 5: Коэффициент прохождения волны для трёхмерной неструктурированной сетки Таблица 1: Ошибка относительно точного решения

Сетка EBR3 EBR5 EBR7 Полиномиальная 5-го порядка

Декартова, h—1 6.59 ■ Ю-3 3.09- 10"4 4.30 • 10~5 2.44 ■ 10'4

Неструкт., h=l 4.79 • 10"3 1.86-Ю-3 1.77- 10"3 9.60 • Ю-5

Неструкт., h=l/2 6.53 • 10~4 5.63 ■ Ю-4 5.79 ■ Ю-4 2.75 • Ю-6

На рисунке 5 приведен коэффициент прохождения плоской акустической волны через область трёхмерной неструктурированной сетки длиной в 5000 характерных шагов. Слева приведено сравнение схем EBR3, EBR5 на барицентрических ячейках и неконсервативной полиномиальной схемы 5-го порядка. Справа добавлены графики схем EBR3 и EBR5 на ортоцентриче-ских ячейках. Видно, что отсутствие точности схемы на линейной функции (при использовании ортоцентрических контрольных объёмов) в случае приводит к очень сильной диссипации схемы на низких частотах.

В таблице 1 приведены результаты расчётов распада 3-мерного гауссо-вого импульса с полушириной Ъ = 6 на момент времени Т = 20. Использовалась неструктурированная тетраэдральная сетка, коэффициент диссипации S = 0.5. Ошибка мерилась в норме С. Видно, что, как и в одномерном случае, схема EBR5 существенно лучше EBR3, но уступает неконсервативной полиномиальной схеме на шаблоне той же ширины.

В параграфе 4.4 предложенная квазиодномерная полиномиальная схема применяется для исследования линейного резонанса резонатора Гельм-гольца, изображённого на рисунке 6. Результаты расчёта по полиномиальной схеме сравнивались с расчётами по схеме EBR5, проведёнными А. Дубенем.

Пятая глава посвящена численному моделированию акустического рупора. Целью численного исследования является оценки коэффициента усиления сигнала акустическим рупором в зависимости от угла раскрыва рупора

Рис. 6: Резонатор Гельмгольца

и угла падения волны относительно оси рупора. Данная задача была поставлена группой мобильных роботов ИПМ им. М. В. Келдыша РАН во главе с д.ф.-м.н. В. Е. Павловским.

М11ж»жЖДиМ1вД

-10 -5 0 5 10

X

Рис. 7: Двумерный рупор

Рупор представляет собой конус с длиной образующей в 10 см и толщиной 3 мм. 2-мерный рупор с углом 2(3 = 7г/2 изображён на рисунке 7, а 3-мерный рупор с тем же плоским углом — на рисунке 8. В качестве математической модели использовалось волновое уравнение для пульсации давления. На стенках рупора ставилось условие непротекания (равенство нулю нормальной производной от пульсации давления), что соответствует предположению абсолютной жёсткости стенок.

В параграфе 5.1 описывается программный код ГЮШЕйе для решения задач газовой динамики и аэроакустики [3], в котором реализованы предлагаемые в настоящей работе численные методики.

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

Рис. 8: Рупор

дифференциальных уравнений 1-го порядка

др „ <3и __

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

В параграфе 5.3 приведены расчёты в двумерной постановке. Рупор моделировался плоским углом.

Численный расчёт: —■— 5° —*— 10°—•— 15°—45°

Аналитика узкого угла:....... 0°

Высокочастотная асимптотика:—■—45°

о 2000 Frequency [Hz] 4000 6000

Рис. 9: Коэффициент усиления волны в двумерном рупоре

На рисунке 9 приведен коэффициент усиления сигнала в двумерном

рупоре при различных углах раскрыва рупора в зависимости от частоты,

в предположении, что рупор направлен на источник звука. Для рупора с

углом 2/3 = 90° проведено сравнение с аналитическим решением, полученным

19

+50dB -,

+40dB-

+OdB

2000

4000

6000

Частота [НЩ

Рис. 10: Коэффициент усиления волны в аксиально-симметричной постановке 2СМВ

15dB

с 10dB о

го

о

'р- 5dB

-5dB -

Расчёт 8 аксиально симметричной геометрии: —— 3=5° —3=15° 3=30° 3=45°

------Аналитическое решение для 3=0°

о Высокочастотная асиптотика для 3=45°

-'-1-'-1-1-1-•-1-1-1

0 500 1000 1500 2000 2500

Frequency [Hz]

Рис. 11: Коэффициент усиления волны в 3-мерной постановке

путём пренебрежения многократным рассеянием. Также приведено сравнение с предельным аналитическим решением при стремлении угла /? —> 0.

В параграфе 5.4 приведены расчёты в аксиально-симметричной постановке для волны, падающей по направлению оси рупора. Были рассмотрены случаи различных углов раскрыва рупора. Также проведено сравнение предельным аналитическим решением при стремлении угла раскрыва рупора к 0. Результаты приведены на рисунке 10.

Исходя из расчётов в 2-мерной и аксиально-симметричной постановках, видно, что качественное поведение кривой коэффициента усиления объ-

ясняется исходя из двух асимптотических решений — асимптотического решения при безконечно узком рупоре (/? = 0) и высокочастотной асимптотики, получаемой путём пренебрежения многократным перерассеянием волн на кромках рупора.

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

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

Основные результаты диссертационной работы.

1. Доказано, что вершинно-центрированные схемы с одномерной реконструкцией переменных имеют произвольно высокий порядок аппроксимации на трансляционно-симметрической сетке, а объёмно-центрированные схемы — на равномерной сетке.

2. Построена объёмно-центрированная схема с квазиодномерной реконструкцией переменых.

3. Построена экономичная формулировка ЕВ11-схемы, позволяющая сократить время счёта до 25%.

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

Благодарности. Автор выражает благодарность:

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

• д.ф.-м.н. Полякову С. В. за безвозмездно предоставленные вычислительные ресурсы;

• библиотеке численного анализа НИВЦ МГУ за предоставленные процедуры вычисления функций Бссссля.

Публикации автора по теме диссертации

1. Абалакин И. В., Бахвалов П. А., Козубская Т. К. Схема на основе реберно-ориентированной реконструкции переменных и ее \¥ЕМО-версия для решения задач аэродинамики и аэроакустики на неструктурированных тетраэдральных сетках // Сборник тезисов Четвертой всероссийской конференции "Вычислительный эксперимент в аэроакустике". — М. : МАКС Пресс, 2012. - С. 5-6.

2. Бахвалов П. А., Козубская Т. К., Корнилина Е. Д. и др. Технология расчетов акустических пульсаций в дальнем поле течения // Математическое моделирование. — 2011. — Т. 11. — С. 33-47.

3. Абалакин И. В., Бахвалов П. А., Горобец А. В. и др. Параллельный программный комплекс ГЮ^ЕИе для крупномасштабных расчетов задач аэродинамики и аэроакустики // Вычислительные методы и программирование. - 2012. — Т. 13. - С. 110-125.

4. Бахвалов П. А. Схема с квазиодномерной реконструкцией переменных на сетках из выпуклых многоугольников для решения задач аэроакустики // Математическое моделирование. — 2013. — Т. 9. — С. 95-108.

5. Павловский В. Е., Бахвалов П. А., Козубская Т. К., Павловский В. В. Слуховой сенсор робота: численная оптимизация и экспериментальное исследование // Тр. Межд. научно-технической конф. "Мехатроника, автоматизация, управление - 2009" (МАУ-2009).— Изд-во ЮФУ, 2009.-С. 261-264.

6. Бахвалов П. А., Козубская Т. К., Павловский В. Е. Численное моделирование приема акустических сигналов в робототехнике //В Тезисах докладов всероссийской открытой конференции по авиационной акустике (октябрь 2009 г.). — М. : Центральный Аэрогидродинамический институт имени проф. Н.Е.Жуковского (ЦАГИ), 2009.- С. 66-67.

7. Бахвалов П. А., Козубская Т. К., Корнилина Е. Д. и др. Технология расчетов акустических пульсаций в дальнем поле течения //В сборнике тезисов Третьей открытой всероссийской конференции "Вычислительный эксперимент в аэроакустике". — М. : МАКС Пресс, 2010.— С. 19-20.

8. Козубская Т. К., Павловский В. Е., Бахвалов П. А. Численное моделирование слухового сенсора робота //В Трудах 53-ей научной конференции МФТИ "Современные проблемы фундаментальных и прикладных наук",

Часть VII, Управление и прикладная математика, Том 3. — М. : МФТИ, 2010.-С. 100-101.

9. Абалакин И. В., Бахвалов П. А., Горобец А. В. и др. Программный комплекс NOlSEtte для расчетов задач газовой динамики и аэроакустики на неструктурированных сетках //В Тезисах докладов Второй всероссийской открытой конференции по авиационной акустике (сентябрь 2011г.).— М. : Центральный Аэрогидродинамический институт имени проф. Н.Е.Жуковского (ЦАГИ), 2011,- С. 100-101.

10. Düben А. P., Abalakin I. V., Bakhvalov Р. А. et al. Exploiting modern supercomputers in a research towards quieter aircrafts // Parallel CFD 2011 Books of abstracts. — Barcelona, Spain, 2011. — P. 30.

11. Abalakin I., Bakhvalov P., Kozubskaya T. EBR schemes: new developments // Book of Abstracts of European Workshop on High Order Nonlinear Numerical Methods for Evolutionary PDEs: Theory and Applications (HONOM 2013), March 2013, Talence, Bordeaux, France.— 2013,— http://honom2013. bordeaux.inria.fr/Bons/Kozubskaya.pdf.

12. Бахвалов П. А. Объемно-центрированная схема с квазиодномерной реконструкцией переменных для решения задач газовой динамики на сетках из выпуклых многоугольников // Сборник тезисов Четвертой всероссийской конференции "Вычислительный эксперимент в аэроакустике". — М. : МАКС Пресс, 2012.-С. 27.

13. Абалакин И. В., Бобков В. Г., ЖдановаН. С. и др. Численное моделирование аэродинамикии и аэроакустики фенестрона в модельной конфигурации //В Тезисах докладов Третьей всероссийской открытой конференции по авиационной акустике (сентябрь 2013г.). — М. : Центральный Аэрогидродинамический институт имени проф. Н.Е.Жуковского (ЦАГИ), 2013. — С. 211-212.

Бахвалов Павел Алексеевич

РАЗВИТИЕ СХЕМ НА ОСНОВЕ КВАЗИОДНОМЕРНОГО ПОДХОДА ДЛЯ РЕШЕНИЯ ЗАДАЧ АЭРОАКУСТИКИ НА НЕСТРУКТУРИРОВАННЫХ СЕТКАХ

АВТОРЕФЕРАТ

Подписано в печать 8.11.2013г. Формат 60х841/1б. Усл. печ. л. 1,25. Тираж 100 экз. Заказ №374. Федеральное государственное автономное образовательное учреждение высшего профессионального образования «Московский физико-технический институт (государственный университет)» Отдел оперативной полиграфии «Физтех-полиграф» 141700, Московская обл., г. Долгопрудный, Институтский пер., 9.

МОСКОВСКИЙ ФИЗИКО-ТЕХНИЧЕСКИЙ ИНСТИТУТ

(ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ) КАФЕДРА МАТЕМАТИЧЕСКОГО МОДЕЛИРОВАНИЯ

На правах рукописи УДК 519.63

04201451770

БАХВАЛОВ ПАВЕЛ АЛЕКСЕЕВИЧ

РАЗВИТИЕ СХЕМ НА ОСНОВЕ КВАЗИОДНОМЕРНОГО ПОДХОДА ДЛЯ РЕШЕНИЯ ЗАДАЧ АЭР О АКУСТИКИ НА НЕСТРУКТУРИРОВАННЫХ СЕТКАХ

Специальность 05.13.18 — «математическое моделирование, численные методы и комплексы программ»

Диссертация на соискание учёной степени кандидата физико-математических наук

Научный руководитель: доктор физико-математических наук, старший научный сотрудник Козубская Т. К.

Москва - 2013

Содержание

Введение ................................................................................5

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

1.1 Схемы для одномерного уравнения переноса на равномерной сетке............11

1.2 Случай нелинейного уравнения и нелинейной системы уравнений............14

1.2.1 Различные типы реконструкции............................................14

1.2.2 Анализ ошибки аппроксимации в конечно-разностном смысле .... 15

1.2.3 Анализ ошибки аппроксимации в конечно-объёмном смысле ..........18

1.2.4 Гибридный подход к реконструкции ......................................21

1.2.5 Объёмно-центрированная схема на прямоугольниках....................22

1.3 Схемы для одномерного уравнения переноса на неравномерной сетке .... 23

1.3.1 Метод перехода к разделённым разностям................................24

1.3.2 Экспериментальное исследование точности построенных схем..........26

1.4 Аналитическое исследование сходимости..........................................27

1.4.1 История вопроса..............................................................27

1.4.2 Доказательство 2-го порядка точности для схемы с 3-точечной реконструкцией ......................................................................29

2 Многомерные схемы с одномерной реконструкцией на равномерной сетке ..................................................................................35

2.1 Равномерная сетка....................................................................35

2.2 Точность схемы с одномерной реконструкцией переменных....................37

2.2.1 Схемы с одномерной реконструкцией переменных........................37

2.2.2 Доказательство произвольного высокого порядка аппроксимации . . 38

2.2.3 Проблема анизотропии на одинаковых многогранниках................40

2.3 Метод контрольных объёмов для построения консервативных вершинно-центрированных схем................................................................42

2.3.1 Двумерный случай ..........................................................42

2.3.2 Трёхмерный случай..........................................................43

2.3.3 Различные виды контрольных объёмов....................................45

2.3.4 Огрублённые контрольные объёмы........................................48

2.3.5 Точность на линейной функции для барицентрических ячеек..........49

2.4 Вершинно-центрированные схемы с одномерной реконструкцией..............51

3 Схемы с квазиодномерной реконструкцией переменных на неструктурированной сетке ..................................................................53

3.1 Определение схемы с квазиодномерной реконструкцией........................53

3.2 Вершинно-центрированные схемы ЕВИ,............................................54

3.3 Экономичная формулировка схемы ЕВИ...........................................59

3.4 Объёмно-центрированная схема с квазиодномерной реконструкцией переменных ................................................................................60

3.5 Сравнение вершинно- и объёмно-центрированных схем..........................02

3.6 Обобщение схемы на уравнения, содержащие вторые производные............64

4 Конечно-разностная полиномиальная схема..................................68

4.1 Построение схемы....................................................................69

4.2 Ресурсоёмкость схемы................................................................70

4.3 Экспериментальное сравнение схем ЕВ11 и полиномиальных схем..............70

4.4 Численное исследование линейного резонанса резонатора Гельмгольца ... 71

5 Численное моделирование акустического рупора ..........................73

5.1 Программный код ГЮГЗЕ^е для решения задач аэродинамики и аэроакустики 73

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

5.3 Решение в двумерной постановке ..................................................78

5.3.1 Аналитическое решение для бесконечно узкого рупора..................78

5.3.2 Двухволновое приближение (высокочастотная асимптитока)..........80

5.3.3 Исследование коэффициента усиления рупора............................83

5.4 Решение в аксиально-симметрической постановке................................84

5.4.1 Анализ сигнала во временном представлении для тестовой постановки 84

5.4.2 Исследование коэффициента усиления рупора............................88

5.5 Решение в трёхмерной постановке..................................................88

Заключение ............................................................................90

Список рисунков ......................................................................91

Список таблиц.........................................................93

Список использованных источников ..............................................94

Введение

Актуальность и мотивация работы.

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

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

Использование перекрывающихся многоблочных сеток было впервые предложено в 1977 году в работе Gören Starius [1] для решения задачи Дирихле для оператора Лапласа в криволинейной области. Рассматривались две сетки, причём искусственная граница каждой из сеток была погружена в другую сетку. Алгоритм решения задачи представлял собой итерационный процесс, на каждой итерации которого выполнялась интерполяция в узлы, лежащие на искусственной границе расчётных сеток, после чего на каждой из сеток выполнялась задача с полученными искусственными граничными условиями.

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

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

Чуть позже перекрывающиеся многоблочные сетки были применены к гиперболической системе уравнений на примере уравнений мелкой воды [2]. Далее, в работах [3], [4] они были применены для решения аэродинамических задач в рамках стационарных уравнений Эйлера, в работе [5] — в рамках стационарных уравнений Навье-Стокса и в работе [6] — для нестационарых задач с использованием неявной схемы интегрирования по времени. Программная реализация метода перекрывающихся сеток подробно описана в работе [7].

Дальнейшее развитие многоблочного подхода заключалось в обобщении его на движущиеся сетки (см., например, [8]), применении более эффективных технологий для решения систем алгебраических уравнений, возникающих при неявном интегрировании по времени (см., например, [9], [10]). Отметим также работу [11], посвящённую вопросу переинтерполяции между сетками на сшивке блоков. Значительная доля исследований в этой области посвящена автоматизации процесса построения многоблочных сеток и определения множества точек, участвующих в процессе переинтерполяции данными между блоками ( [12], [13] и др.) Несмотря на значительный прогресс в этом направлении, проблема построения мпогоблочной структурированной сетки в трёхмерном случае остаётся до конца не решённой.

Более существенной трудностью использования многоблочных сеток для решения широкого класса аэродинамических задач является проблема консервативности на сшивке различных блоков. Решение этой проблемы впервые было предложено в работе [14], однако, следуя [15], можно сказать, что вычисления по этой схеме никогда не были опубликованы. В работе [15] для обеспечения консервативности предлагается отказаться от наложения сеток и использовать стыкующиеся (patched) сетки, получаемые из перекрывающихся сеток на этапе предварительной обработки. Однако без дополнительной геометрической коррекции такой подход может приводить к появлению ячеек сколь угодно малого размера, и поэтому его универсальность до конца не ясна. Во многих работах (см., например, [16]) для задач с разрывными решениями продолжают использоваться неконсервативные схемы. Таким образом, задача построения оператора интерполяции для сжимаемых уравнений, обеспечивающего консервативность, остаётся до конца не решённой, см. также обзор в [17]. Это существенно ограничивает применение схем для решения сжимаемых уравнений Эйлера и Навье-Стокса при наличии разрывов и высоких градиентов решения на перекрывающихся многоблочных сетках, в особенности

при необходимости моделирования течения с высокой точностью.

Альтернативой схемам на многоблочных сетках являются схемы на неструктурированных сетках. Среди схем высокого порядка точности на неструктурированных сетках, пригодных для решения аэроакустических задач, активно развиваются два главных направления: конечно-элементные схемы, главным образом — метод Галёркина с разрывными базисными функциями (DG), и конечно-объёмные схемы, основанные на полиномиальной реконструкции неременных.

В отличие от схем на перекрывающихся многоблочных сетках, разработка методов счёта на неструктурированных сетках изначально велась для задач с разрывными решениями. Повышение порядка аппроксимации за счёт линейной реконструкции переменных, предложенное в одномерном случае в работе В. П. Колганом [18] и независимо от него Ван Лиром [19] (см. также [20], [21], [22]), было применено к неструктурированной треугольной сетке в работе F. Angrand и A. Dervieux [23] и к тетраэдральной сетке в [24]. Однако для задач, требующих высокой точности, использование лишь линейной реконструкции переменных было явно недостаточным.

Метод полиномиальной реконструкции переменных на неструктурированной сетке, точной на полиномах произвольно высокого порядка, был впервые описан в работе Barth [25], будучи на тот момент уже хорошо известным. В работе Abgrall [26] была построена схема методом ENO [27] [28], совмещающая высокую точность полиномиальной реконструкции с возможностью счёта разрывных задач. Первые применения к неструктурированным сеткам взвешивания полиномов по аналогии с WENO были опубликованы в работах [29] и [30], однако полноценная WENO-схема на неструктурированных сетках была предложена в работе Ни и Shu [31]. Для эффективного применения WENO-схем на неструктурированных сетках понадобились дополнительные исследования конечно-объёмных схем на многомерной декартовой сетке [32] [33] [34]. Дальнейшее развитие и применение ENO и WENO схем на неструктурированных сетках связано с работами [35], [36], [37] и др.

Метод Галёркина с разрывными базисными функциями (Discontinious Galerkin method, DG) был впервые применён в 1973 году для уравнения переноса нейтронов [38]. В отличие от конечно-объёмных схем, в которых для повышения порядка точности используются значения функции на широком шаблоне, конечно-элементные схемы подразумевают повышение порядка точности за счёт задания нескольких значений на одной ячейке. Первые применения метода DG к задачам газовой динамики связаны с работами В. Cockburn и

C.-W. Shu [39], [40], [41], [42]. За последние 20 лет он претерпел существенное развитие как минимум в трёх направлениях: выбор оптимального набора базисных функций [43], уточнение оценок точности и исследование супер- и сверхсходимости [44] [45] [46] [47] [48] и построение лимитеров [49] [50]. Отметим также работы [51] и [52], в которых используется анализатор разрыва, и при обнаружении последего производится переключение с DG на конечно-объёмную WENO-схему.

Схемы, построенные методом Галёркина с разрывными базисными функциями, и конечно-объёмные схемы с полиномиальной реконструкцией переменных близки по своим свойствам [53]. Они характеризуются возможностью построения схем, точных на полиномах сколь угодно высокого порядка, и низкой чувствительностью к качеству сетки. Недостатком обоих методов является их большая вычислительная стоимость для нелинейных задач, особенно при необходимости воспроизведения разрывов.

История работ по теме диссертации.

В 1998 году A. Dervieux и С. Debiez предложили группу разностных схем под названием «Mixed-element-volume MUSCL methods with weak viscosity» (смешанные конечно-элементные-конечно-объёмные низкодиссипативные методы типа MUSCL) [54] [55]. Их идея заключалась в построении консервативной схемы 2-го порядка на произвольной неструктурированной сетке, которая в случае декартовой сетки вырождалась бы в схему 5-го порядка на основе 6-точечных разностей для дискретизации пространственной производной. В той же работе были предложены TVD-модификации данной схемы. Методика была применена к 2-мерному уравнению переноса. Детальный анализ одномерной схемы для линейных задач был проведён в работе С. Debiez [56]. К системе линеаризованных уравнений Эйлера эта методика была применена в работе [57] для 2-мерных задач и в работе [58] для 3-мерных.

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

Обобщение данной методики на нелинейные задачи было предложено в работе А. Dervieux, И. А. Абалакина и Т. К. Козубской в 2006 году в работе «High Accuracy Finite Volume Method for Solving Nonlinear Aeroacoustics Problems on Unstructured Meshes» (высокоточный конечно-объёмный метод для решения нелинейных аэроакустических задач на неструктурированных сетках) [60]. Как и для линейных задач, на декартовой сетке разра-

батываемая схема должна была вырождаться в схему 5-го порядка. Поскольку известно в рамках конечно-объёмного подхода (построения схем относительно средних значений по ячейке) таких методов не существует [61] (см. также главу 1 диссертационной работы), то было предложено строить конечно-разностную схему (т. е. записываемую относительно точечных значений). Таким образом, можно сказать, что эта работа открыла новый класс разностных схем, а именно, консервативных конечно-разностных схем повышенной точности на неструктурированных сетках. Эта схема также была опубликована в книгах [62].

Среди применений этого класса схем для расчёта промышленных задач (и тестовых задач, позволяющих показать применимость схем для промышленных задач) можно отметить зарубежные [63] и российские работы [64] [65] [66] [67] [68] [69] [70] [71] [72] [73] [74] [75].

Хотя в ранее упомянутых работах [59], [54] и др. монотонизация рассматриваемого класса схем проводилась путём построения TVD схем, для этой цели также возможно использование WENO подхода, см. [76] [77]. В этих работах рассматриваемые схемы были названы схемами с реконструкцией переменных вдоль направления ребра (Edge-Based Reconstruction, EBR).

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

Краткое содержание диссертации.

Диссертация состоит из введения, пяти глав, заключения и списка литературы.

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

Основное внимание уделяется методу перехода к разделённым разностям для построения схем на неравномерной сетке, позволяющий сохранить консервативность, обеспечить точность на линейной фу�