О построении wеnо схем для гиперболических систем уравнений на неструктурированных сетках

УДК 519.6

В. А. Горюнов, Р. В. Жалнин, Е. Е. Пескова, В. Ф. Тишкин

О ПОСТРОЕНИИ WENO СХЕМ ДЛЯ ГИПЕРБОЛИЧЕСКИХ СИСТЕМ УРАВНЕНИЙ НА НЕСТРУКТУРИРОВАННЫХ СЕТКАХ

Аннотация.

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

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

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

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

Ключевые слова: WENO схема, неструктурированная сетка, высокий порядок точности.

V. A. Goryunov, R V. Zhalnin, E. E. Peskova, V. F. Tishkin

ON CONSTRUCTION OF WENO SCHEMES FOR HYPERBOLIC SYSTEMS ON UNSTRUCTURED MESHES

Abstract.

Background. Mathematical simulation of fluid and gas flows is reduced to solving the equations of the Euler system in areas of complicated geometry. Real flows are characterized by the appearance of gas dynamic discontinuities. It leads to the usage of numerical methods of high order accuracy. The purpose of this paper is to construct essentially non oscillatory high order scheme (WENO scheme) on unstructured meshes for gas dynamics equations; and to compare the obtained results with the numerical results using first-order accuracy.

Materials and methods. The basic idea of WENO scheme is linear combination of polynomials constructed by the ENO scheme. Weights in the linear combination depend on the smoothness of solution in each set. For treatment of negative weights the authors carried out disintegration thereof.

Results. The authors developed an essentially non oscillatory third-order scheme (WENO scheme) on unstructured meshes for gas dynamics equations. The obtained results were compared with the numerical results using first-order accuracy.

Conclusions. The researchers developed a third order scheme using a combination of linear polynomials. A series of test calculations for the Riemann problem using this scheme was performed. It is concluded that the proposed scheme smears the solutions on discontinuities less than the first order scheme.

Key words: WENO scheme, unstructured mesh, high-order accuracy.

Введение

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

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

1. Математическая модель и разностная схема

Рассмотрим систему уравнений газовой динамики в переменных Эйлера:

эи Э^

dt

- + -

дх

- + -

дУ

= 0,

где

U =

Р

v v

Р

pu pv

( 2 +2 ^ и + v £ +--------

/у

f\u ) =

pu

2+

puv

2

pu + p

pu

£ + -

+ pu

; f2(U) =

pv

£ + -

pv

puv

pv2 + p

2 , 2 ^ u + v

v

+ pv

Система дополнена уравнением состояния

р =

Здесь р,р, є - плотность, давление и удельная внутренняя энергия соответственно; и, V - компоненты вектора скорости; у - показатель адиабаты.

Введем треугольную сетку ю = {Р;- = (х^, уі ), і = 1,2,..., N} , содержащую все внутренние и граничные точки расчетной области. На ю построим триангуляцию Л(ю) ={лот =Д(ри,Р]т,Ркт),Ріт,Р]т,Ркт є ю т = 1,2,...,М} .

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

И 1 (•

—и

Иі Д.-1 з

'ЭА

Здесь ^ ^(1), ^ (2))Т , п - внешняя нормаль к грани ячейки ЭЛг-.

Интеграл в (3) рассчитываем по квадратурной формуле Гаусса:

ЭЛг- к=1 у=1

где ^(и + (о(к),t), и (о(к),t))• Пк - дискретные потоки, которые рассчитываются по схеме распада разрыва [1]; и + (о(к\t), и-(о(к),t) - «левые» и

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

Используем двухточечную квадратуру Гаусса д = 2. Для ребра

треугольника с координатами /1 и р точки Гаусса определяем следующим

1 73 1

образом: 01 = с/1 + (1 - с)/, О2 = с/2 + (1 - с)р, где с = — +-, Ш = ш = _.

2 6 2

2. Линейная реконструкция

Пусть дана триангуляция области {Л1, Л 2,..., Л N} и заданы средние значения некоторой функции и (х, у) для каждой ячейки сетки

Л (/ = 1,2,...,К):

— 1 Г

иI = П""Г I и

N Л

Для каждой ячейки Лг- (I = 1,2,...,N) построим полином р(х,у) степени, не больше к, который интерполирует функцию и(х,у) с порядком точности к + 1 такой, что

1 Г —

I «'I лг

Для построения полиномов Р1 (х, у) к -й степени для ячейки Л| будем

гг О < Л Л 1

рассматривать шаблоны = {—1,...,Лк}, где К =-- - количество

ячеек в каждом из шаблонов. Используя выражение (6) для каждого треугольника из шаблона и решив линейную систему К X К уравнений, найдем искомый полином р^ (х, у).

Для построения линейной схемы третьего порядка точности [2] для нахождения газодинамических параметров для ячейки —о используем шаблон, включающий в себя два ряда соседних ячеек ^ = {—о, Лг-, Л|а, Лгь, Л j, Л ja,

—Ь, Лк, Лка, ЛкЬ } (рис. 1).

Рис. 1

Используя данный шаблон, строим полином второй степени р2 (х,у,г) = ао + а1х + а2у + азх2 + а4у2 + а5ху методом наименьших квадратов. Для каждой точки Гаусса ,у0) находим набор коэффициентов

{с/ }/=!, которые зависят только от геометрии ячеек:

N

Р 2

N

l=l

Здесь N - число треугольников в шаблоне; и - среднее значение и(х,у) в ячейке.

Идея WENO схемы заключается в следующем [2, 3]. Строим линейные полиномы р, (х, у), взвешенная сумма которых дает тот же результат, что и

квадратичный полином р2 (х, у).

Используя шаблон S, строим линейные полиномы

Рі (х,у) = а) + о!х + а2у, і = 1,..., 9 , на следующих шаблонах-кандидатах [2]: 51 = (А),Аj,А^}, 52 = (Ао,А,,Ак},

£3 = (Ао, Аі , А j}, £4 = (Ао, Аі , Аіа}, £5 = (А), А,-, АіЬ}, 56 = (А), А}-, А^},

£7 = (А0, А} , А}Ъ Ь £8 = (А0, Ак, Ака }, £9 = (^ Ак, АкЪ }.

Используя выражение (6) для каждого треугольника из шаблона и решив линейную систему 3 X 3, найдем искомый полином р, (х, у). Например,

для шаблона-кандидата 51 = (А), А j, Ак } получим следующую систему:

1 Г —

|—г J Pl(x, y)dxdy = uo,

0 A A0

1 r _

I—I I Pl(x, y)dxdy = u j, A

Iа jA;

1 r _

7—7 I Pl(x, y)dxdy = uk.

M—

Лк

Для каждой точки Гаусса (х0,у0) находим набор коэффициентов

{ (О}3

|с/ | , которые зависят только от геометрии ячеек:

р

/=1

Здесь и/ - среднее значение и в ячейке.

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

9

К

/'=1

который удовлетворяет условию

к

Из равенства (10) получаем линейную систему уравнений вида

М у = с ,

Т

где с = (с1, с2,..., сN) - вектор, координатами которого являются коэффици-

енты в (7) для полинома второй степени. Каждый столбец матрицы М состоит из коэффициентов в (8) для полиномов первой степени.

3. WENO реконструкция

Алгоритм WENO построения интерполяционного полинома заключается в использовании комбинации всех возможных полиномов для данной ячейки [2]:

N

Рweno

т=1

где N - количество полиномов; wш - нелинейные весовые коэффициенты, сохраняющие все свойства существенно не осциллирующих схем:

N ~ У

2 w = 1 w =-------ш— — =--------ш-----

/ . ш ’ ш N ’ ш 2 ' '

ш=1 X-' -г.

/ , »г

ш=1

Здесь е - малая положительная величина, введенная, чтобы избежать

_з

деления на ноль. В расчетах принимаем е = 10 ; Шш - индикатор гладкости

для полинома Рш (х, у):

Шш = 2 11—011

фи—0

где к - степень полинома Рш (х, у).

Линейные весовые коэффициенты {уг- }^=1 зависят только от геометрии ячейки, поэтому могут быть отрицательными. Если шт(у1, У2,.., У N) < 0, то применяем технологию расщепления весовых коэффициентов [4]. Разбиваем линейные коэффициенты на две группы:

Г+= ^(Уг + 3 |ъ 1^ Т_=Т+_Уг,

N

а±= 2, у±=у± /, /=l,...,N.

/=1

Используя линейные коэффициенты у±, вычисляем нелинейные весовые коэффициенты . Затем проводим WENO реконструкцию Р±епо (х°,у0) для каждой группы нелинейных весовых коэффициентов. В результате WENO реконструкция представляется следующим выражением:

Pweno Iх ,y ) ) Pweno (x ,y ) ) Pweno (x ,y

4. Тестовые расчеты

Была проведена серия тестовых расчетов для задачи Римана (задача о распаде разрыва). Расчеты были выполнены в двумерной постановке для начальных данных Сода [5]:

(р1,ыьVI,Рх) = (1,0,0,1), (р2,«2,V2,Р2) = (0.125,0,0,0.1).

Здесь Р1, М1, Vl, Р1 и Р2, «2, V2, Р2 - плотность, компоненты вектора скорости и давление в канале слева и справа соответственно.

Задача рассматривалась в следующей области: [-1,1] х[-1,1]. Разрыв располагался вдоль прямой х = 0. Расчет велся до времени ^ = 0,3.

На рис. 2-4 представлены значения плотности, скорости и давления. Полученные результаты сравнивались с результатами расчетов с использованием схемы первого порядка точности.

Рис. 2. График плотности

Рис. 3. График давления

Рис. 4. График скорости Заключение

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

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

Список литературы

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

2. Hu, C. Weighted Essentially Non-oscillatory Schemes on Triangular Meshes / C. Hu, C.-W. Shu // Journal of Computational Physics. - 1999. - Vol. 150, № 1. - P. 97-127.

3. Zhang, Y.-T. Third order WENO schemes on three dimensional tetrahedral meshes / Y.-T. Zhang, C.-W. Shu // Communications in Computational Physics. - 2009. -Vol. 5. - P. 836-848.

4. Shi, J. A technique of treating negative weights in WENO schemes / J. Shi, C. Hu,

C.-W. Shu // Journal of Computational Physics. - 2002. - Vol. 175. - P. 108-127.

5. Sod, G. A. A survey of several finite difference methods for systems of nonlinear hyperbolic conservation laws / G. A. Sod // J. Comput. Phys. - 1978. - Vol. 27. -P. 1-31.

References

1. Godunov S. K., ZabrodinA. V., M. Ivanov Ya., Krayko A. N., Prokopov G. P. Chislen-noe reshenie mnogomernykh zadach gazovoy dinamiki [Numerical solution of multidimensional problems of gas dynamics]. Moscow: Nauka, 1976, p. 400.

2. Hu C., Shu C.-W. Journal of Computational Physics. 1999, vol. 150, no. 1, pp. 97-127.

3. Zhang Y.-T., Shu C.-W. Communications in Computational Physics. 2009, vol. 5, pp. 836-848.

4. Shi J., Hu C., Shu C.-W. Journal of Computational Physics. 2002, vol. 175, рр. 108127.

5. Sod G. A. J. Comput. Phys. 1978, vol. 27 рр. 1-31.

Горюнов Владимир Александрович

доктор физико-математических наук, профессор, кафедра экспериментальной физики, Мордовский государственный университет имени Н. П. Огарева

E-mail: gorval1934@mail.ru

Жалнин Руслан Викторович

кандидат физико-математических наук, доцент, кафедра прикладной математики, Мордовский государственный университет имени Н. П. Огарева

E-mail: zhalnin@gmail.com

Пескова Елизавета Евгеньевна ассистент, кафедра прикладной математики, Мордовский государственный университет имени Н. П. Огарева

E-mail: lizanika@mail.ru

Тишкин Владимир Федорович

доктор физико-математических наук, профессор, заместитель директора по научной работе, Институт прикладной математики имени М. В. Келдыша Российской академии наук

E-mail: office@keldysh.ru

Goryunov Vladimir Aleksandrovich Doctor of physical and mathematical sciences, professor, sub-department of experimental physics, Ogaryev Mordovia State University

Zhalnin Ruslan Viktorovich Candidate of physical and mathematical sciences, associate professor, sub-department of applied mathematics, Ogaryev Mordovia State University

Peskova Elizaveta Evgen'evna Assistant, sub-department of applied mathematics, Ogaryev Mordovia State University

Tishkin Vladimir Fedorovich Doctor of physical and mathematical sciences, professor, deputy director for reseach, Keldysh Institute of Applied Mathematics Russian Academy of Sciences

УДК 519.6 Горюнов, В. А.

О построении WENO схем для гиперболических систем уравнений на неструктурированных сетках / В. А. Горюнов, Р. В. Жалнин, Е. Е. Пескова, В. Ф. Тишкин // Известия высших учебных заведений. Поволжский регион. Физико-математические науки. - 2014. - № 1 (29). - С. 79-87.