Том ХЫН
УЧЕНЫЕ ЗАПИСКИ ЦАГИ
2012
№ 4
УДК 533.6
МОДИФИКАЦИЯ ЧИСЛЕННОГО МЕТОДА ГОДУНОВА ДЛЯ УРАВНЕНИЙ ГАЗОДИНАМИКИ
Ю. М. ЛИПНИЦКИЙ, А. В. САФРОНОВ
Представлены модификации метода С. К. Годунова для уравнений газодинамики на основе новых приближенных решений задачи о распаде разрыва, обеспечивающих выполнение условия неубывания энтропии в численных расчетах. Преимуществом подходов является также отсутствие необходимости вычисления скорости звука на границах ячеек сетки, упрощающее применение предложенных схем для случая сложного уравнения состояния с переменными свойствами газа. Предложенные приближенные решения задачи распада разрыва могут быть применены для численного решения широкого круга задач аэрогазодинамики методом Годунова.
Ключевые слова: распад газодинамического разрыва, проблемы энтропии в численном решении гиперболических уравнений, схемная вязкость, TVD схема.
ВВЕДЕНИЕ
Условие неубывания энтропии наряду с выполнением интегральных законов сохранения является одним из необходимых условий применимости численных методов для расчета течений газа с разрывами. Решение энтропийной проблемы в численных расчетах впервые предложено Дж. фон Нейманом [1] путем введения искусственной вязкости, которая проявляется в градиентных зонах течения, «размазывая» разрывы.
При численном решении уравнений газодинамики методом Годунова эту роль играет схемная вязкость [2], причем метод Годунова 1-го порядка, основанный на точном итерационном решении задачи Римана о распаде газодинамического разрыва, обладает минимальной схемной вязкостью, обеспечивающей неубывание энтропии при сквозном расчете сложных течений с разрывами и волнами разрежения.
Схематически начальный разрыв газа с различными состояниями в левом и правом полупространстве распадается на три волны: левую волну, контактный разрыв и правую волну. Обозначим их скорости И^, ф и Ид соответственно. Левые и правые волны могут быть в зависимости от начального перепада давления как волнами разрежения, так и волнами сжатия (скачками). Точное решение задачи Римана единственно, зависит от автомодельной переменной, скорости волн удовлетворяют неравенству ф < И < Заданной функции потока (компонентам вектора потока) отвечают два решения: одно соответствует дозвуковому случаю, другое — сверхзвуковому. Поэтому можно представить, что функция
ЛИПНИЦКИЙ Юрий Михайлович
доктор технических наук, профессор, начальник центра теплообмена и аэродинамики ЦНИИмаш
САФРОНОВ Александр Викторович
кандидат физико-математических наук,
начальник отдела газодинамики старта ЦНИИмаш
потока уравнений газодинамики при переходе из одного состояния в другое в автомодельной задаче о распаде разрыва обладает свойством выпуклости.
Точное решение задачи распада разрыва приводит к трудоемким вычислениям, кроме того, в ряде физических случаев теоретическое решение этой задачи затруднительно. В этой связи широкое распространение получили методы на основе приближенного решения задачи Римана. К данному классу можно отнести известные схемы LxF (Lax-Friedrichs) [3], Русанова [4], HLL (Harten, Lax, Leer) [5], Roe [6], EO (Engquist — Osher) [7], HLLC (Toro, Spruce, Speares) [8].
Для выполнения энтропийного условия схемная вязкость методов на основе приближенного решения задачи Римана должна быть выше, чем метода на основе точного решения. Доказательство этого положения представляет проблему в связи с нелинейностью уравнений газодинамики.
Решением проблемы энтропии в численном решении гиперболических уравнений является также применение кинетического метода релаксации, например [9 — 13]. В этом варианте исходная нелинейная гиперболическая система уравнений преобразуется в линейную систему с релаксационным источником. Причем кинетическая интерпретация численных методов упрощает их энтропийный анализ.
Сравнительный анализ численных методов на основе решения задачи Римана и других подходов, применяемых для уравнений газодинамики можно найти, например, в работах [14, 15].
Интенсивное развитие методов Годунова связано с повышением порядка вычислений на основе процедур интерполяции параметров к граням ячейки сетки с ограничителями, обеспечивающими монотонность схемы. Первый такой ограничитель предложен Колганом [16], изложение теории и обзор вариантов можно найти, например, в работах [17 — 19]. После интерполяции к граням ячейки (реконструкции) параметры на границе ячейки в схемах повышенного порядка вычисляются из задачи распада разрыва аналогично схеме 1-го порядка. Поэтому повышение порядка не устраняет принципиальные недостатки и в первую очередь необходимо выполнение условия неубывания энтропии в базовых численных схемах 1-го порядка, которые анализируются в настоящей работе.
Отметим, что повышение порядка аппроксимации уменьшает различие между методами на основе точного и приближенных решений задачи Римана.
Из рассматриваемых методов доказанным является выполнение условия неубывания энтропии в газодинамических расчетах по схемам LxF, Годунова, Русанова, HLL, EO, а также HLLC в кинетической интерпретации [12, 13] и в массовых переменных [20].
Ниже в разделе 1 на примере скалярного закона сохранения с выпуклой функцией потока, с помощью подхода E-flux [21], проведен энтропийный анализ рассматриваемых схем и предложены модификации алгоритмов. В разделе 2 приведено обобщение подходов для численного решения уравнений газодинамики.
В качестве примера, на котором можно проиллюстрировать характер модификаций разностных схем, рассмотрим задачу Коши для скалярного закона сохранения:
с выпуклой функцией потока/(и). Для определенности примем /"(и) > 0 .
Уравнение (1) может иметь разрывные решения (скачки), скорость распространения которых 5 определяется соотношением Ренкина — Гюгонио: /(и+) - /(и_) = s(u+-и_). В случае гладких решений функция и(х, ¿) постоянна вдоль характеристик Жх/Ж = /' (и).
Энтропийное условие допустимости разрыва Олейник [22] заключается в том, что, например, в случае 5 > 0, характеристическая скорость за разрывом должна быть больше скорости распространения разрыва, а характеристическая скорость перед разрывом — меньше:
1. ЭНТРОПИЙНЫЙ АНАЛИЗ РАЗНОСТНЫХ СХЕМ
Щ + fx = 0
u (x,0) = u0(x), 0 < x < 1, 0 < t < T,
(1)
f'(u- ) > 5 > f'(u+ ).
(2)
Запишем разностную схему для уравнения (1) в консервативном виде:
Mf+1 = un — с( fi+1/2 - fi-1/2), (3)
At л ■
где с = — ; n — номер шага по времени с интервалом At; i — номер ячеики сетки по оси х c раз-Ax
биением Ax.
В методах типа Годунова поток на границе ячеек fi + 1/2 вычисляется из решения задачи распада разрыва с начальными данными, соответствующими параметрам в соседних ячейках сетки (i, i + 1). При этом схема устойчива при числе Куранта CFL = smaxc < 1, где smax — максимальная скорость распространения возмущений в расчетной области.
Как упомянуто, кроме выполнения законов сохранения для единственности численного решения необходимо выполнение энтропийного закона. Условие, обеспечивающее неубывание энтропии в ячейках сетки при численных расчетах, являющееся следствием теоремы Олейник, введено Ошером [21]. Это условие (E-flux), которое необходимо выполнять при аппроксимации потока на границе ячеек (fi+1j2) имеет следующий вид:
f+1/2 — f (u) < 0, Vu с [u¡,u;+1]. (4)
иг+1 — щ
Поток в методе Годунова на границе ячеек (1/,+°д,унов ), который вычисляется из точного решения задачи распада разрыва с начальными данными, соответствующими параметрам в соседних ячейках сетки, может быть представлен в виде [21]:
Я
i+1/2
Годунов =
uc[ui ,ui+1]
min f (u), u < u+1,
uc[u. ,u i+1 ]
' i+1 (5)
max (f (u)),u > u,+1.
Поток в методе Годунова удовлетворяет условию Е-/их, причем
/Е * /Годунов, щ < Щ+1,
/Е * /Годунов, Щ > иг+1. (6)
С целью сравнительного анализа диссипативных и энтропийных свойств численных методов разностное уравнение (3) представляется в следующем виде:
иГ1=ип _-/- /л+2 а+1/2(и+1- и,)-1 ег_1/2(иг - им), (7)
где Q —схемная вязкость
Й+1/2 =с
^ fi+1 + fi 2fi+1/2 ^
u,-,i — u,
(8)
Н+1
В этом случае аппроксимация потока на границе ячеек выражается следующим образом:
¿+1/2 = /Ш2+/- _ 2- е1+1/2(и,+1 _ и, ). (9)
2 2-
Схемная вязкость в методе Годунова имеет вид [21]:
Q
Годунов _
+1/2
= с max
uc[uj ,Uj+1
( f+1 + f - 2f (u) > ui+1 - ui
(10)
Энтропийному условию удовлетворяют методы, в которых при аппроксимации потока на границе ячеек схемная вязкость больше схемной вязкости метода Годунова [11]:
П > пГодун°в Q+1/2 г q+1/2 ,
(11)
что аналогично выполнению условия E-flux (4).
Согласно [11] рассматриваемые численные схемы монотонны и устойчивы при выполнении условия CFL и в случае
0+1/2 < Q+1/2 * QLF/2 - 1.
(12)
Здесь =Сfi+1/2 — схемная вязкость в методе Roe [7], Пил/2 — схемная вязкость
в методе Лакса LxF [3], принимающая максимальное значение.
Приведем сравнительный анализ величин схемной вязкости методов на основе приближенного и точного решений задачи распада разрыва. Ниже рассматриваются параметры на границе ячейки сетки. Индекс «L» относится к i-й ячейке сетки, индекс <R» — к i + 1-й, индекс «i + 1/2» — к параметрам на границе этих ячеек.
Схемная вязкость в методе EO [7, 21] имеет вид:
uR
Л f '(u)| du
qeo Qi+1/2
= с
uL
u p u :
(13)
В работе [21] показано, что поток EO удовлетворяет условию E-flux.
«Двухволновые» схемы HLL, Русанова и LxF представляются из соотношений Ренкина — Гюгонио, записанных для двух скачков (рис. 1):
(f - fL ) = sL (u*~ uL X (f - fR ) = sR (u - uR ).
(14)
Здесь SL, — скорости левой и правой волн соответственно. Причем для рассматриваемой выпуклой функции потока < ¿к.
Аппроксимация потока на границе ячеек из соотношений (14) имеет вид:
f = f * = fLsR - fRsL + sLsR (uR - uL)
fi+1/^ — J -
при этом
i+1/2
= u =
sR - sL
uRsR uLsL fR + fL
sd st
(15)
(16)
Для определения скоростей волн используются следующие варианты:
•л =-$ь = 1/о, (17)
^ = шах(| /|,|/ |), (18)
зь = шт(0М(/[, /к)), ^ = шах(0,8ир/, /к )).(19)
В случае (17) поток (15) представляет схему ЬхГ.
Рис. 1. Геометрия энтропийного E-flux анализа
В этой схеме используется максимальная оценка скоростей волн по всей расчетной области. Заметим, что в этом случае поток f (u*) из (16) отвечает схеме 2-го порядка LxW [23].
При вычислении скоростей волн согласно (18) поток (15) соответствует схеме Русанова, в которой применяется локальный максимум скоростей волн. При выборе варианта (19) аппроксимации потока на границе ячеек (14) соответствуют схеме HLL.
Рассмотрим случай смены знака производной функции f (u) при ur > ul. Этот случай является аналогом волны разрежения и характерен для исследования энтропийных проблем. Геометрическая интерпретация потоков рассматриваемых схем представлена на рис. 1.
В соответствии с (5), (6) поток на границе ячеек схемы Годунова соответствует минимуму функции f (u) на рассматриваемом интервале, а схемы, удовлетворяющие условию E-flux, находятся ниже этого уровня, как показано на рис. 1. В данном случае поток EO (13) совпадает с потоком Годунова. Поток HLL находится в зоне E-flux, на пересечении касательных к функции f (u), проходящих через точки (ul, fL) и (ur, JR). Потоки Русанова и LxF находятся также в зоне E-flux ниже потока HLL, поскольку линии, проведенные из точек (ul, fL) и (ur, JR) в этих схемах ближе к вертикали, чем в методе HLL. Очевидно, выбор скоростей волн (17) — (19) обеспечивает выполнение условия неубывания энтропии в двухволновых схемах LxF, Русанова и HLL, для любой выпуклой функции f(u).
В случае изменения знака характеристик (fL < 0, fR > 0), поток Roe, как показано на рис. 1,
не удовлетворяет энтропийному условию. На рис. 1 показан поток f (u*), которой при определении u по формуле (16) с выбором скоростей волн в виде (17) представляет собой схему LxW [23], при этом поток f (u*) больше потока схемы Годунова, но меньше потока Roe. В рассматриваемом случае чем меньше поток на границе ячеек, тем больше схемная вязкость (8), поэтому при смене знака характеристик схема LxW 2-го порядка предпочтительнее схемы Roe.
Согласно [11] и в соответствии с данными рис. 1 известные потоки по величине схемной вязкости располагаются следующим образом:
s^Roe ^ .^Годунов ^ глЕО ^ s^HLL ^ ^Русанов ^ ,^¡LxF _ 1 Qi+1/2 < Qi+1/2 < Qi+1/2 < Qi+1/2 < Qi+1/2 < Qi+1/2 = 1
На рис. 1 показан способ учета контактного разрыва в задаче распада разрыва с выполнением условия E-flux. Алгоритм заключается в следующем: по параметрам u , полученным с помощью формулы (16), находится наклон касательной в этой точке s* = f' (u*), затем вычисляется
поток f(u ), и далее определяется искомое решение, находящееся на пересечении данной касательной линии с линиями, обозначающими левую и правую волны. Потоки, отвечающие этому
обозначены f и Jr , а переменные ul и ur соответственно. Для этого случая можно записать следующие соотношения на разрывах:
(fL — fL ) = SL (uL — uL X (fR — L ) = S*(uR — uL X (fR — fR ) = SR (uR — uR ). (20)
При этом в силу выпуклости функции f(u) выполняются неравенства
Sl < s* = f' (u*) < Sr . (21)
Выбор аппроксимации на границе ячеек зависит от знака s* :
fi+1/2 = fL при s* > 0 и¿+1/2 = fR при s* < 0.
Решение задачи на основе соотношений (20) можно представить в виде:
fi+1/2 = f* [f (u*) — f *], (22)
-Sr
У S* S l SR S * j
(23)
На рис. 1 данные потоки обозначены XContact, X — означает способ выбора скоростей волн sl и Sr, соответствующий схемам LxF, Русанова, HLL (17) — (19). Из рис. 1 видно, что одно из
пересечений касательной к функции f (u) в точке (u*, f (u*)) с линиями, обозначающими левую и правую волны, при выпуклой функции f (u) в случае f < 0, fR > 0 всегда будет лежать ниже минимума этой функции. Поэтому потоки XContact удовлетворяют энтропийному условию E-flux.
Из анализа соотношений (22),(23) и данных рис. 1 выявлена следующая возможность модернизации двухволновых схем. Если в качестве скорости средней волны s* принять предельные
значения, т. е. из точки (u*, f(u*)) провести линии, параллельные левой и правой волне, то получим следующее значение параметра ю в уравнении (22):
ю=ю
GFORCE
= min
-Sr
\ sr s
L JR
Sd - S
(24)
L
™=™contact= min
При выборе sr = -sl = smax, полученная схема (15), (16), (22), (24) соответствует схеме GFORCE [24], поэтому в соотношении (24) параметр ю обозначен &gforce.
Новые схемы (15), (16), (22), (24) с выбором скоростей волн (18) или (19) являются обобщением схемы GFORCE. На рис. 1 положение данных потоков обозначено XGFORCE. Из геометрии видно, что полученные схемы удовлетворяют условию E-flux. Потоки XGFORCE, основанные на максимальной оценке скорости средней волны (рис. 1), которая всегда выше скорости контактного разрыва, «размазывают» контактный разрыв, но их схемная вязкость (8) меньше, чем соответствующих потоков при ю = 0 (LxF, Русанов, HLL).
Исходя из вышеизложенного, предложенные модификации численных методов (XContact и XGFORCE) по величине схемной вязкости можно расположить следующим образом:
ÓRoe ^ ^Годунов ^ r>EO ^ ^XContact ^ ^XGFORCE ^ s^HLL ^ ^Русанов ^ s^LxF _ л i+1/2 < Qi+1/2 < Qi+1/2 < Qi+1/2 < Qi+1/2 < Qi+1/2 < Qi+1/2 < Qi+1/2 = 1
В соответствии с (12) схемы XContact и XGFORCE устойчивы и монотонны.
В общем случае при выборе параметра ю в диапазоне 0 <ю <ю contact схемы (15), (16), (22) с выбором скоростей (17) — (19) удовлетворяет условию E-flux, поэтому возможны и иные варианты, с различными переключателями.
2. УРАВНЕНИЯ ГАЗОДИНАМИКИ
Рассмотрим двумерный случай в расщеплении по оси х. Уравнения газодинамики представляются в следующем виде:
и + Ех = 0, (25)
где и — вектор консервативных переменных и Е — вектор потока:
р pw
U = pw , F = pw2 + P
pv pwv
PE _ pEw + Pw
Здесь р — плотность газа; V, V — нормальная и касательная к грани ячейки компоненты скорости; Е = е + м> 2/2 + V 2/2 — полная энергия на единицу объема; е = р/р/(у - 1) — внутренняя энергия; р — давление; у — показатель адиабаты.
Методы, описанные для скалярного закона сохранения в разделе 1 (3), (15), (16), (22), (23), (24), обобщаются для системы (25), путем замены величины и на вектор консервативных переменных и, величины/— на вектор потока Р, а значения скоростей волн SL, я* и sя соответственно на Жь, Ж* и Жя:
тП+1
где
иг1 = ип _-(Р,+1/2 _Р-1/2), Р* = № _ + ЖЩя (ик _ иь)
Жя _ Жь ияЖя _ иьЖь _ Ря + рг
и* = -я"я
Жя _ Жь
Р,+1/2 = Р [Р (и ) _ Р ],
_ЖЬ
Жя
Ж* _ Жь Жя _ Ж *
(26)
(27)
(28)
(29)
или
ю = юОРОяСЕ = т1П I
_ЖГ
Жя
Жя _ Жь Жя _ Жь ,
(30)
где
Поток Р (и *) в формуле (28) вычисляется следующим образом:
Р (и *) = Ж,и * + 0*,
е* = {0,/,0,/Ж*}7 .
(31)
Скорость (Ж*) и давление (р*) на контактном разрыве определяются из соотношений на разрывах для уравнений газодинамики в лагранжевых переменных [20] с введением массовых скоростей через поверхности левой ть = рь(^ь - Жь) и правой тя = ря(Жя - ^я) волн:
Ж* = = wя =
* шя^я + т1Уь _Ря + Рь
ть + Шя
Р = Рь = Ря =
тьРя + тяРь _ тьтя(^я _)
т
т
(32)
(33)
Уравнения в универсальном виде (28), полученные из соотношений на разрывах, в зависимости от выбора скоростей волн Жя, Жь и способа вычисления параметра ю включают в себя девять численных методов: классические схемы ьхР, Ньь, Русанова (ю = 0), их варианты с учетом контактного разрыва Х-СоШаМ (29) и модификаций Х-ОРОяСЕ (30), где «X» обозначает соответствующий способ выбора скоростей волн.
Выбор скоростей волн:
1) Жя = _Жь = 1/ - , соответствует схеме ьхР;
2) Жя = _Жь = max(abs(wя) + ся, abs(wL) + сь), соответствует схеме Русанова;
ю = Юсоп(ас1 = т1П
3) Ж = тт(м>ь - е^, - ер), Жр = шах(^ + е^,, ^^ + ер), соответствует схеме ИЬЬ.
Следует подчеркнуть, что все перечисленные подходы энтропийно обоснованы.
Заметим, что при аппроксимации потока на границе ячеек сетки по формулам (26) — (33) не требуется привлечения уравнения состояния. Это упрощает применение подходов для сложного уравнения состояния [25, 26] и расчетов многокомпонентных течений с переменными теплофизическими свойствами [27].
2 4 6 8 10 12 14 16 18
Рис. 2. Картина заполненных контуров плотности (результаты расчетов по схеме ЯЬЬ-СоШае^
2 1.6 1.2 0.8 0.4
Рре
х /Ra
rms -1 -
-2
-3
-4
-5
-6
12
16
20
RK3, KorerS
----HLL-Cortact
Годунов — ■ — HLL-GFORCE
■ Л .
Рис. 3. Изменения давления на оси струи
III liter
4000 8000 12000 16000 20000
Рис. 4. Сходимость расчетов
Обобщение изложенных методов на трехмерный пространственный случай аналогично методу Годунова [2].
В качестве иллюстраций на рис. 2 — 4 даны результаты численных расчетов сверхзвуковой осесимметричной струи воздуха, истекающей в затопленное пространство. Уравнения Эйлера решались методом установления. Параметры на срезе сопла: число Маха Ма = 2, отношение давления на срезе сопла к внешнему давлению ра/ре = 1.5, температура торможения Т00 = 300 К. Струя истекает справа налево, срез сопла расположен от нуля до 1 по вертикальной оси. На срезе сопла задавался сверхзвуковой поток, на оси струи — условия симметрии, на гранях ячеек, примыкающих к свободным границам, потоки задавались из решения задачи Римана с параметрами, соответствующими состояниям в примыкающей ячейке и невозмущенной внешней среде. Расчетная область занимала 3 х 20 радиусов выходного сечения сопла (Ра = 1). Использовалась квадратная сетка с количеством разбиений на радиус среза сопла Та = 20.
Расчеты проводились по ТУБ схеме Рунге — Кутта 3-го порядка по времени с ограничителем Корена 3-го порядка по пространству (РК3-Когеп) [19] с решением задачи распада разрыва по схемам Годунова: с точным решением задачи Римана, а также представленными ЯЬЬ-СоМаМ и ЯЬЬ-ОЕОРСЕ.
В потоке образуются скачки, волны разрежения и контактный разрыв (рис. 2). Разрывы проходят поперек ячеек сетки и рассчитываются сквозным образом (без выделения).
0
4
8
0
На рис. 3 приведены изменения давления на оси струи, обозначения кривых соответствуют рис. 4.
На рис. 4 показана сходимость по времени в виде изменения логарифма максимального (в расчетной области) приращения плотности за шаг счета в зависимости от количества шагов счета по времени (итераций).
Из рис. 3, 4 видно, что результаты расчетов с применением предложенных способов решения задачи распада разрыва близки к схеме Годунова.
ЗАКЛЮЧЕНИЕ
Приведены модификации численных методов решения уравнений газодинамики на основе приближенного решения задачи Римана из соотношений на разрывах.
Предложены новые способы учета контактного разрыва и алгоритм с максимальной скоростью средней волны в методах с верхней оценкой скоростей левой и правой волн в задаче распада разрыва.
На примере скалярного закона сохранения с выпуклой функцией потока, а также на основе кинетической интерпретации для уравнений газодинамики дано энтропийное обоснование предложенных модификаций.
Результаты работы могут быть применены для решения широкого круга задач аэрогазодинамики.
Работа выполнена при поддержке грантов РФФИ № 10-08-00501, № 10-01-00677.
ЛИТЕРАТУРА
1. Neumann J., Richtmyer R. A method for the numerical calculation of hydrody-namic shocks // J. Appl. Phys. 1950, 21, № 3, р. 232 — 237.
2. Численное решение многомерных задач газовой динамики / Под ред. С. К. Годунова. — М.: Наука, 1976.
3. Lax P. D. Weak solutions of nonlinear hyperbolic equations and their numerical computation // Comm. Pure Appl. Math. 7, 159 (1954).
4. Русанов В. В. Расчет взаимодействия нестационарных ударных волн с препятствиями // ЖВМ и МФ.1961, 1, N 2, р. 267 — 279.
5. Harten A., Lax P. D., B. Van Leer. On upstream diffrencing and Godunov-type schemes for hyperbolic conservation laws // SIAM Review. 1981. 25, N 1, р. 35 — 61.
6. R o e P. L. Approximate Riemann solvers, parameter vectors, and difference schemes // J. Comput. Phis. 1981. 43, N 2, р. 357 — 372.
7. Enguist B.,Osher S. One-sided difference approximation for nonlinear conservation laws // Math. Comp. 36(1981), р. 321 — 351.
8. Toro E. F., Spruce M., Speares S. Restoration of the contact surface in the HLL Riemann solver // Shock Saves. 1994, 4, р. 25 — 34.
9. Jin S., X i n Z. P. The relaxation schemes for systems of conservation laws in arbitrary spacedimensions // Comm. Pure Appl. Math. 1995. 48 : 235 — 276.
10. L e V e q u e R. J., P e l a n t i M. A class of approximate riemann solvers and their relation to relaxation schemes // J. Comput. Phys. 2001, 172, р. 573 — 591.
11. Tadmor E. Entropy stability theory for difference approximations of nonlinear conservation laws and related time-dependent problems // Acta Numerica (2003), p. 451 — 512. Cambridge University Press, 2003.
12. Bouchut F. Entropy satisfying flux vector splittings and kinetic BGK models// Nu-mer. Math. 2003, 94, р. 623 — 672.
13. Сафронов А. В. Кинетические схемы для уравнений газодинамики // Вычислительные методы и программирование. 2009. Т. 10, с. 62 — 74.
14. Toro E. F. Riemann solvers and numerical methods for fluid dynamics. — Berlin. Heidelberg: Springer-Verlag, 2009.
15. Воронич И. В. Сравнительный анализ группы численных методов газовой динамики. — М.: МФТИ, 2007, 151 с.
16. Колган В. П. Применение принципа минимальных производных к построению конечно-разностных схем для расчета разрывных решений газовой динамики // Ученые записки ЦАГИ. 1972, Т. III, № 6, с. 68 — 77.
17. Bram van Leer. Upwind and high-resolution methods for compressible flow: From Donor Cell to Residual-Distribution Schemes // Commun. Comput. Phys. 2006. V. 1, N 2, p. 192 — 206.
18. B e r g e r M., A ft o s m i s M. J. Analysis of slope limiters on irregular grids // AIAA Paper 2005-0490. 2005.
19. Сафронов А. В. Оценка точности и сравнительный анализ разностных схем сквозного счета повышенного порядка // Вычислительные методы и программирование. 2010, 11.
20. Сафронов А. В. Разностный метод для уравнений газодинамики из соотношений на разрывах // Математическое моделирование. 2008. 20, № 2, с. 76 — 84.
21. O s h e r S. Riemann solvers, the entropy condition, and difference approximations // SIAM J. Numer. Anal. 1984; 21(2) : 217 — 235.
22. Олейник О. А. О единственности и устойчивости обобщенного решения задачи Коши для квазилинейного уравнения // Успехи мат. наук. 1959. 14, № 2 (86), с. 165 — 170.
23. Lax P., Wendro ff B. Systems of conservation laws // Comm. Pure Appl. Math. 13 (1960), р. 217 — 237.
24. Toro E. F., Titarev V. A. MUSTA fluxes for systems of conservation laws // J. of Comp. Phys. 216 (2006), р. 403 — 429.
25. Сафронов А. В. Разностный метод решения уравнений гидродинамики с двучленным уравнением состояния на основе энтропийно-согласованного приближенного решения задачи Римана. Высокопроизводительные вычисления в задачах механики и физики: Сб. научных трудов, посвященный 75-летию со дня рождения А. В. Забродина. — М.: ИПМ им. Келдыша, 2009, с. 143 — 149.
26. Липницкий Ю. М. , Сафронов А. В. Численные методы Годунова типа для уравнений газодинамики (Риман — солверы). — Тезисы докладов // Третья всероссийская конференция «Вычислительный эксперимент в аэроакустике». — Светлогорск, Калининградской области. 2010.
27. Сафронов А. В. Численный метод расчета струй продуктов сгорания при старте ракет // Космонавтика и ракетостроение. 2007. № 1(46), с. 72 — 79.
Рукопись поступила 9/III2011 г.