CC BY
9
УДК 519.6
DOI: 10.14529/mmp210212
МЕТОД ИССЛЕДОВАНИЯ ДИССИПАТИВНЫХ СВОЙСТВ РАЗНОСТНЫХ СХЕМ В ЭЙЛЕРОВЫХ КООРДИНАТАХ
Е.С. Шестаковская1, Я.Е. Стариков1, Н.Л. Клиначева1
1 Южно-Уральский государственный университет, г. Челябинск, Российская Федерация
В настоящее время численные методы расчета ударно-волновых течений жидкости и газа в эйлеровых координатах получили широкое распространение, поэтому исследование их характеристик является актуальной задачей. В работе представлен подход к оценке диссипативных свойств таких разностных схем на сильных разрывах. Идея метода заключается в построении уравнения производства энтропии, погрешность аппроксимации которого может быть выражена комбинацией погрешностей аппроксимации уравнений, составляющих разностную схему. В качестве критерия диссипативно-сти разностной схемы используется уравнение производства энтропии на слабой ударной волне. В работе проведена оценка диссипативных свойств метода крупных частиц с использованием предложенного метода.
Ключевые слова: разностные схемы; эйлеровы координаты; диссипативность; производство энтропии.
Введение
Общепринятыми характеристиками разностных схем являются аппроксимация и устойчивость. Методы исследования этих свойств разностных схем подробно изучены в [1]. Однако в разностных схемах, описывающих сильные разрывы, должен присутствовать механизм диссипации энергии, обеспечивающий возрастание энтропии. В результате его действия может быть создано не только необходимое изменение энтропии, но и излишнее, обусловенное свойствами конкретной разностной схемы. Возникает вопрос, какое изменение энтропии можно считать допустимым. Анализ вклада аппроксимационной вязкости в энтропию для ряда разностных схем приведен в [2]. В качестве критерия диссипативности разностной схемы предложено использовать уравнение производства энтропии на слабой ударной волне [3, 4]
1 /82 Р\
АУ3 + 0(АУ4,А52).
Таким образом, в зависимости от скорости роста энтропии на ударных волнах, схемы можно разделить на два класса: сильно диссипативные и слабо диссипативные. Физический смысл этого критерия весьма прост: разностная схема является приемлемой, если изменение энтропии из-за погрешностей аппроксимации не превосходит ее изменения в характерных физических процессах, какими являются слабые ударные волны. В сильно диссипативных разностных схемах слабые ударные волны неразличимы на фоне погрешностей. Иными словами, сильно диссипативные разностные схемы имеют низкую «разрешающую способность:», и наоборот, слабо диссипативные РС - высокую [5]. Метод исследования диссипативных свойств разностных схем в лагранжевых координатах [6], заключается в построении уравнения производства
энтропии, погрешность аппроксимации которого может быть выражена комбинацией погрешностей аппроксимации уравнений, составляющих разностную схему. В настоящее время разностные схемы, использующие эйлеров подход для описания сплошной среды все более активно применяются, особенно для решения многомерных задач. В этих разностных схемах, также как и в лагранжевых, используются различные механизмы диссипации энергии на сильных разрывах. Таким образом, исследование диссипативный свойств разностных схем в эйлеровых переменных является важной, но практически не изученной задачей. В настоящей работе показана возможность применения метода для схем в эйлеровых координатах на примере широко используемого метода крупных частиц.
1. Метод исследования
Запишем уравнение скорости изменения энтропии
дЛ, <Ы <Ы '
где Т - температура, Б - энтропия, е - удельная внутренняя энергия, Р - давление, V - удельный объем. Поскольку вдоль траектории каждой частицы сохраняются еч лагранжевы координаты, то полная производная какой-либо величины совпадает с еч частной производной в этих координатах. В соответствии с [6] погрешность аппроксимации уравнения производства энтропии имеет вид
аэ т'чд2р\ (ву\3
Правая часть тождества (1) может быть выражена различными способами через комбинацию погрешностей аппроксимации уравнений из исходной системы, либо же их следствий. Выберем вариант разностной схемы, содержащий погрешности уравнений неразрывности, движения и энергии в лагранжевых координатах
дУ ди _ ди дР _ <9£ дтп 1' дt дтп 2'
Тогда погрешность аппроксимации уравнения производства энтропии ш3 имеет вид
дБ
Т— = ш8 = Рил + ш3 - 11ш2, (3)
где и - скорость.
Запишем систему уравнений в эйлеровых переменных
др дРи
= (4)
дрЦ дрЦЦ (5) _<9£_дх дх Ш2'_
Вестник ^ЭУрГУ. Серия «Математическое моделирование 109
и программирование» (Вестник ЮУрГУ ММП). 2021. Т. 14, № 2. С. 108-116
дрЕ дрЕЦ дРЦ _ <9£ дх дх Шз' где Е = е + 0, 5и2 - удельная полная энегия.
Покажем, что для системы (4) - (6) выражение (3) имеет тот же вид. Для уравнения неразрывность проведем преобразования, доказывающие возможность перехода между двумя формами. Начнем с его записи в эйлеровых координатах
9р дрЦ_ = о
дЬ дх
Преобразуем второе слагаемое для выделения выражения субстанциональной производной
др др 817 дЬ дх ^ дх
Запишем первые два слагаемые как полную производную и сделаем замену р
¿р д£ _ Ц 1_д£ _
сИ ^ дх ' сИ V дх
После этого остается привести уравнение к требуемому виду путем несложных действий
С учетом т = рх, дт = рдх, уравнение (7) примет вид
¿У ди
V
dt dm
Окончательно получим вид уравнения в лагранжевых координатах
dV dU
дЬ дт
Более кратко изложим аналогичные преобразования уравнения движения и энергии. Рассмотрим уравнение движения в эйлеровых координатах
дрЦ дрЦЦ дР_ _ 0 дЬ дх дх
После преобразований получим
— + и— + с/— + и— + — - о
^ дЬ дЬ дх ^ дх дх
Слагаемые 2 и 3 образуют уравнение неразрывности, домноженное на и. Разделим уравнение на р для выделения субстанциональной производной
д£ ид£ дЬ дх р дх
Окончательно получим
dU дР dt dm
0
Закончим рассмотрение эквивалентности уравнением энергии. Проведем аналогичные действия. Запишем его в эйлеровых координатах
дрЕ дрЕЦ дРЦ _ 0 дЬ дх дх '
— + + идЕ | дри - О
^ дЬ дЬ дх ^ дх дх дЕ дРи_ дЬ дт
Таким образом, подстановка выражений (4) - (6) в (3) приводит к получению уравнения (1). Это показывает справедливость нахождения погрешности аппроксимации уравнения производства энтропии через комбинацию погрешностей исходных уравнений разностной схемы для методов в эйлеровых координатах.
2. Исследование одной разностной схемы в эйлеровых координатах
Реализуем данный подход применительно к методу крупных частиц [7]. Для упрощения ограничимся рассмотрением одномерного случая, а также примем направления потоков слева направо с их вычислением с использованием формул первого порядка точности. Рассмотрим разностные уравнения неразрывности, движения и энергии
Рг +1 ~ Рг Р?Ф? + Щг) ~ Р1-Ж-1 + У?) = д т 21г
рП+гцп+г _ рпщ рп(1пфп + (1п+1)_рп_1(1п_1фп_1 + (1п) рга^^га Т 21г 2к
рп+1Еп+1_рпЕп р^п + - р^ф^ + й'ПЕ^! ,
т 21г
, (РГ + РГ+гЖ + Щг) ~ 1 + РГЖ-1 + Щ) 0
Ак
Идея МКЧ подразумевает использование предварительных значений скорости и и энергии Р.
Пп — 11п__— (Рп — Рп )
иг ~ иг о „га I, ^ *+1 г—1/ '
2рП к
т
+ 1/2^+1/2 ~~ 1 1-\12и1-\12)—^Т
Рг к
рп _ рп ( рп ттп рп ттп )
рг _ рг — (Рг+1/2тг+1/2 — Рг—1/2тг—1/2)
1 1 т
Е? - (-ЛР„г + р^ж + и:+1) - -(р^ + рщиг-г + и:)у
г г+1/ г 1 ~ г+и ^^ г—1 1 ~ г г—1 < ~ г )) ^п^
Отметим, что в настоящей работе рассматривается метод крупных частиц без использования дополнительной искусственной вязкости, поэтому значения Р™+ х и и™ ± вычисляются как 2 2
Р?Ц = \(Р? + Р?+1),
Теперь приступим к разложению разностных уравнений в ряды Тейлора в точке (г, п) и выделению погрешностей аппроксимации. Для упрощения данного процесса, начнем с уравнений предварительной скорости и энергии. Ограничимся слагаемыми первого порядка малости по времени и пространству
т дП ~ дР
Ег = Е'г ~ ТГ^^Р'г^ + ^'г —К) + 0(Т2, к2). г г г дх г дх '
Вместе с этим потребуются первые и вторые производные предварительных величин по пространству. В результате получим
ди ди тд2Р 1 т дрдР , ^, 2
+ + (7)
dx дх р дх2 р2 дх дх d2U _ d2U т д3Р + ^дрд^Р + ^dj^dP_ + ^дрд^Р _ 2тдРдрдр + Q^2 ^
дх2 дх2 р дх3 р2 дх дх2 р2 дх2 дх р2 дх дх2 р3 дх дх дх
дЁ _дЕ 2 dPdU Р d2U Р dU dp U д2Р U дР dp
дх дх р дх дх р дх2 р2 дх дх р дх2 р2 дх дх
д2 E д2 E 3 д2 P ди 3 дР д2и 4 дрдидР P д3 U ид3 P
~ -т---^^г^гТ + ^г^—Т^^^т------7Г-ГТ+
дх2 дх2 р дх2 дх р дх дх2 р2 дх дх дх р дх3 р дх3
Р д2U U д2Р Р dp d2U U др д2Р Р dU д2р UdPd2p
р2 дх2 р2 дх2 р2 дх дх2 р2 дх дх2 р2 дх дх2 р2 дх дх2
Р_ Ш др_ _ U_ дР_ др_ + 2 р4 дх дх дх р4 дх дх дх '
Далее перейдем к нахождению погрешностей самих разностных уравнений. Продемонстрируем используемый подход на уравнении неразрывности
рТ1 - Рг , Р?(й? + йщ)-Р?-Лй?-1 + й?) др 1д2р 1 ~ ди
+-2h-= m + 2к{рЛЩ
1 д2 II12. , n др, 1 д2р12..п ~ ди, 1 д2 ¿7, 2чч 3
1 (On"d^h Д- OT>9ph 9й dPh2 find2Ph2, , n(2 h2N
Теперь заменим производные предварительной скорости по пространству ^ на ранее полученное выражение (7). В результате получим
РГ1 - Р? РгФ'г + U?+l) ~ Ptlfel + U?) = т 2 h
_др ndU_ dp 1сРр _ д^Р_ _ 1 dUdp _ lrjn&Ph , n(ji h2, dt+Pl дх+ г дх + 2dt2 дх2 2дхдх 2гдх2+[,)'
Ц2 Bulletin of the South Ural State University. Ser. Mathematical Modelling, Programming
& Computer Software (Bulletin SUSU MMCS), 2021, vol. 14, no. 2, pp. 108-116
Заметим, что первые три слагаемые в правой части являются исходным дифференциальным уравнением, следовательно, ошибку аппроксимации составляют оставшиеся члены. Окончательно заключим
Проведем аналогичные действия с уравнением движения, которое имеет следующую запись в разностном виде
р7+1ц?+1~рж р?йпй? + щг) - + и?) рг+1~рг-1=0
т 21г 2к '
Разложение в ряд Тейлора в точке (г, п) будет иметь вид
ди др 1 д2р 1 д2и дидр ди др
I ттп г , гт ^^ | ^ ^ | ^ I П^пттп^^ | „.пцп г
- К - V—— Л " ^Е/»—* -
1 дх2 <9ж 2 <9ж <9ж 2^г г &г2 2 г <9ж
Таким образом, выделяя лишние слагаемые, получим погрешность аппроксимации уравнения движения
1 тп д2 р 1 п д2 и дидр ттп д2Р ди дР 1 п ди ди ^ Ш2 = ~2 ~ 2^ " ¥¥Т + Щ ~д^Т + + 2Л , 1 „пТТпд2ии . 3 ттпди др 1 ттпттпд2р . 2 т,2\
Остается лишь провести вышеприведенные действия в отношении уравнения энергии. Рассмотрим его разностный вид
р?+1е:+1-Р?Е.г + +
г 2Н
, + ¿У+Ж^Г + РЙ-1) ~ №! + + Е?) = п
4/1
Запишем итоговое выражение для погрешности этого уравнения
3 = ~2 -д*Т ~ 2 " ¥¥Т + + ^ + 17 + + + рпд_ид_ит + „пдР^ + ^пди^ + цпЕпдр^
дх дх г дх2 г дх дх г дх дх 2 г дх дх г г дх
+ 2Р* 2 ' ' &г2 + 2Рг * &г2 + р дх2 Т р2 <9ж <9ж Т+
1 др др 2 2 р дх дх
Теперь с помощью уравнения (3) найдем итоговую погрешность ш3, объединив результаты для отдельных уравнений. Результатом является следующее выражение
ш3 = Риг+— Пи2 = 2 г дг2 д2 х 2 г дх дх 2 г г дх2 2 г дг2 2 г дг2
№ _итита_Хт+эрэЕ^р
дЬдЬ г дх дх г г дх2 г г дх2 дх дх г дх2 +рпди дирпди др др 1 ди дЕ 1 д2р 1 д2Е (8)
Рд^Рт _ РдрдР 1дРдР 1 ипип^Рт+1 ип0п^ит+ипдидр р дх2 р2дхдх рдхдх 2 г 1 Ы2 2 г Рг Ы2 1 дЬ дЬ 1 дП дП 1 д2Ц 3 дП дп 1 д2п
„пТТп Ъ „пТТпТТп Ъ "гтгт I Ъ ТТпТТпТТп 1 Ъ
~2Рг и* ах^ ~ 2Рг и* и* д^*1 - 2 и* д^^ ~ 2 и* и* д.
Для возможности сравнения полученной погрешности и допустимого физически обоснованного изменения энтропии на слабой ударной волне, рассмотрим преобразование выражения (2) к виду, содержащему только производные по времени. Для этого используются следующие выражения [4]
дР (дР\ дУ
dt \dV J S dt ' д2Р_(д2Р\ (dV\2 (дР\ d2V
дЬ2 \дУ2)3\дг) \dVJsdt2 После ряда действий, рассматриваемое уравнение принимает вид
12 V dV2) Д dt )
д2P др т2 дРд2р т2 дР (дрх 2
р3 dt2 dt 12р2 dt dt2 6р3 dt \дЬ
+ I тгг I +•••
(9)
Сравнение выражений (8) и (9) показывает превышение физически обоснованного изменения энтропии в методе крупных частиц. Это обусловлено наличием в итоговой погрешности слагаемых первого порядка по времени и пространству, в то время как выражение (9) ограничивается членами второго порядка. Следовательно, схема метода крупных частиц является сильно диссипативной на сильных разрывах.
Таким образом, в работе показана возможность применения рассматриваемого метода анализа диссипативных свойств численных методов в эйлеровых координатах, суть которого заключается в построении уравнения производства энтропии и нахождении его погрешности аппроксимации через комбинацию погрешностей аппроксимации исходных разностных уравнений.
Авторы выражают благодарность И.Р. Макеевой за значимые замечания и дискуссии в ходе написания научной статьи.
2
Литература
1. Яненко, Н.Н. О первом дифференциальном приближении разностных схем для гиперболических систем уравнений / Н.Н. Яненко, Ю.И. Шокин // Сибирский математический журнал. - 1969. - Т. 10, № 5. - С. 1173-1202.
2. Рождественский, Б.Л. Системы квазилинейных уравнений / Б.Л. Рождественский, Н.Н. Яненко. - М.: Наука, 1968.
3. Зельдович, Я.Б. Физика ударных волн / Я.Б. Зельдович, Ю.П. Райзер. - М.: Физматгиз, 1963.
4. Забабахин, Е.И. Некоторые вопросы газодинамики взрыва / Е.И. Забабахин. - Сне-жинск: РФЯЦ-ВНИИТФ, 1997.
5. Куропатенко, В.Ф. Основы численных методов механики сплошной среды / В.Ф. Куро-патенко, Е.С. Шестаковская. - Челябинск: Изд-во ЮУрГУ, 2017.
6. Куропатенко В.Ф. Локальная консервативность разностных схем для уравнений газовой динамики / В.Ф. Куропатенко // Журнал вычислительной математики и математической физики. - 1985. - Т. 25, № 8. - C. 1176-1188.
7. Белоцерковский, С.М. Метод крупных частиц в газовой динамике / С.М. Белоцерков-ский, Ю.М. Давыдов. - М.: Наука, 1982.
Елена Сергеевна Шестаковская, кандидат физико-математических наук, доцент, кафедра «Вычислительная механика:», Южно-Уральский государственный университет (г. Челябинск, Российская Федерация), [email protected].
Ярослав Евгеньевич Стариков, студент, кафедра «Вычислительная механика», Южно-Уральский государственный университет (г. Челябинск, Российская Федерация), [email protected].
Наталия Леонидовна Клиначева, кандидат физико-математических наук, доцент, кафедра «Вычислительная механика», Южно-Уральский государственный университет (г. Челябинск, Российская Федерация), [email protected].
Поступила в редакцию 26 марта 2021 г.
MSC 76M20 DOI: 10.14529/mmp210212
A METHOD FOR STUDYING THE DISSIPATIVE PROPERTIES
OF DIFFERENCE SCHEMES IN EULER COORDINATES
E.S. Shestakovskaya1, YA.E. Starikov1, N.L. Klinacheva1
1 South Ural State University, Chelyabinsk, Russian Federation
E-mails: [email protected], [email protected], [email protected]
At present, numerical methods for calculating shock-wave flows of liquid and gas in Eulerian coordinates have become widespread; therefore, the study of their characteristics is an urgent task. The paper presents an approach to assessing the dissipative properties of such difference schemes on strong discontinuities. The idea of the method is to construct the entropy production equation, the approximation error of which can be expressed by a combination of the approximation errors of the equations that make up the difference scheme. The equation of entropy production on a weak shock wave is used as a criterion for the dissipativity of the difference scheme. The paper evaluates the dissipative properties of the large particle method using the proposed method.
Keywords: difference schemes; Euler coordinates; dissipativity; entropy production.
References
1. Yanenko N.N., Shokin Y.I. The First Differential Approximation to Finite-Difference Schemes for Hyperbolic Systems of Equations. Siberian Mathematical Journal, 1969, vol. 10, no. 5, pp. 868-880. DOI: 10.1007/BF00971662
2. Rozhdestvenskii B.L., Yanenko N.N. Sistemy kvazilineinykh uravnenii [Systems of Quasilinear Equations]. Moscow, Nauka, 1968.
3. Zel'dovich YA.B., Raizer YU.P. Fizika udarnykh voln [Physics of Shock Waves]. Moscow, Fizmatgiz, 1963.
4. Zababakhin E.I. Nekotorye voprosy gazodinamiki vzryva [Some Questions of the Gas Dynamics of the Explosion]. Snezhinsk, RFYaTs-VNIITF, 1997.
5. Kuropatenko V.F., Shestakovskaya E.S. Osnovy chislennykh metodov mekhaniki sploshnoi sredy [Fundamentals of Numerical Methods of Continuum Mechanics]. Chelyabinsk, Publisher Center of South Ural State University, 2017.
6. Kuropatenko V.F. Local Conservativeness of Difference Schemes for the Equations of Gas Dynamics. Computational Mathematics and Mathematical Physics, 1985, vol. 25, no. 8, pp. 1176-1188. DOI: 10.1016/0041-5553(85)90157-0
7. Belotserkovskii S.M., Davydov YU.M. Metod krupnykh chastits v gazovoi dinamike [The Method of Large Particles in Gas Dynamics]. Moscow, Nauka, 1982.
Received March 26, 2021
