ЖУРНАЛ ВЫЧИСЛИТЕЛЬНОЙ МАТЕМАТИКИ И МАТЕМАТИЧЕСКОЙ ФИЗИКИ, 2007, том 47, < 3, с. 445-459
УДК 519.63
О РАЗНОСТНЫХ АППРОКСИМАЦИЯХ ПЕРЕОПРЕДЕЛЕННЫХ ГИПЕРБОЛИЧЕСКИХ УРАВНЕНИЙ КЛАССИЧЕСКОЙ
МАТЕМАТИЧЕСКОЙ ФИЗИКИ1}
© 2007 г. Д. П. Бабий*, С. К. Годунов*, В. Т. Жуков**, О. Б. Феодоритова**
(* 630090 Новосибирск, пр-т Коптюга, 4, Ин-т матем. СО РАН;
** 125047 Москва, Миусская пл., 4, Ин-т прикл. матем. РАН) e-mail: godunov@math.nsc.ru; zhukov@kiam.ru; feodor@kiam.ru Поступила в редакцию 17.10.2006 г.
Рассматриваются вопросы разностной аппроксимации переопределенных систем гиперболических уравнений. Для системы уравнений гидродинамики, магнитной гидродинамики, Максвелла и упругости приведены формулировки расширенных переопределенных систем. Обсуждаются подходы к построению разностных схем для таких систем. Библ. 9. Фиг. 5.
Ключевые слова: переопределенные гиперболические уравнения, разностные аппроксимации, уравнения газовой динамики, уравнения акустики, уравнения Максвелла.
1. ВВЕДЕНИЕ
Мы остановимся на трудностях, которые приходится преодолевать при разработках алгоритмов решения переопределенных систем уравнений, называемых гиперболическими. Вообще говоря, гиперболическими называют системы, состоящие из числа уравнений, совпадающего с числом неизвестных, и либо имеющими определенную структуру характеристического коноида (предложение И.Г. Петровского), либо записывающиеся с помощью коэффициентов, образующих симметрические (эрмитовы) матрицы (предложение Фридрихса). Одна из этих матриц должна быть положительно-определенной. Полезно отметить, что системы, гиперболические по Петровскому, совсем не обязательно допускают запись в виде симметрических гиперболических по Фридрихсу.
Пример, демонстрирующий это обстоятельство, был построен в 1984 году В.В. Ивановым, однако его нельзя считать общеизвестным, так как он опубликован только на русском языке в препринте (см. [1]).
Будем придерживаться определения Фридрихса и называть систему из N уравнений для N неизвестных функций, образующих N-мерный вектор u, гиперболической, если она записана в виде
.ди „ du г A -Т- + Bk = f, д t k д xk J
где квадратные матрицы A = A* > 0, Bk = B* либо постоянны, либо являются гладкими функциями от u (случай квазилинейных уравнений).
Простейшая система уравнений газовой динамики в лагранжевых координатах, разностную схему для которой один из авторов разрабатывал зимой 1953-1954 гг. (она была опубликована лишь в 1959 г., см. [2]), а именно
дulдt + dp/дm = 0, E = E(v, s), -dv/dt + дu/дm = 0, p = -Ev(v, s), (1.1)
дs Ш = 0,
^ Работа является расширенным изложением пленарного доклада на XI международной конференции по гиперболическим уравнениям (17 июля 2006 года, Франция, Лион) и выполнена при финансовой поддержке грантов 047.016.003 РФФИ-NWO (Нидерланды), "Ведущие научные школы" НШ-9019.2006.1 РФ, междисциплинарного проекта < 5 Президиума СО РАН, программы < 3.1 ОМН РАН.
после введения гамильтониана H(u,p, s) = E+pv + u2/2, dH = udu - vdp + Tds может быть записана в виде, предлагаемом Фридрихсом:
/ н ** uu H up Hus \ д д t / \ u ( 0 1 0 \ --д--д m / \ u
н pu Hpp Hps p + 1 0 0 p
н ( ** su Hsp Hss ( s 0 0 0 V ( s /
Нетрудно видеть, что эта запись при обозначении Н = ри (с учетом равенства И6 = 0) эквивалентна следующей системе:
и дНи + дНи = 0
д t дт
р дНр + дН = о, (1.2)
дt дт
Т д-Н + дН = 0. д t д т
Умножая последние уравнения на выписанные слева от них множители и,р, Ти суммируя, мы приходим к дополнительному равенству, совместному с системой (1.2):
д( Е + — У 2 ) + дрри = о
дt дт ' (1.3)
2
Е + и- = иНи + рНр + ^ - Н.
Таким образом, добавляя уравнение (1.3) к (1.1) или (1.2), получаем совместную переопределенную систему из четырех уравнений.
Интересно отметить, что из этой переопределенной системы симметрическую гиперболическую подсистему можно выделить иначе. Пусть эта подсистема составлена из уравнений
u
du + dp = 0
Т д t д т
р дV+ ди =о
-Т - Э7 дт = ' (1.4)
д( Е + и2) Л
1 У 2 ] + дри = о Т д t д т
слева от которых выписаны "интегрирующие множители". Линейная комбинация уравнений с этими множителями в виде коэффициентов оказывается "законом сохранения энтропии" дs/дt = 0.
Если воспользоваться термодинамическим тождеством
ds = - Т^и - (- V) + Т^Е + 7)'
d
u2л
E + p v - Ts--
F_2
T
T) - vi-T) + (E + |)d-1
и ввести обозначения
Т 42
= Л. Т
Г
Е + р V - Т - — Ь = --- = Ь(40, 4,, 42), Ь =
то уравнения (1.4) примут вид
Так как
Т
40
42
Ыа„ д ь4
4о
+
4о
= 0,
д t д т
д Ь4 д Ь4
4 + = о,
дt дт
дь4 дь4
42
+
42
д t д т
= 0.
(1.5)
4142
(1.6)
4о Ь40 + 4! Ь4, + 42 Ь42 Ь = 5,
42
4о Ь40 + 4! Ь4, + 42 Ь42 Ь = 0,
то равенство = 0 можно получить из уравнений (1.4) путем взятия их линейной комбинации с выписанными слева множителями 40, 41, 42. Квазилинейная форма выписанных уравнений (1.6)
Т д-к. Ь44 дt
Ь = 0 Т44 д т
при выпуклом производящем потенциале Т(40, 41, 42) свидетельствует об их гиперболичности.
Именно система (1.4), составленная из законов сохранения количества движения, объема и энергии, была положена в основу разработанной схемы, на решениях которой безусловно выполнялось требование неубывания энтропии. По этой причине схему можно было использовать при расчете разрывных решений - ударных волн. На ударных волнах энтропия растет и переопределенную систему, составленную из (1.1) и (1.3), уже нельзя считать совместной.
Размышления над описанными обстоятельствами подтолкнули к выделению класса "термодинамически согласованных" уравнений, образующих симметрические гиперболические системы, составленные из дивергенций (законов сохранения)
дь4 дМ]4
-г-4 + 4 = 0
д t д х
и совместных (на гладких решениях) с дополнительным законом сохранения из [3], [4]
д(4—4 - Ь) , д[4М{, -М1 ]
д t
д х1
= 0.
В дальнейшем, однако, выяснилось, что далеко не все уравнения классической математической физики к такому виду приводятся. Вот, например, аналогичная формализация уравнений магнитной гидродинамики (недивергентная):
д t
дЬ4о , д(икь)4о
дt
ди д( ЫкЪ )ъ
д Хк
= 0,
д п.
д хк
" - ^ эХ- = °,
(1.7)
4
2
0
4
Ыи, d(ukL )и
д к,
\ _ £ _- = 0
д г дхк Нк дхк
Хотя система (1.7) без труда переписывается в симметрическом гиперболическом виде, она не содержит законов сохранения количества движения и энергии. Эти законы (дивергентные равенства)
д Lu. д[ ( ukL) и , ~ к , LK] dt dxk
= 0,
д( Я0 Lq0 + h,Lki + U,Lu, - L ) , д[ Uk ( ^0 Lq0 + h,Lki + U,Lui) - U,k,Lkt ]
9o
t
dxt
= 0
справедливы не на всех решениях системы (1.7), а только на решениях, удовлетворяющих дополнительному уравнению
дЬк /д xi = 0,
совместному с (1.7), как это вытекает из ее следствия
Э_
t
J_
L q 0
дЦр
x
кдх.
J_
Lq 0
д^П
x
= 0.
Аналогично обстоит дело и ] натах.
уравнениях нелинейной теории упругости в эйлеровых коорди-
Отмеченное обстоятельство вызывает трудности при численном решении уравнений магнитной гидродинамики. Обзор этих сложностей и описание всевозможных технологических приемов их преодоления составляет значительную часть в монографии [5], с которой мы знакомились после того, как профессор Ван-Леер обратил наше внимание на затруднения, причиной которых является накопление вычислительных погрешностей. Эти погрешности, по-видимому, накапливаются в основном при долговременном расчете больших гладких полей, тогда как в работах, реферируемых в [5], при построении расчетных схем упор делается на их адаптации к расчету зон, содержащих ударные переходы. Размышляя об этом, мы решили, что стоит разобраться в проблеме расчета гладких решений у совместных переопределенных систем, на некоторое время забыв о том, что у этих систем существуют и разрывные решения. К вопросу о расчете разрывных решений, на наш взгляд, следует вернуться лишь после того, как ситуация с гладкими решениями станет совершенно ясной. Так же как и 53 года тому назад, когда один из авторов готовил свой вариант расчета разрывных решений газовой динамики к пуску первой русской серийной ЭВМ "Стрела", мы решили сначала детально исследовать задачи, описываемые линейными уравнениями с постоянными коэффициентами. Вот об этих исследованиях пойдет речь ниже. Хотя они пока находятся на начальном этапе, мы в процессе их проведения натолкнулись на новые, неожиданные для нас приемы построения вычислительных алгоритмов.
2. СХЕМЫ С КОРРЕКТИРОВКОЙ ДЛЯ УРАВНЕНИЙ АКУСТИКИ Начнем с разбора простейшего примера двумерных уравнений акустики:
дм , др , ,, д^ . др и,
д + дХ = а(х' *г}, э7 + дУ = ь (х'*г},
др ди д 7 . .
д + дХ + ду = с (х'*г},
0 < х <п, 0 < у <п, г > 0.
(2.1)
Ограничимся для простоты изучением решений, удовлетворяющих граничному условию р = 0 на границе области и заданным при г = 0 начальным условиям.
_д t дх ду_
(2.4)
Дополнительные уравнения, совместные с этой системой, имеют вид
д u2 + v2 + p2 , dpu , дpv , ,
-г---—— + + — = au + bv + cp,
д t 2 дх ду (2 2)
д (u v ) = дa(х, у, t) дЬ (х, у, t)
a? uy- Vx) дУ дХ •
Основная трудность, как известно, состоит в выборе разностной аппроксимации, обеспечивающей достаточно точное соблюдение уравнения (2.2). Расширим нашу систему, добавив в нее новые неизвестные потенциалы ю1, ю2, ю и новые уравнения
дю1 ( а дю2 и( дю ( А
— = a(х, у, t), = b (х, у, t), "ду = p(х, у, t), (2 3)
MjÍх, у, 0) = u(х, у, 0), Ю2(х, у, 0) = v(х, у, 0), ю(х, у, 0) = 0. Легко убедиться, что ю удовлетворяет волновому уравнению
д2 ю - д^ю - д!ю = д
дt2 дх2 ду2 д tL и что u^, у, t), у(х, у, t) представимы в виде
u (х, у, t) = ю j (х, у, t) - юх( х, у, t), v(х, у, t) = ю2 (х, у, t) - юу( х, у, t) .
Расчет одного временн0го шага t —»-1 + т разобьем на два этапа. На первом этапе используется аппроксимация симметрической гиперболической системы
д и дp , .. д v дp ,, ..
д + дх = a(х,у, у), э7 + ду = b (х,у, у),
дp , д^ . д v . , ír.
dt + д! + ду = c (х,у, у ), (2.5)
дю! дю2 дю
= a, = b, т— = p.
Э/1 Э/1 э^
Несложный пример такой достаточно точной аппроксимации, основанный на симметрической гиперболичности, будет приведен ниже. После того как первый этап закончен, на втором этапе производим
Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.