Моделирование газодинамических течений на основе лагранжево-эйлеровой схемы les-asg Текст научной статьи по специальности «Физика»

www.volsu.ru

КОМПЬЮТЕРНОЕ МОДЕЛИРОВАНИЕ

DOI: http://dx.doi.org/10.15688/jvolsu1.2015.5.3

УДК 519.6 ББК 22.193

МОДЕЛИРОВАНИЕ ГАЗОДИНАМИЧЕСКИХ ТЕЧЕНИЙ

_ ____о____о ____л

НА ОСНОВЕ ЛАГРАНЖЕВО-ЭЙЛЕРОВОЙ СХЕМЫ LES-ASG1

Антон Владимирович Белоусов

Студент Института математики и информационных технологий, Волгоградский государственный университет anton.belousov.v@mail.ru, math@volsu.ru

просп. Университетский, 100, 400062 г. Волгоград, Российская Федерация

Сергей Сергеевич Храпов

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

доцент кафедры информационных систем и компьютерного моделирования, Волгоградский государственный университет xss-ip@mail.ru, infomod@volsu.ru

просп. Университетский, 100, 400062 г. Волгоград, Российская Федерация

Аннотация. Приведены формулы для численного моделирования газодинамических течений на основе лагранжево-эйлеровой схемы LES-ASG в двумерном случае. Приводится объяснение отличительной особенности численной схемы LES-ASG от схемы cSPH-TVD [1; 2]. Подробно рассматрива-g ется как лагранжев, так и эйлеров этап. Рассмотрены условия устойчивости

^ численной схемы. Проведен анализ и сравнение с численной схемой MUSCL

в одномерном случае.

m Ключевые слова: численные схемы, LES, ASG, cSPH-TVD, лагранжево-

g эйлерова схема.

а

* Введение

m

< В данной работе мы подробно рассмотрим все этапы и формулы численной схемы

о LES-ASG, которая является модификацией численной схемы cSPH-TVD. Проведем срав-

о нение схемы LES-ASG и схемы MUSCL для одномерного случая. В дальнейшем данная

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

1. Основные уравнения

Будем исходить из интегральных законов сохранения массы для однородной несжимаемой жидкости, а также законов сохранения энергии и импульса для «жидкой частицы» с «объемом» £(£), деформирующейся произвольным образом в процессе движения, в двумерном случае:

— JJ р dxdy = 0,

т

(1)

^ jj pudx dy = - JJ (^X - P f^ dx dy,

т

т

(2)

II P^ dxdy = - II (jy - P dxdy,

m m

itUedxdy=-U

m m

плотность среды; u и v —

(

d(up) д (vp)

x

+

- pufx - pvfy) dxdy,

)

(3)

(4)

где р — плотность среды; и и V — скорость газа вдоль ординат х и у соответственно;

д 4

е — полная энергия единицы объема; р — изотропное давление; /х = — и /„ =

д х

д4 , „

= —— удельная потенциальная внешняя сила; 4 — гравитационный потенциал.

Система уравнений (1)-(4) замыкается уравнением состояния идеального газа р = (у — — 1)ре, где у — показатель адиабаты, а е — удельная внутренняя энергия.

Введем характерное значение плотности р0, а также пространственную /0, и временной ¿0 масштабы задачи. Далее будем использовать только безразмерные величины (р, и, V, р, е, 4), переопределив их следующим образом

Ро

Р = — , u

u о

,

v

vt о 1о '

4=40 4 /2 •

pt0 Pol2;

et0_ Pol 2;

2. Численная схема ЬЕБ-ЛБв

Численная схема ЬЕБ-АБО (лагранжево-эйлерова схема с антисимметричной сеткой) — это модифицированный еБРИ-ТУЭ подход [1;2;7]. Основной отличительной особенностью схемы ЬЕБ-АБО от схемы еБРИ-ТУЭ является различие в расчете лагран-жевого этапа. В еБРИ-ТУЭ методе на лагранжевом этапе для расчета каждой ячейки задействуются все близлежащие ячейки, как это показано на рисунке 1. Из-за этого возникают нежелательные осцилляции. А в схеме ЬЕБ-АБО при расчетах не задей-ствуются диагонально расположенные соседние ячейки, как это делается на эйлеровом этапе. Пример показан на рисунке 2.

Рис. 1. Нахождение значений в ячейке на лагранжевом этапе с использованием cSPH-TVD схемы

Рис. 2. Нахождение значений в ячейке на лагранжевом этапе с использованием LES-ASG схемы

Воспользуемся стандартными процедурами дискретизации сплошной среды, применяемыми в численных схемах, основанных на эйлеровом и лагранжевом подходах [3;4]. Покроем расчетную область равномерной эйлеровой (неподвижной) сеткой с пространственным шагом h, где h^- = h = const (г и j — пространственные индексы) и х0+1 = х® + hi, у®+1 = у0 + hj. Совместим, в начальный момент времени, частицы с ячейками эйлеровой сетки. Число ячеек и частиц равно N.

Введем временные слои tn+1 = tn + тп с неравномерным шагом тп. На первом лагранжевом этапе, используя модифицированный в работе [6] SPH-подход, рассчитываем изменения интегральных характеристик и положений частиц, обусловленные действием газодинамических и внешних сил.

Частицы будем характеризовать массой М, импульсом Р и энергией Е соответственно: Mij(t) = //s..(t) pdxdy; Pid(t) = //s..(t) puvdxdy; Eid(t) = //s..(t) edxdy. После лагранжева этапа необходимо вернуть частицы в исходное состояние Ах™+1 — х0, — у®, вычислив при этом изменение интегральных характеристик частиц (M^f1, P^f1, ^jf1), вызванное таким перемещением. Именно на эйлеровом этапе возникает нужда использования неподвижной сетки для расчета потоков массы, импульса и энергии на границах ячеек в момент времени tnf1/2 = tn+тп/2, обусловленных перетеканием вещества через границы ячеек. Соответствующее изменение интегральных характеристик частиц пропорционально разности втекающих и вытекающих в ячейку потоков. Для расчетов потоков применяется модифицированный в [6] TVD-подход и приближенное решение задачи Римана.

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

чет в области течения нестационарных границ «вещество — вакуум».

2.1. Лагранжев этап

Введем дополнительную вспомогательную функцию ф, которая удовлетворяет соотношению р= ф2/2. Посредством введенной функции преобразуем в уравнениях (2)-(4) величины:

дф д дф д( и ) иф дф ф д( иф)

д

х

ф

дх' ду Фду, дх

2 дх 2 дх

д(ор) Vф дф + ф д(vф)

(5)

д 2 д 2

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

^ 0 = -1

(

\

А(х, Ь)йхйу

(6)

/

являющихся аналогом средних значений функции А = (р, ри, рv,e,р, ф,4) в ячейках, при конечно-объемной аппроксимации на неподвижной сетке.

Подставляя в систему (1)-(4) соотношения (5) и учитывая (6), преобразовав при этом интегральные характеристики частиц с учетом (5), получим

¿И

г,3

Ц

г,3>

(7)

где

и.

г,3

( рг,3 \ (ри)г,3

\ ег,3 )

Р = _ & =_ 1г'] дх, 1г'] ду,

Ц

г,3

0

дф

дхг дф

фг,3 ^ — рг,^

ф

г,3

и

дф д( иф)

г,3

д хг

+

хг

V

— (ри)3 & — (рп\з £

дф | д(vф),

% .

ду.

г,3

/

и

(р и)

г,3

(ро)г,з

— и иг,3 — скорости центра масс частиц. значение величины фг,^

рг,з рг,3

можно выразить через компоненты вектора консервативных переменных Ига- в виде

и Ог

— скорости центра масс частиц. Значение величины фг консерв

фг,з = л/2р г,з, (8)

где

,

(у — 1)

(е -

I с г,]

рг,3 2

2

Кз рг,з 2

)

(9)

2

Для аппроксимации пространственных производных, входящих в уравнение (7), будем использовать модифицированный БРИ-подход со сглаживающим ядром Ш [17; 18]. В качестве сглаживающего ядра будем использовать кубический сплайн Монагана:

Ш (д)

2 —

3

1 — 3 я2 + У3, 0 <д< 1;

4(2—^ 0,

1 < < 2; 2 < д.

(10)

Ш(д)

где д можно представить как дх (10) и (11), следует, что

2 —

3

9 2

—3д + 7 д2,

3

4

—4(2 — ^, 0,

1хг хк \ к

или ду

0 < д < 1;

1 <д< 2;

2 < ,

\ Уг — Ук \ к

(11)

. Исходя из соотношения

дШ г,к дШ (д)дд ' 8%дп(хг — хк)

Ш (д)—

д хг

д д х

к

дШ3,к _ дШ(1) д1 _ ш' (д)^1дп(Уз — Ук)

дУз

д д

где

Ш г,к = Ш (\хг — хк\,к)

= Ш (\у, — Ук\,к)

к

Ш г,к = Ш (\х? — хк\,к),

Ш,к = Ш(\у] — ук\, к).

(12)

(13)

(14)

(15)

Заменяя входящие в уравнение (7) пространственные производные конечными суммами, содержащими аналитически вычисляемые производные от сглаживающего ядра Ш, получим

Ц

,

г+1

ф , ф к,

к=г-1

3 + 1

0

дШ г,к + Рг,,, дШкк

хг

ф ,

к= -1

+ф ,

ф ,

к

+1

М

= -1

+1

М

к= -1

ф , к

ф ,

дШ^к , рг,з

4к- "Ж

к= -1

ф к,

дУ3

+

ф ,

4

дШ

, к

, к

дУз

иг,о + ик,] дШг,к . (ри)г^.,. г,к

д х

+

ф ,

4

дШ ^

к,

+0г,кдШ-1,к (ро) ^ , +- 4 г ,к

д хг

дШ

}

+

, к

2

%

ф

,

%

}

(16)

На шаге «предиктор» находим методом ломаных промежуточные значения и , в

2

о

*

момент времени tn+i = tn + т„

U.

г,3

(

Ura. + — + 2h

0

(фг+ij - Vi-I,j) - Pi,j(4i+i,j - 4i-i,j) -<Vi,j(<$i,j+i - Vi,j-i) - Pi,j(4ij+i - 4i,i-i)

ф

i,3

Фг+i,

Ui,j + Ui+i,j

- I

фг,3 Vi,3+i

2

Vi,j + Vi,j+i

фг-i,

Щ, +Ui-1,

<Pi,i-i

2

Vi,j + Vi,j-i

(17)

2 ^ 1 2 V -(Рикз - и ) - Мч СФм+1 - 1) У

знак «~» над вектором консервативных переменных И^- означает, что центр масс частиц смещен относительно исходного положения. На шаге «корректор» по наклону интегральной кривой в точке И^- и вычисляем приращение значений И^- и на временном слое гп+1 = + тп:

U

,n+i

Uj + U

1,3

1,3 + f Q« (U,x.

к, xk, Ук)

Ax.

1,3

A

1,3

т

T

Ui,j + Uk,j

Vn•+ V* ■ ъ,3 г,3

(18)

2

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

*

*

*

2

2

2.2. Эйлеров этап

На эйлеровом этапе вычисляются потоки массы, импульса и энергии, обусловленные перемещением жидкости через границы ячеек, в момент времени 1п+у2 = tn + тп/2. На данном этапе находится приближенное решение задачи Римана. Разность втекающих и вытекающих в ячейки потоков позволяет определить изменения характеристических «жидких» частиц, рассчитанных на предыдущем этапе в момент времени Ьп+1 = 1п + тп [14]. Представляя (1)-(4) в дифференциальной форме, а затем в консервативной эйлеровой форме при отсутствии сил газодинамического давления, получим:

дЬ дх ду

Применив стандартную процедуру конечно-разностной аппроксимации к уравнению (19), для г,]-й ячейки получим соотношение:

т тга+1 _ "Гт Тп (т?п+У2 Т71га+1/Л _ (г^п+1/2 г,п+1/2\ (20)

ЦЬЗ = ИЬЗ ^ {Fi+1/2 Г г-1/2 ) ^ ^+1/2 ^3-1/2) , (20)

Г п+1

где значения И^- вычисляются на лагранжевом этапе по формуле (17), а значения

тлП+1/2 -гл,т-тП+У2 ч г^п+1/2 ^,-¡^+1/2 ч

потоков на границах ячеек = с(Иi±1|2j) и = ^(И^^±1/2) находятся

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

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

Запишем выражения для и(х,у,1) слева (Ь) и справа ( К) от границы хс0+у2- и

1 /° •

У 1,3+1/2-

и 1+1/2,3

^+1/2,^ — иг+1, | 2 +

и., +12 -

+1/2 (ь Дх£Л

2 I ,

+1/2 (2 + дхз^«»

^2+ 2 ¡^х^1,

■

(21)

где

и

г,3+1/2

и

п+1/2 ■ ^ - Д^м^ е

(2+ДТ)«

1,3

2 -

и^З+1/2 — иг^+1 [ 2 +

,га+1/2

иу г,3+1,

Дх* ■

Дхп+1 — ДХг,з + Дхг,3 2 + 2

тп и™ + и.

2

п+1

1 2 ' 2 2

,п+1 _

+

~п+1 Т Уп + V ■

(22)

(23)

Наклоны кусочно-линейного распределения (21) и (22) должны удовлетворять ТУЭ-условию [12]. Для этого воспользуемся функцией-ограничителем [20]

@п — £ I 'и1+1,3 - и1з и1з - иг-1,3

®Пг ? — С

у

2

2

4,3+1

2

2

(24)

(25)

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

С(а, Ь) — 1 [Б1дп(а) + 81§п(6)] шт(|а|,

ограничитель ван Лира [15; 16]

2аЬ

С{а,Ь) — { ^ГЬ, аЬ> 0 0, аЬ < 0,

(26)

(27)

ь

ограничитель ван Альбады [8]

функция вирегЬее

_ (а2 + Qb+jb2 + Qa .28. С(а, Ь)_-а2 + Ь2 + 2Z-' (28)

С(а, b) _ max(a, b) max[min(1, 2 к), fmin(2, к)],

min(a,b) (29)

к

max(a, b)'

где а и Ь — вектора наклонов Q распределения величины и внутри ячейки; С — малая константа. Перечисленные ограничители (26)-(29) удовлетворяют ТУЭ-условию, ограничитель штшс^ (26), кроме того, сохраняет монотонность реконструируемой функции.

С учетом (21), (22) и методов приближенного решения задачи Римана (ЬЕ, ИЬЬ, ИЬЬС) [21], можно вычислить значения потоков и С™_+11/|, входящих в уравнение

(20). Рассмотрим значение потоков Г^/^2, введем следующие обозначения

иь,к — и±/2, — гсад,

тогда для потоков на границах ячеек, используя метод Лакса — Фридрихса (ЬЕ), имеем

[9]:

^п+1/2 Гь + ¥К иь — ип Ь ±1/2 — -2- + -2-, (30)

где

( РцвЦцв. \

F

L,R

2

PL,R ui,R + PL,R

PL,R UL,RVL,R

V uL}R(el,R + PL,R) J

Б* — шах(|5'ь|, |) — скорость распространения единственного разрыва, разделяющего две области с постоянными значениями и^, Для минимальной Бь и максимальной Бн скоростей распространения волн внутри ячейки справедливы оценки [10; 11]

— ш\п(иь — сь,ик — сЕ), — шах(иь + сь,ик + ск),

(ри)ь^ I УРцв. л где ULR — -, CLR — л/- — адиабатическая скорость звука в газе; pLR —

Рь^ у Рь, в.

, -,4 / и2ькрь,К У2ькрЬ>и\ г^п+1/2

— (у — 1) I е ь^-------- I. Аналогично находим значения потоков 2,

используя вместо основного пространственного индекса индекс .

В методе ИЬЬ значения потоков для уравнения (20) определяются по системе, представленной ниже [13]

F

Fl, 0 < SL;

<n+\/2 _ SRFL — SLFR + SLSK(UL — Uß) c c ,Qn

i±i/2 _ \-^-^-, sl s 0 s sR; (31)

I SR — SL

1

Fr, 0 > SR.

Значения Бь и Бк вычисляются аналогично методу ЬЕ.

Метод НЬЬС предпологает расчет потоков исходя из системы

Г

,п+1/2 г±1/2

П ) п )

Бь > 0; Бк < 0; Бс > 0; < 0.

Значение П (^ПИ/^ находим через

где

Вл

Бс (БьРь - Гь)

Бь Б<

Вз

В В2

Вз

Бс(Бь(ру)ь - Гь)

с

Бь Б{

с

В2

Бс(Бь(ри)ь - Гь) + Бь(рь + Рь(Бь - иь)(Бс - иь))

Бь Б<

с

В

Бс(Бьеь - Гь) + Бь(рь + Рь(Бь - иь)(Бс - иь))Бс

4 —

Бь Б

С

Значение П (тП+т/^ находим через

п (*п+:г?) —

В2 Оз

\»4)

где

О

2—

О,

п _ Бс(БкРк - Гк) _ Бс(БкМк - ¥к) —-^-^-, —-^-^-,

Бк - БС Бк - БС

Бс (Бк(ри)к - Г к) + Бк(рк + Рк(Бк - ик )(БС - ик))

Бк - Бс

Бс(Бкек - Гк) + Бк(рк + Рк(Бк - ик)(Бс - ик))Бс

Бк - Бс

В отличие от НЬЬ нам потребуется еще найти Бс

Б

с

Рк - Рь + Р ьиь(Бь - иь) - Ркик(Бк - ик) р( Бь - иь) - Рк(Бк - ик)

где Бь и Бк вычисляются по формулам

Бь — иь - сь<х(рь,рс), Б к — ик - ска(р к, Рс), Ь < а;

а(а, Ь)

* - О

Ь > а.

(32)

(33)

(34)

(35)

(36)

Акустическое приближение давления рс находится как

рс = max ^0,1 (pL + pR - 4(ur - uL)(pL + Pr)(cL + cr)^ .

(37)

2.3. Условия устойчивости

Для устойчивости численной схемы LES-ASG необходимо, чтобы за время интегрирования тга:

1) на лагранжевом этапе центр масс частиц смещался на расстояние, не превышающее h/2 относительно начального положения;

2) на эйлеровом этапе возмещения распространялись на расстояние, меньшее размера ячейки h.

С учетом этих условий, временной шаг тга для алгоритма LES-ASG должен определяться из условий:

т„ = К min

n Ii-1

' Vkjl

h

h

/

(38)

и VI | + я/ 1'°г,3 I + °г,з. где 0 < К < 1 — число Куранта; с — адиабатическая скорость звука в газе.

3. Проведение вычислительных экспериментов

Проведем несколько различных тестов. В качестве ограничителя используем min-mod, а для приближенного решения задачи Римана будем использовать метод ЬЕ.

Рассмотрим тест (А) для сгустка плотности в виде окружности на квадратной расчетной области с размерностью N — 102. В качестве начальных значений: показатель адиабаты у — 1,4; в областях задается однородная плотность р(х, у) — 1; и(х, у) — 0; у(х, у) — 0; х и у Е [0,1], а энергия:

Р

e(x, у)

10у Р

. Y

2 N

R2 <Г-

2 N R2 ä N •

где Я — это радиус отклонения от центра расчетной области. На рисунках 3-8 представлены результаты эксперимента.

Рассмотрим тест (В) о взаимодействии двух взрывных волн, который был применен Вудвартом [22] и Колеллой для изучения свойств численной схемы РРМЬР годуновского типа, реализующий счет давления на лагранжевом этапе [22]. Начальное состояние состоит из трех областей, с показателем адиабаты у — 1, 4 и ограниченных твердыми стенками. В областях задается однородная плотность р(х, у) — 1, и(х, у) — 0, у(х, у) — — 0, х и у Е [0,1], а разрыв по давлению:

(1000, 0 <д< 1;

р(х, у) — \ 0, 01, 0,1 < х < 0, 9; [100, 0, 9 < х < 1.

В результате формируются две сильные сходящиеся ударные волны. На рисунке 9 показана динамика взаимодействия двух взрывных волн для расчетной области N — 400.

Рис. 5. Тест А, распределение плотности после 10-й расчетной итерации

Рис. 6. Тест А, распределение плотности после 30-й расчетной итерации

х

Рис. 7. Тест А, распределение энергии после Рис. 8. Тест А, поле скоростей после

30-й расчетной итерации 30-й расчетной итерации

Рассмотрим тест (С), в котором смоделируем задачу Римана о распаде произвольного разрыва. Зададим начальные условия для квадратной расчетной области размером

N — 103, с левой границы от разрыва рь — 1, рь — 1 ис правой границы от разрыва рк — 1, рк — 0,1, распределение скоростей и и V нулевые на всей расчетной области. Показатель адиабаты у — 1,4. Данный тест применяется для оценки устойчивости численной схемы при моделировании сильных ударных волн. Решение состоит из сильной ударной волны, контактной волны и левой волны разряжения. На рисунках 10 и 11 представлены профильные разрезы давления и плотности. Функция-ограничитель тттоА позволяет получить профиль наиболее приближенным к монотонному, при этом масштабы численной размазки получаются более существенными.

Рис. 9. Тест В, динамика взаимодействия двух взрывных волн для N = 400: графики распределения давления р(х,Ы/2) в разрезе по центру расчетной области, вдоль движения волн

Рис. 10. Тест С, распределение давления при моделировании распада произвольного разрыва

Рис. 11. Тест С, распределение плотности при моделировании распада произвольного разрыва

В работе [1] проводилось сравнение погрешностей численных схем ЬЕБ-АБО и МиБСЬ с использованием кусочно-линейной и кусочно-постоянной реконструкции. В качестве основного теста рассматривалось решение уравнения переноса с неоднородным распределением скорости. В таблице приведены результаты эксперимента.

Погрешность вычислений схем МУБСЬ и ЬББ-ЛБО [1]

N Кусочно-постоянная реконструкция Кусочно-линейная реконструкция

МиБСЬ ЬЕБ МиБСЬ ЬЕБ

300 7,5455 х 10-2 7,5451 х 10-2 2,0163 х 10-2 1,0191 х 10-2

600 4,2171 х 10-2 4,2169 х 10-2 5,7162 х 10-3 5,7275 х 10-3

1200 2,2633 х 10-2 2,2633 х 10-2 2,6335 х 10-3 2,6309 х 10-3

2400 1,1818 х 10-2 1,1818 х 10-2 1,2402 х 10-3 1,2411 х 10-3

4800 6,0794 х 10-3 6,0794 х 10-3 5,9244 х 10-4 5,9223 х 10-4

9600 3,0903 х 10-3 3,0903 х 10-3 2,7158 х 10-4 2,7163 х 10-4

Заключение

По результатам численных экспериментов с использованием схемы ЬЕБ-АБО, для моделирования газодинамических течений, можно сделать вывод, что модификация схемы еБРИ-ТУЭ проведена успешно. Возмущения, которые возникали в угловых ячейках при использовании еБРИ-ТУЭ-метода, значительно уменьшены и профиль сечения имеет более гладкое распределение. Сравнивая результаты вычислений с использованием метода ЬЕБ-АБО и МиБСЬ, можно с уверенностью говорить о схожести точности численных схем [1]. В дальнейшем численная схема ЬЕБ-АБО будет использована для создания программного комплекса, нацеленного на моделирование газодинамических течений, который в будущем можно задействовать как основу для развития и создания более сложных и точных модулей, нацеленных на моделирование газодинамических течений.

ПРИМЕЧАНИЕ

1 Работа выполнена при поддержке грантов РФФИ № 15-45-02655, № 15-47-02642, № 15-02-06204, № 14-07-97030.

СПИСОК ЛИТЕРАТУРЫ

1. Белоусов, А. В. Разработка программы для численного газодинамического моделирования на основе лагранжево-эйлеровой схемы ЬЕБ-АБО / А. В. Белоусов, С. С. Храпов // Вестник Волгоградского государственного университета. Серия 1, Математика. Физика. — 2015. — № 1 (26). — С. 30-39.

2. Жумалиев, А. Г. Численная схема еБРИ-ТУО: моделирование фронта ударной волны / А. Г. Жумалиев, С. С. Храпов // Вестник Волгоградского государственного университета. Серия 1, Математика. Физика. — 2012. — № 16. — С. 24-27.

3. Еремин, М. А. Конечно-объемная схема интегрирования уравнений гидродинамики / М. А. Еремин, А. В. Хоперсков, С. А. Хоперсков // Изв. Волгогр. гос. техн. ун-та. —

2010. - № 6:8. - C. 24-27.

4. Кузьмин, Н. М. Численное моделирование эволюции неустойчивых мод джетов, выходящих из молодых звездных объектов / Н. М. Кузьмин, В. В. Мусцевой, С. С. Храпов // Астрон. журн. - 2007. - № 84:12. - C. 1089-1098.

5. Куликовский, А. Г. Математические вопросы численного решения гиперболических систем уравнений / А. Г. Куликовский, Н. В. Погорелов, А. Ю. Семенов. - М. : ФИЗМАТ-ЛИТ, 2001. - 656 с.

6. Численная схема для моделирования динамики поверхностных вод на основе комбинированного SPH-TVD-подхода / C. C. Храпов, А. В. Хоперсков, Н. М. Кузьмин, А. В. Писарев, И. А. Кобелев // Вычислительные методы и программирование. - 2011. -Т. 12, № 1. - C. 282-297.

7. Численная схема CSPH-TVD: исследование влияния ограничителей наклонов / Н. М. Кузьмин, А. В. Белоусов, Т. С. Шушкевич, С. С. Храпов // Вестник Волгоградского государственного университета. Серия 1, Математика. Физика. - 2014. - № 1 (20). -C. 22-34.

8. Albada, G. D. van. A comparative study of computational methods in cosmic gas dynamics / G. D. van Albada, B. van Leer, W. W. Roberts // Astron. Astrophys. - 1982. -P. 76-84.

9. Courant, R. On the solution of nonlinear hyperbolic differential equations by finite differences / R. Courant, E. Isaacson, M. Rees // Comm. Pure. - 1952. - P. 243-255.

10. Davis, S. F. Simplified Second-Order Godunov-Type Methods / S. F. Davis // SIAM J. Sci. Stat. Comput. - 1988. - P. 445-473.

11. Einfeldt, B. On Godunov-Type Methods for Gas Dynamics / B. Einfeldt // SIAM J. Numer. Anal. - 1988. - P. 294-318.

12. Harten, A. High Resolution Schemes for Hyperbolic Conservation Laws / A. Harten // J. Comput. Phys. - 1983. - P. 357-393.

13. Harten, A. On upstream differencing and Godunov type methods for hyperbolic conservation laws / A. Harten, P. Lax, B. van Leer // SIAM review. - 1983. - P. 35-61.

14. Hudson, J. Numerical techniques for conservation laws with source terms / J. Hudson // PhD Thesis. - Reading : University of Reading, 1998. - P. 1-118.

15. Leer, B. van. Towards the Ultimate Conservative Difference Scheme II. Monotonicity and Conservation Combined in a Second Order Scheme / B. van Leer // J. Comput. Phys. -1974. - P. 361-370.

16. Leer, B. van. Towards the Ultimate Conservation Difference Scheme V. A Second Order Sequel to Godunov's Method / B. van Leer // J. Comput. Phys. - 1979. - P. 110-136.

17. Monaghan, J. J. Simulating free surface flows with SPH / J. J. Monaghan // Comput. Phys. - 1994. - P. 399-406.

18. Monaghan, J. J. Smoothed Particle Hydrodynamics / J. J. Monaghan // Annual Review of Astronomy and Astrophysics. - 1992. - P. 543-574.

19. Roe, P. L. Some Contributions to the Modelling of Discontinuous Flows / P. L. Roe // Proceedings of the SIAM/AMS Seminar. - 1983. - P. 163-193.

20. Sweby, P. K. High Resolution Schemes Using Flux Limiters for Hyperbolic Conservation Laws / P. K. Sweby // SIAM J. - 1984. - P. 995-1011.

21. Toro, E. F. Restoration of the Contact Surface in the HLL / E. F. Toro, M. Spruce, W. Speares // Shock Waves. - 1994. - P. 25-34.

22. Woodward, P. The Numerical Simulation of Two-Dimensional Fluid Flow with Strong Shocks / P. Woodward, P. Colella // J. Comput. Phys. - 1984. - P. 115-173.

REFERENCES

1. Belousov A.V., Khrapov S.S. Razrabotka programmy dlya chislennogo gazodinamicheskogo modelirovaniya na osnove lagranzhevo-eylerovoy skhemy LES-ASG [Development of the Program for Numerical Gasdynamic Modeling on the Basis of Lagrange —

Euler Scheme the LES — ASG]. Vеstnik Volgogradskogo gosudarstvеnnogo univеrsitеta. Sеriya 1, Matеmatika. Fizika [Science Journal of Volgograd State University. Mathematics. Physics], 2015, no. 1 (26), pp. 30-39.

2. Zhumaliev A.G., Khrapov S.S. Chislennaya skhema cSPH-TVD: modelirovanie fronta udarnoy volny [A Numerical Scheme Based on the Combined SPH-TVD Approach: Simulation of the Shock Front]. Vеstnik Volgogradskogo gosudarstvеnnogo univеrsitеta. Sеriya 1, Matеmatika. Fizika [Science Journal of Volgograd State University. Mathematics. Physics], 2012, no. 16, pp. 24-27.

3. Eremin M.A., Khoperskov A.V., Khoperskov S.A. Konechno-obyemnaya skhema integrirovaniya uravneniy gidrodinamiki [Finite Volume Sheme of Integration for Hydrodynamics Equations]. Izv. Volgogr. gos. tеkhn. un-ta, 2010, no. 6:8, pp. 24-27.

4. Kuzmin N.M., Mustsevoy V.V., Khrapov S.S. Chislennoe modelirovanie evolyutsii neustoychivykh mod dzhetov, vykhodyashchikh iz molodykh zvezdnykh obyektov [Numerical Modeling of the Evolution of Unstable Modes of Jets From Young Stellar Objects]. Astron. zhurn. [Astronomy Reports], 2007, no. 84:12, pp. 1089-1098.

5. Kulikovskiy A.G., Pogorelov N.V., Semenov A.Yu. Matеmatichеskiе voprosy chishnnogo rеshеniya gipеrbolichеskikh sistеm uravmniy [Mathematical Problems in the Numerical Solution of Hyperbolic Systems]. Moscow, FIZMATLIT Publ., 2001. 656 p.

6. Khrapov C.C., Khoperskov A.V., Kuzmin N.M., Pisarev A.V., Kobelev I.A. Chislennaya skhema dlya modelirovaniya dinamiki poverkhnostnykh vod na osnove kombinirovannogo SPH-TVD-podkhoda [Numerical Scheme for Modeling the Dynamics of Surface Water Based on the Combined SPH-TVD-Approach]. Vychislitеlnyе mеtody i programmirovaniе, 2011, vol. 12, no. 1, pp. 282-297.

7. Kuzmin N.M., Belousov A.V., Shushkevich T.S., Khrapov S.S. Chislennaya skhema CSPH-TVD: issledovanie vliyaniya ogranichiteley naklonov [Numerical Scheme CSPH — TVD: Investigation of Influence Slope Limiters]. Vеstnik Volgogradskogo gosudars^nnogo univеrsitеta. Sеriya 1, Matеmatika. Fizika [Science Journal of Volgograd State University. Mathematics. Physics], 2014, no. 1 (20), pp. 22-34.

8. Albada G.D. van., Leer B. van., Roberts W.W. A Comparative Study of Computational Methods in Cosmic Gas Dynamics. Astron. Astrophys., 1982, pp. 76-84.

9. Courant R., Isaacson E., Rees M. On the Solution of Nonlinear Hyperbolic Differential Equations By Finite Differences. Comm. Pure., 1952, pp. 243-255.

10. Davis S.F. Simplified Second-Order Godunov-Type Methods. SIAM J. Sci. Stat. Comput., 1988, pp. 445-473.

11. Einfeldt B. On Godunov-Type Methods for Gas Dynamics. SIAM J. Numer. Anal., 1988, pp. 294-318.

12. Harten A. High Resolution Schemes for Hyperbolic Conservation Laws. J. Comput. Phys, 1983, pp. 357-393.

13. Harten A., Lax P., Leer B. van. On Upstream Differencing and Godunov Type Methods for Hyperbolic Conservation Laws. SIAM review, 1983, pp. 35-61.

14. Hudson J. Numerical techniques for conservation laws with source terms. PhD Thesis, Reading, University of Reading, 1998, pp. 1-118.

15. Leer B. van. Towards the Ultimate Conservative Difference Scheme II. Monotonicity and Conservation Combined in a Second Order Scheme. J. Comput. Phys., 1974, pp. 361-370.

16. Leer B. van. Towards the Ultimate Conservation Difference Scheme V. A Second Order Sequel to Godunov's Method. J. Comput. Phys., 1979, pp. 110-136.

17. Monaghan J.J. Simulating Free Surface Flows with SPH. Comput. Phys., 1994, pp. 399406.

18. Monaghan J.J. Smoothed Particle Hydrodynamics. Annual Review of Astronomy and Astrophysics, 1992, pp. 543-574.

19. Roe P.L. Some Contributions to the Modelling of Discontinuous Flows. Proceedings of the SIAM/AMS Seminar, 1983, pp. 163-193.

20. Sweby P.K. High Resolution Schemes Using Flux Limiters for Hyperbolic Conservation Laws. SIAM J, 1984, pp. 995-1011.

21. Toro E.F., Spruce M., Speares W. Restoration of the Contact Surface in the HLL. Shock

Waves, 1994, pp. 25-34.

22. Woodward P., Colella P. The Numerical Simulation of Two-Dimensional Fluid Flow with Strong Shocks. J. Comput. Phys., 1984, pp. 115-173.

GASDYNAMIC MODELING ON THE BASIS OF THE LAGRANGIAN AND EULER SCHEME LES-ASG

Anton Vladimirovich Belousov

Student, Institute of Mathematics and IT, Volgograd State University anton.belousov.v@mail.ru, math@volsu.ru

Prosp. Universitetsky, 100, 400062 Volgograd, Russian Federation

Sergey Sergeevich Khrapov

Candidate of Physical and Mathematical Sciences, Associate Professor, Department of Information Systems and Computer Simulation, Volgograd State University xss-ip@mail.ru, infomod@volsu.ru

Prosp. Universitetsky, 100, 400062 Volgograd, Russian Federation

Abstract. In this work the numerical scheme LES-ASG which is modification of CSPH-TVD is considered in detail. Formulas of all stages of calculation with use of LF, HLL, HLLC methods for the solution of a task of Riemann are presented. Using LES-ASG it was succeeded to achieve more smooth distribution in comparison with CSPH-TVD. By results of comparison of LES-ASG and MUSCL it is possible to draw a conclusion on similarity of accuracy of numerical schemes. In this paper, compare the two types of solving the transport equation with inhomogeneous distribution of velocity, namely LES (Lagrange — Euler scheme) and MUSCL (Monotonic Upstream-Centered Scheme). According to the research we can assume that the scheme LES and MUSCL is equally well applicable for modeling the transport equation. The numerical scheme LES-ASG will be used as a basis for creation of the program complex aimed at modeling of gasdynamic currents with use of the OpenMP and CUDA technologies.

Key words: numerical schemes, LES, ASG, cSPH-TVD, MUSCL, Lagran-gian and Euler scheme.