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

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

Развитие схем

на основе квазиодномерного подхода

для решения задач аэроакустики

на неструктурированных сетках

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

численные методы и комплексы программ

Автореферат

диссертации на соискание учёной степени

кандидата физико-математических наук

Москва 2013

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

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

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

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

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

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

Автореферат разослан 2013 года.

Ученый секретарь диссертационного совета Д 212.156.05 Федько О.С.

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

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

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

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

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

В 1998 году A. Dervieux и C. 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, B. Koobus, F. Alauzet, И.

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

Koobus, M.-V. Salvetti, S.Camarri, И. Абалакина, А. Горобца, А. Дубеня, Т.

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

• семинар Математическое моделирование ИПМ им. М. В. Келдыша РАН (2013), а также на следующих конференциях:

• Международная научно-техническая конференция "Мехатроника, автоматизация, управление - 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).

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

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

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

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

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

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

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

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

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

Доказательство данного утверждения основано на введении вспомогательного оператора проецирования функции на сетку: значение сеточной функции определяется как где коэффициенты i находятся из условия точности вычисления пространственной производной на полиномах 2-го порядка. Это условие представляет собой систему уравнений с 4-диагональной матрицей, сумма элементов в каждом столбце которой (кроме одного) равна 0. Домножение этой системы уравнений на матрицу T с элементами Tij = {1, i j; 0, i > j} приводит к 3-диагональной матрице (с дополнительным заполненным первым столбцом) с элементами порядка O(1), которая при отношении соседних шагов сетки не выше kcrit = 1.148 имеет диагональное преобладание. В негативной норме правая часть исходной системы имеет порядок O(h2 ), а, следоваmax тельно, правая часть домноженной системы имеет тот же порядок в норме C. Доказывается, что число обусловленности матрицы домноженной системы ограничено, и поэтому i = O(h2 ). Исходя из этого легко показывается, что в смысле нового оператора проецирования схема имеет 2-й порядок аппроксимации, а, следовательно, в предположении устойчивости схемы, и 2-й порядок сходимости. Если же мерить ошибку в смысле точечного значения, то возникает дополнительная ошибка величиной O(h2 ), что не ухудшает общей оценки.

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

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

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

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

Рис. 2: Схема 5-точечной реконструкции для равномерной двумерной сетки В параграфе 2.2 доказывается следующее утверждение. Объёмноцентрированные схемы с одномерной реконструкцией переменных на равномерной сетке имеют порядок аппроксимации O(h2M +1 ), где M полуширина шаблона реконструкции.

Доказательство проводится в 2 шага. На первом шаге показывается, что для схем с одномерной реконструкцией на равномерной сетке условие произвольного высокого порядка аппроксимации равносильно точности на линейной функции формулы Гаусса-Грина Это следует из того, что граница ячеки равномерной сетки состоит из пар противолежащих сегментов gk и g k, и разделённая разность потоков (Fgk Fgk )/(rg rk ) через каждую пару противолежащих сегментов с порядком 2M +1 аппроксимирует производную по соответствующему направлению (rk rg )/|rk rg |.

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

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

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

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

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

Вершинно-центрированные схемы с одномерной реконструкцией переменных на трансляционно-симметрической сетке имеют порядок аппроксимации O(h2M +1 ), где M полуширина шаблона реконструкции.

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

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

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

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

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

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

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

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

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

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

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

1. Определим прямую, проведённую через центр масс ячейки i и центр масс сегмента ij. Будем говорить, что центр масс сегмента ij находится Рис. 4: Схема двух 5-точечных реконструкций для определения потока на сегменте ij правее центра масс ячейки i, и обозначать единчный, направленный вправо направляющий вектор этой прямой за eij.

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