Научная статья на тему 'Моделирование газодинамических течений на основе лагранжево-эйлеровой схемы les-asg'

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

CC BY
261
51
i Надоели баннеры? Вы всегда можете отключить рекламу.
Ключевые слова
ЧИСЛЕННЫЕ СХЕМЫ / LES / ASG / СSРН-ТVD / ЛАГРАНЖЕВO-ЭЙЛЕРOВА СХЕМА / NUMERICAL SCHEMES / CSPH-TVD / MUSCL / LAGRANGIAN AND EULER SCHEME

Аннотация научной статьи по физике, автор научной работы — Белоусов Антон Владимирович, Храпов Сергей Сергеевич

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

i Надоели баннеры? Вы всегда можете отключить рекламу.

Похожие темы научных работ по физике , автор научной работы — Белоусов Антон Владимирович, Храпов Сергей Сергеевич

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.
i Надоели баннеры? Вы всегда можете отключить рекламу.

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

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.

Текст научной работы на тему «Моделирование газодинамических течений на основе лагранжево-эйлеровой схемы 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)

где

и.

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

г,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

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

хг

ф ,

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

к= -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

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

где значения И^- вычисляются на лагранжевом этапе по формуле (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

П ) п )

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

Бь > 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

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

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.

i Надоели баннеры? Вы всегда можете отключить рекламу.