Об одном методе сквозного счета ударных волн Текст научной статьи по специальности «Математика»

УДК 519.63

БОТ: 10.14529/шшр140106

ОБ ОДНОМ МЕТОДЕ СКВОЗНОГО СЧЕТА УДАРНЫХ ВОЛН

В.Ф. Куропатенко

Сильные разрывы - ударные волны возникают в сплошной среде при динамических внешних воздействиях. На поверхности сильных разрывов законы сохранения принимают вид нелинейных алгебраических уравнений, связывающих скачки величин по обе стороны разрыва. На сильном разрыве энтропия терпит скачок. В этом заключается принципиальное различие между ударными волнами и волнами с непрерывным изменением величин. В однородных разностных методах сильный разрыв заменяется слоем конечной ширины, сравнимой с размером сеточной ячейки. Такое свойство разностных схем получило название дистракции. Поскольку состояние за разрывом связано ударной адиабатой с состоянием перед разрывом, то в области дистракции сильного разрыва должен действовать механизм, обеспечивающий возрастание энтропии. Физическая вязкость и теплопроводность в уравнениях механики сплошной среды не устраняют необходимости введения поверхности сильного разрыва и, следовательно, не могут обеспечить величину дистракции, сравнимую, с несколькими ячейками разностной сетки. В работе рассмотрены несколько разностных схем, в которых диссипация энергии в слое дистракции определяется уравнениями, справедливыми на поверхности сильного разрыва.

Ключевые слова: ударная волна; разностный метод; дистракция; диссипация энергии; законы сохранения.

1. Идея метода

Ударная волна - поверхность сильного разрыва заменяется слоем конечной ширины, содержащим несколько ячеек сетки. Параметры вещества, находящегося в одной сеточной ячейке, за несколько шагов по времени изменяются от состояния перед сильным разрывом до состояния за сильным разрывом. Эти состояния в случае идеальной среды связаны законами сохранения массы, количества движения и энергии в виде системы нелинейных алгебраических уравнений

р - р0 - W(и - ио) = 0, (1)

и - ио + W(VI - Уо) = 0, (2)

Ег + 1 и2 - Ео - 1 ио2) W - (Ргиг - РоЩ) = 0, (3)

где Р - давление, У - удельный объем, и - скорость, Е - удельная внутренняя энергия,

W = —— - скорость сильного разрыва. Состояние вещества перед разрывом, отмеченное

аЬ

индексом «0>, считается заданным, а состояние за разрывом с индексом «1> - текущим. К уравнениям (1)-(3) добавляется уравнение состояния

Рг = Р(Уг, Ег). (4)

Система четырех уравнений (1)-(4) содержит пять величин. Таким образом, чтобы определить конкретное состояние за разрывом, нужно задать одну из текущих величин, т.е. выбрать ее в качестве параметра.

В разностных методах сквозного счета ударных волн сильный разрыв заменяется ударным слоем шириной в несколько ячеек сетки [1—6]. Происходит дистракция сильного разрыва [7]. Уравнения (1)—(3) связывают состояния вещества на границах этого слоя. Внутри же слоя действуют различные механизмы диссипации энергии [2-5]. Рассмотрим разностный метод счета ударных волн, в котором сильный разрыв заменяется пакетом следующих друг за другом ударных волн меньшей амплитуды. Таким образом внутри ударного слоя в каждый момент времени Ьп в каждой ячейке сетки создается сеточная ударная волна, которая после перехода на следующий шаг по времени ¿п+1 = ¿п + т исчезает. Сеточная ударная волна в ячейке сетки с индексом г + 0, 5 схематически изображена на рисунке.

t 11

Схематическое изображение «:сеточной> ударной волны АВ и состояний перед ней

0 и за ней 1

Чтобы не усложнять понимания идеи метода, рассмотрим разностную схему в простейшей постановке. Одномерное течение идеальной среды с плоской симметрией в лагранжевых координатах и без теплопроводности описывается системой законов сохранения

дУ ди ди дР дє дРи _ 0

дЬ дт ’ дЬ + дт ’ дЬ + дт '

Умножим каждое из уравнений (5) на йтйЬ и проинтегрируем по площади сеточной ячей-

ки. Применив к полученным интегралам теорему о среднем значении, получим систему разностных уравнений

У+0,5 — Уі+0,5 и;*+1 — иі 0 (6)

-------Т--------------— _ 0 (6)

ип+1 — ТТп Р* Р *

иі+0,5 иі+0,5 , Рі+1 — Рі _ 0 (7)

т + Н _0 (7)

ЄП+01Б — єП+0,5 , (Ри)*+і — (Ри)* „ (8)

Т + Н _0 (8)

где величины с индексами п, п + 1 и нижним индексом і + 0, 5 являются средними на промежутке ті < т < ті+і в моменты времени Ьп и Ьп+1, величины с верхним индексом * и нижними индексами і и і + 1 являются средними на промежутке Ьп < Ь < Ьп+1 при значениях ті и ті+1. Иными словами ¡і _ / (ті) , /* _ / (Ь*) , /п _ / (Ьп). Масса сеточного интервала Н связана с координатами его границ Н _ (жп+1 — Х?) /Уі+о 5 • Величина Н от времени не зависит и, таким образом, сохраняется при переходе от одного момента времени к другому. Удельная полная энергия є есть сумма удельной внутренней и удельной кинетической 1 - 2

энергии є _ Е + ^ и- ■

Вообще говоря, уравнения (6)—(8) являются точными до тех пор, пока не конкретизированы координаты точек, в которых определены величины, входящие в эти уравнения. Уравнения (6)—(8) являются также общими до тех пор, пока не указаны уравнения для определения вспомогательных величин и*, Р*, (Ри)*. Эти величины называются вспомогательными величинами, т.к. после завершения перехода к решению в момент Ьп+1 они забываются. Поскольку при нахождении решения в момент Ьп+1 используются законы сохранения на сильном разрыве в виде (1)—(3), то тем самым на каждом временном шаге диссипация энергии в отличие от [2-4] определяется единственным физически обоснованным механизмом - законами сохранения на поверхности сильного разрыва.

Для пояснения идеи рассмотрим процедуру перехода с момента Ьп, где все сеточные функции известны, к моменту Ьп+1 в одной сеточной ячейке, изображенной на рисунке. Величина Ш может быть и положительной и отрицательной. Не умаляя общности, будем рассматривать случай Ш > 0.

Термодинамические величины Р ¿+о 5 , Р п+0 5 , Е п+о 5 , С п+о 5 характеризуют состояние вещества в момент Ьп в ячейке с номером і + 0, 5. Поскольку они являются результатом применения теоремы о среднем значении, то будем считать их постоянными на промежутке ті < т < ті+1. Это предположение означает, что на границах сеточных интервалов возникли разрывы. Будем рассматривать эти разрывы как сеточные ударные волны. Состояние с одной стороны разрыва - это состояние перед разрывом, а одна из величин в соседнем интервале - это величина за разрывом. Условие, что разрыв является сеточной ударной волной, имеет вид

ип — ип

иі+0 ,5 иі-0 , 5 < 0

ті+о , 5 — ті-0 , 5

Сравнение давлений Р+о 5 и Р'(—о 5 позволяет определить знак Ші. Если Рп—о 5 > Рп+о 5, то Ші > 0 и ударная волна будет при Ь > Ьп распространяться в интервале с номером і + 0, 5 (как показано на рисунке). В качестве величин перед разрывом берутся сеточные значения в момент Ьп

и0 _ иі+0,5, Р0 _ Рі+0,5, р0 _ Рп+0,5, Е0 _ Е!+0,5, У0 _ 1/р0,

а в качестве величин за разрывом берутся либо скорость Ц—о 5, либо давление р—о 5 ■

Далее решается система уравнений (1)-(4), в результате чего определяются величины и1,Р1,У1,Е1 и Ш и, соответственно, вспомогательные значения и* и Р*

и* _ Ц—0,5, Р* “ Р+0,5 + Ші (иі—0,5 — Щщ) , если задано Ь\,

(Р"-0.5 — Р+0,5) (9)

Р* _ Рі—0,5, и* _ иі+0,5 + ------------------------------------Ш-, еслизадано Р1.

Аналогично определяются вспомогательные значения и*+1, Р*+1

иі*+1 _ иі+0,5, р*+1 _ Рі+1,5 + Ші+1 (иі+0,5 — С/і+1,^ , если задано и1 и Ші +1 > 0

(р^+0,5 — Р+1,5)

р*+1 _ Рі+0,5, иі*+1 _ и¿+1,5 + ------Ш+1------, если задано Р1 и Ші+1 > °-

Покажем теперь, что уравнения (6)-(8) со вспомогательными значениями и *, Р * строго соответствуют мгновенным законам сохранения. Поскольку ячейка сетки состоит из двух частей, разделенных траекторией сеточного сильного разрыва - линией АВ, то в момент Ьп+1 масса вещества, находящегося за разрывом, равна Шт, масса вещества перед разрывом Н—Шт. Термодинамические величины, средние в массе Н в момент Ьп+1 находятся с помощью

мгновенных законов сохранения массы, количества движения и энергии. В рассматриваемом

тШ

случае усреднение идет по двум массам с массовыми концентрациями -------- для величин за

п тШ

разрывом и 1---------для величин перед разрывом.

п

Рассмотрим получение удельного объема У++5. Первый из мгновенных законов сохранения имеет вид

УПЙ = + (1 - тр) ^+1- (10)

Значение У™+1 это постоянное в интервале Шт значение удельного объема за фронтом сеточной ударной волны. Величины Уо, Цз, ^1, входящие в уравнение (2) совпадают с сеточными величинами Уо = У+0,5, Ц = Ц*, Ц = Ц+0,5

УГ+1 = У+0,5 - Ш (и* - Ц+0,5) .

Значение У0га+1 определяются из разностного уравнения

т

ут+1 _ ^0,5 + (ц*+1 — Ц+«) ■

Подставив У1п+1 и Уо"+1 в уравнение (10), получим разностное уравнение (6).

Рассмотрим теперь мгновенный закон сохранения количества движения

иЙад _ Тщ-и"* + (1 — ТШ)и°“+1- (11)

Значение ип+1 - это постоянное в интервале тШ значение скорости за фронтом сеточной ударной волны. Оно связано со значением ЦЩ _ Ц+о 5 уравнением на разрыве (1) в виде

иг+1 _ Ц+0,5 + Ш (Р1 — Рі+0,5) ■ (12)

Поскольку вспомогательное давление Р* в силу (9) совпадает с давлением за фронтом сеточной ударной волны Р1 _ Рі*, то уравнение (12) принимает вид

ггп+1 _ ттп + 1 ( р* рп )

Ц1 _ Ці+0,5 + Ш \Рі — Рі+0,5,1 ■

В интервале перед сеточным разрывом Р*+1 _ Р^+о 5. Поэтому за промежуток времени т скорость вещества в интервале Н — Шт изменится в соответствии с разностным уравнением

Цп+1 _ Ц+0 5 —

'п (т>* туп \

і+0,5 — Н — Шт ^Р+1 — Р+0,^ ■

Подставив Цп+1 и ЦЦ+1 в уравнение (11), получим разностное уравнение (7).

Наконец проделаем аналогичную процедуру с мгновенным законом сохранения энергии

4Й5 _ ІШ^Т+1 + (1 — тНШ) 4+ (13)

Значение удельной полной энергии єп+1 является средним в интервале за фронтом сеточной ударной волны, значение ¿п+1 - среднее в интервале перед фронтом сеточной ударной волны.

Значения єп+1 и ¿п+1 выражаются через основные и вспомогательные значения сеточных функций с помощью уравнений (3) и (8) в виде

^ _ ^+0,5 + Ш (Рг*иг* — Р+0,5иГ+0,5) , (14)

Єп + 1 _ Єі+0,5 — Н — Шт (Рі*+1Ці*+1 — Р+0,5иі+0,5) ■ (15)

Подставим (14) и (15) в (13). В результате получим разностное уравнение (8).

Подчеркнем, что при получении разностных уравнений (6)-(8) было сделано несколько предположений:

- В слое «:размазанной> ударной волны все функции кусочно постоянны.

- Разрывы на границах сеточных интервалов являются сильными разрывами (сеточными ударными волнами).

- Функции перед сеточным разрывом постоянны. Они выбираются в момент Ьп в рассматриваемом интервале.

- Функции Ц или Р за сеточным разрывом постоянны. Одна из них выбирается в момент Ьп в соседнем интервале, а вторая рассчитывается из уравнений (1)-(3).

- Решение в момент Ьп+1 получается путем применения мгновенных законов сохранения.

2. Дивергентная разностная схема

Рассмотрим с небольшими изменениями разностную схему из [9]. Все термодинамические величины и скорости определены в серединах сеточных интервалов, узлы сетки имеют координаты Ьп, Шг. Разностные уравнения имеют вид (6)—(8). Вспомогательные величины Р*, и* определяются с помощью двух различных алгоритмов в зависимости от того, разрежение или сжатие происходит на вспомогательном промежутке Шг-о,5 < ш < Шг+о,5.

Если внутри вспомогательной ячейки ЦП+о 5 — иП—о 5 — 0, то решение в указанном сеточном интервале является непрерывным, и Р*, и* определяются разностными уравнениями в виде

т та2

77* — ТТп (туп рп Л р* _ рп ' ^ ¡тт-п ттп \

иг = иг — 2п \Рг+о,5 — Рг—о,^ , Рг = Рг *+о>5 — иг—°>5/ '

Значения ип и рп находятся интерполяциями по и?+о5, и—о 5 и Ргп+о5, Рп—о5.

Вспомогательные величины Р*, и* используются только в уравнениях (6) и (7) для нахождения У+О.5. Ц++о15- Вместо уравнения энергии (8) на волне разрежения используется следствие из законов сохранения в виде

уп + 1

Еп+1 — Еп + І Р (У, Е) (1У _ 0. (16)

V п

Интегрирование вдоль изэнтропы успешно применялось в [10-14] для определения параметров вещества на волнах разрежения. Этот метод обеспечивает любую наперед заданную точность определения энтропии и устраняет ложную диссипацию энергии. Уравнение (16) может быть решено разными способами. Один из возможных способов интегрирования вдоль изэнтропы основан на использовании структуры УРС вещества [15]

р = Рх (У) + Рт (У, Ет), Е = Ех (У) + Ет'

Поскольку зависимости Рх (У) и Ех (У) заданы, то с их помощью находятся значения

Рхп+1 (Уп+^

и Еп+1 (Уп+1) ' В качестве зависимости между тепловым давлением Рт , тепловой энергией Ет и удельным объемом У возьмем уравнение, являющееся определением Рт

_ , дЕт

Рт — —

дУ / s

Согласно [15] зависимость Рт от У и Ет чаще всего представляется в виде

Рт _ Г (У) Ет/У,

где

Г(У)_ — ,

у ’ й 1п У,

а функция в (У) есть аналог характеристической температуры Дебая [15]. Из этих трех уравнений после интегрирования вдоль изэнтропы получается уравнение

Е_п+1 _ ЕТпвп+1 (Уп+1) /вп (Уп),

позволяющее определить Ет(Уп+1)' Особо следует отметить, что эти уравнения справедливы только при изменении У вдоль изэнтропы Б=еопв1;. Такой расчет внутренней энергии и давления обеспечивает любую необходимую точность определения энтропии, а уравнение производства энтропии принимает вид

т ( * ). = о

Рассмотрим дивергентную разностную схему на волне сжатия. Вспомогательные величины вычисляются из уравнений на поверхности сильного разрыва (1)-(3) и уравнения состояния (4). Величины по обе стороны разрыва задаются следующим образом.

Если Ц+о,5 - ип_о,5 < о, то

1. Ц1 _ ип-0,5, (Р, У, Е, и)о _ (Р, У, Е, и)п+о,5, Ш > 0 при рп-0,5 > Р+о,5,

г+о,5 и г—о,!

-о,5, (РЛ,Е,и) о — \Р,У,Е,и) г+о,5, Ш >и при Р г—о,!

2. и = ип+о,5, (Р, У, Е, и)о = (Р, У, Е, и)п—о,5, Ш < 0 при Р-—о,5 < Р+о,5.

Остальные величины с индексом «1> находятся из (1)-(4). Если ограничиться рассмотрением только случая Ш > 0, то Р*, и* определятся уравнениями (9).

Для исследования диссипации энергии на сеточной ударной волне запишем в соответствии с [10-13] разностные законы сохранения (6)-(8) вместе со вспомогательными величинами (9) в дифференциальной форме

дУ ди ди дР де дРи

¿1, — + ——: = ¿2, ТТГ + = ¿3'

дг дм ’ дг дм ’ дг дм

Погрешности аппроксимации ¿1, ¿2, ¿3 имеют вид

тд2У Н д2и 2 тд2и д2и Нд2Р 2

Ш1 2 дг2 2 дш2 + , ¿2 2 дг2 + дш2 2 дш2 + ,

тд2е н д ( дР) н д ( ди) д ( ди) 2

2 дг2 2 дш \ дш) + 2 дш\ дш) + дш \ дш) + '

Согласно [10, 12], уравнение производства энтропии этой разностной схемы таково

дБ Т= ¿3 + Рш1 - и^2' дг

Подставив сюда выражения для ¿1, ¿2, ¿3, получим уравнение производства энтропии на сеточной ударной волне

тдБ = нш (ди)2 - та2 (дди)2 - т (дР)2 + о».

дг \ дш) 2 \ дш) 2 \ дш

Для дальнейшего упрощения этого уравнения воспользуемся уравнениями (9) для вспомогательных величин. Представим все входящие в них величины в виде рядов Тейлора в точке

п дР ди

гп, шг+о 5 . В результате получим связь производных —— и ——, которая на слабой ударной

ш ш

волне при Ш ~ а имеет вид

Р и 2

Я-- = а Я--+ 0 '

ш ш

Подставив эту связь в уравнение производства энтропии на сеточной ударной волне, получим уравнение

тж =,ш (1 - К >( £)2+02'

Рассмотрим монотонность этой разностной схемы на акустической (Ш = а) волне сжатия. Основные уравнения вместе со вспомогательными значениями (9) примут вид

та2

рп+1 _ рп ' ^ ¡ттп тти \

Рг+о,5 = Рг+о,5 1иг+о,5 - и г—о,5^ ,

^4+5 = игП+о,5 - Н (Р+1,5 - Р+о,5 - а (^¿+1,5 - 2иН-о,5 + и1—о,5)) '

Заменим Р и и их выражениями через инварианты а и в

аП+15 + вП+15 = аП+о,5 (1 - К) + вп+о,5 (1 + К) + К (аП—о,5 - вП—о,5) ,

аП+о15 - вП++о15 = аП+о,5 (1 - К) - вп+о,5 (1 - 3К) - 2Квп+1,5 + К (аП—о,5 - ^—о^) ' Сложив эти два уравнения, получим

аП+о,5 = аП+о,5 (1 - К) + КаП—о,5 - Квп+1,5 + 2Квп+о,5 - Квп—о,5' (17)

В случае вп=°опв1 уравнение (17) принимает вид

аПЙ5 = аП+о,5 (1 - К)+ аП—о,5К'

Оба коэффициента положительны при 0 < К < 1 и, таким образом, дивергентная разностная схема Куропатенко на волне сжатия монотонна.

Далее, аналогично [7, 9], исследуем дистракцию ударной волны. Перейдем к автомодельной переменной £ = ш - Шг. Тогда уравнения (6)-(8) вместе со вспомогательными значениями и*, Р*, определяемыми уравнениями (9) и погрешностями аппроксимации ¿1, ¿2 и ¿3, примут вид

ШУ' + и'------— У" - ^ и'' + О2 =0, (18)

тШ2 Н

Ши' - Р------— и'' - 2Р'' + НШи'' + О2 = 0, (19)

Ше' - (Ри)' - тЩ-Ри'' - Н (иР')' + Н (Ри')' - НШ (ии')' + О2 = 0' (20)

Проинтегрировав эти уравнения по £, получим

т^К2 Н

ШУ + и--------— У' - 2и' = ШУо + ио + О2, (21)

тШ2 Н

Ши - Р--------—и' - ^р' + НШи' = ШУо - Ро + о2, (22)

Ше - Ри - тЩ- (Ри)' - Нир' + Нри' - НШии' = Шео - Роио + О2' (23)

С помощью (18)-(20) заменим в (21)-(23) и' и Р' на У'. Затем с помощью (21)-(23) заменим и и Р на У' В результате для идеального газа получим уравнение, описывающее профиль

У (£)'

т+П . £ + (У - Уо)(У - У1) + О (т 2, Н2) = 0' (24)

(7 + 1) ¿£ У

Решение этого уравнения имеет вид

£= (7 +Н1()1(-¡, ’%<) (У‘1п (У - У1) -Уо 1п (У" - У»' (25)

Из (25) следует, что

£ = £о = +то при У = Уо, £ = £. = -то при У = У1'

Таким образом, дистракция разрыва в дивергентной разностной схеме метода Куропатенко при К < 1 бесконечна, а при К ^ 1 дистракция Ок1 ^ 0.

Для определения эффективной дистракции ОК 1 находятся точки пересечения прямой линии У(£), имеющей максимальный наклон УМ, со значениями Уо и У.. Эффективная ширина ударного слоя определяется уравнением

А£ = УУ-У' (26)

УМ

Для определения УМ продифференцируем (24) и в полученном уравнении положим нулю У''. В результате получим значение У = Ум, при котором У'' = 0, и выражение производной У' при У = УМ

Ум = у/УУъ ум = 221{1+-1К)(^- Ту)

Подставив УМ в (26), получим после деления на Н выражение для эффективной дис-тракции в дивергентной разностной схеме Куропатенко

чЭ 2 /л т^\ (+ Л/У1

Ок. = 7----тт (1 - К И ' (27)

К1 (7 + 1Г Лу/Уо-\/У1/ )

Рассмотренная дивергентная разностная схема устойчива при соотношении шагов

та2 2 (дР

ПТ < 1 где а = -(Ш, ,

3. Недивергентная разностная схема

Сетки для скорости и для термодинамических величин различаются. Значения Р, V, Е определяются в серединах сеточных интервалов по массе, значения скорости - в узлах сетки Ьп, Шг.

В случае волны разрежения при Ц+_1 — ЦП > 0 разностные уравнения имеют вид

хп+1 = хп + тЦ, (28)

Хп+1 — Хп+1

О! = Х+1 ь Х . (29)

Удельная внутренняя энергии Е^)^ и давление Рг”++15 находятся из уравнения (16) вместе с уравнением состояния (4), один из способов решения которого изложен в §2. После определения Рг”+15 и значения Рг”1+15 в соседнем интервале определяется скорость Цп+1

ггп+1 тт Рп+1 — Рп+1

и — Ц + Р‘+°’5 к Р-°’5 =0. (30)

Разностные уравнения (28)-(30) в дифференциальной форме имеют вид

ди дР дх ТТ дх тг

ТТГ + ^— = ^2, ^- — и = ^ —-----V = Ш5.

дЬ дш дЬ дш

Таким образом, независимыми погрешностями аппроксимации являются Ш2, ^4, ^5

т д2и т2 д3и Ь2 д3и т ди т2 д2и ^3 Ь2 д3х

^2 = — 2 — ¥ — 24 дШ3 + ° ,Ш4 =2 Ж — ¥ ° ,Ш5 = — 24 дШ3 + ° '

Для волн сжатия выполняется условие и?+1 — Цп < 0. В этом случае в соответствии с основной идеей метода считается, что в сеточном интервале находится сеточная ударная волна. Ее местоположение, знак скорости Ш и значения Ц и Цо определяются профилем Р(ш) в соответствии с условиями,

если Рп+15 >Рп-015, то Ш< 0, Шр = Шг+1, Ц = ип, Ц1 = ип+1,

если Р™++15 < Рп+15, то Ш > 0, Шр = Шг, Цо = Ц+1, Ц = Цп.

В качестве величин перед сильным разрывом берутся величины в сеточном интервале в момент Ьп

Ро = Рп+о,5, V) = У%,5, Ео = еп+0,5.

Поскольку скорости определены на концах интервала, для определения полной энергии +0,5 нужно доопределить скорость в середине интервала. Из предположения, что фон перед сеточной ударной волной постоянный, следует

ип [ иг+1, если Ш > 0,

иг+0>5 \ ип, если Ш< 0.

На следующем этапе определяются значения за сеточной ударной волной. Решая систему уравнений (1)-(4) с заданным значением Ц1, получим Р1, Vl, Е1, Ш.

Для определения У++5 применим мгновенный закон сохранения массы (10). Величина тШ это тот промежуток по ш, который сеточная ударная волна прошла со скоростью Ш за время т = Ьп+1 — Ьп. Сеточная ударная волна распространяется с постоянной скоростью Ш и постоянными величинами перед и за разрывом. Из этого следует, что Vl не изменяется со

временем и, таким образом, П]^1 = Пь Поскольку скорость перед разрывом постоянная, то и среднее на промежутке Н _ тШ значение У0га+1 должно быть равно У0 = У-+0 5. Из

уравнения (2) выразим У, и в полученном уравнении заменим У0, П и П значениями

У+0 5, Ц+і и П™. В результате получим

УЇЙ = У+о,5 + Н (П+ — ип). (31)

Это уравнение совпадает с уравнением (6), если в качестве вспомогательных значений П*+1, и* взять скорости концов интервала і + 0, 5

П * = ип и * = ип

иі+1 = иі+1, иі = Пі •

Такое определение вспомогательных значений не противоречит требованию инвариантности вспомогательных величин к изменению индекса і.

При переходе с одного временного шага ¿п на следующий шаг іп+1 величина Н = Ші+1 — Ші сохраняется

хп — хп хп+1 — хп+1 _ Хі+1 хі _ хі+1 хі /о9ч

Н = У п = Уп+1 • (32)

Уі+0,5 Уі+0,5

Подставив Уп+5 из (31) в уравнение (32), получим

Хп+1 _ Хп+1 = хп х і+1 хі х і+1

_ Хп + т (Пі+1 _ иі) •

Траектории частиц не зависят от индекса і. Следовательно

хп+1 = хп + тПп.

Рассмотрим теперь мгновенный закон сохранения количества движения (11). Значения ип+1 и и)п+1 это средние значения скорости в промежутках тШ и Ь — тШ. За ударной волной все величины постоянны и ип+1 = Ць Для среднего значения и?'+1 напишем разностное уравнение

и0‘+1 = Ц+015 — (Р+1 — Рж).^ • (33)

Будем считать, что в интервале г + 1,5 тоже распространяется сеточная ударная волна с Шг+1.5 > 0. Тогда

Р*+1 = (Р1)г+1.5 = Рп+1.5 — Шг+1.5 (Ц™+2 — Ц +О •

Очевидно, что

р* / рп

Р г+1 = Р г+0.5

и таким образом среднее значение Ц^+1 = Ц+о 5. Выразим Ц из (1) и подставим его вместе с ип+1 из (33) в (11). В результате получим разностное уравнение (7).

Мгновенный закон сохранения энергии имеет вид (13). Величину ^п+1 выразим из (3). С учетом введенных выше обозначений получим

^п+1 = ^п+0.5 + Ш (Р**и** — Р"+о.5ЦГ+1) • (34)

Величина еп+1 определяется из разностного уравнения

Рп+1 = Ґ-

^ л.

4+0,5 _ Н _ тШ (Р*+1 Пі*+1 _ Рі?+0,5иі+0 • (35)

Подставив (34), (35) в (13) и проведя простые преобразования, получим уравнение (8).

Переход от мгновенных законов сохранения (10), (11), (13) к разностным уравнениям (6)—(8) возможен при вполне конкретном выборе вспомогательных значений Р*, П*.

Из всех искомых величин остались неопределенными скорости ип+1 и Пі+1. В рассмотренном нами случае, когда Ш > 0, в момент Ьп скорость в интервале і + 0, 5 считалась постоянной и равной П?+1, т.е. Пі+0,5 = П+1. Завершить расчет величин в момент Ьп + 1 следует так, чтобы такая же ситуация была бы и в момент ¿п+1, т.е.

П++1 = идй,

Величина удельной внутренней энергии определяется уравнением

1 Л_ + ! Л2

En+1 _ Fn+i __ 1 (ттп+1 \

Ei+0,5 _ i+0,5 2 V i+0’5)

2 у г+°.5/

та

Недивергентная разностная схема устойчива при соотношении шагов —- < 1. Дистракция

Ь

разрывов определяется уравнением (27).

Разностные уравнения обладают двумя важными свойствами:

1. При т = 0 профили величин, заданных в момент Ьп, не изменяются.

тШ

2. При соотношении шагов —— = 1 ударная волна не размывается и при постоянном

Ь

фоне стационарная ударная волна «:прыгает> из точки в точку. Действительно, при тШ

—— = 1 из уравнений (10), (11), (13) следует, что

У+Й = V!, Ц++о15 = Ц1, 6"+о1б = 61, Ц+11 = Ц1.

Заключение

Метод определения вспомогательных величин Рг*, Ц*, основанный на использовании законов сохранения на сильном разрыве для расчета параметров сеточных ударных волн и на применении мгновенных законов сохранения, позволяет строить и дивергентные, и недивергентные разностные схемы. В [16, 17] описан богатый опыт их применения для моделирования ударных волн и волн разрежения.

Работа проводилась при финансовой поддержке РФФИ. Грант 13-01-00072.

Литература

1. Куропатенко, В.Ф. О разностных методах для уравнений гидродинамики / В.Ф. Куропатенко // Труды матем. инст. им. В.А.Стеклова. - 1966. - Т. 74, вып. 1. - С. 107-137.

2. Neumann, J. A Method for the Numerical Calculation of Hydrodynamical Shocks / J. Neumann, R. Richtmayer // J. Appl. Phys. - 1950. - V. 21, № 3. - P. 232-237.

3. Lax, P.D. Weak Solution of Nonlinear Hyperbolic Equations and Their Numerical Computations / P.D. Lax // Comn. Pure and Appl. Math. - 1954. - V. 7. - P. 159-193.

4. Годунов, С.К. Разностный метод расчета ударных волн / С.К. Годунов // Успехи математических наук. - 1957. - № 12, вып. 1. - С. 176-177.

5. Куропатенко, В.Ф. Метод расчета ударных волн /В.Ф. Куропатенко // ДАН СССР. -1960. - Т. 3, № 4. - С. 771-772.

6. Рождественский, Б.Л. Системы квазилинейных уравнений и их приложения к газовой динамике j Б.Л. Рождественский, Н.Н. Яненко. - М.: Наука, 1968.

7. Куропатенко, В.Ф. Исследование дистракции разрывов в методах расчета ударных волн j В.Ф. Куропатенко, И.Р. Макеева j j Математическое моделирование. - 2006. - Т. 18, № З. - С. 120-128.

8. Ступоченко, Е.В. Релаксационные процессы в ударных волнах j Е.В. Ступоченко.

С.А. Лосев, А.И. Осипов. - М.: Наука, 1965.

9. Куропатенко, В.Ф. Разностный метод расчета ударных волн с повышенными свойствами монотонности j В.Ф. Куропатенко, И.Р. Макеева j j Препринт ВНИИТФ. - 1997. - № 120.

10. Куропатенко, В.Ф. Локальная консервативность разностных схем для уравнений газовой динамики j В.Ф. Куропатенко jj Журнал выч. матем. и матем. физики. - 1985. -Т. 25, № 8. - С. 1176-1188.

11. Куропатенко, В.Ф. О полной консервативности разностных законов сохранения j В.Ф. Куропатенко jj Вопросы атомной науки и техники. Серия: Численные методы решения задач математической физики. - 1982. - Вып. З(11). - С. З-5.

12. Куропатенко, В.Ф. О точности вычисления энтропии в разностных схемах для уравнений газовой динамики j В.Ф. Куропатенко jj Численные методы механики сплошной среды: сб. - 1978. - Т. 9, № 7. - С. 49-59.

13. Куропатенко, В.Ф. Связь дивергентности с консервативностью разностных схем для уравнений газовой динамики / В.Ф. Куропатенко // Вопросы атомной науки и техники. Серия: Математическое моделирование физ. процессов. - 1990. - Вып. 2. - С. 6З-69.

14. Куропатенко, В.Ф. Методы расчета ударных волн j В.Ф. Куропатенко j j Энциклопедия низкотемпературной плазмы. Серия Б. Математическое моделирование в низкотемпературной плазме. Часть 2. - 2008. - Т. VII-I. - С. 496-506.

15. Куропатенко, В.Ф. Модели механики сплошной среды j В.Ф. Куропатенко. - Челябинск: ЧелГУ, 2007.

16. Куропатенко, В.Ф. О влиянии свойств разностных схем на математическое моделирование динамических процессов j В.Ф. Куропатенко, И.А. Доровских, И.Р. Макеева jj Вычислительные технологии. - 2006. - Т. 11, часть 2. - С. 9-11.

17. Комплекс программ ВОЛНА и неоднородный разностный метод для расчета неустано-вившихся движений сжимаемых сплошных сред j В.Ф. Куропатенко, Г.В. Коваленко. В.И. Кузнецова, Г.И. Михайлова jj Вопросы атомной науки и техники. Серия: Математическое моделирование физ. процессов. - 1989. - Вып. 2. - С. 9-25.

Валентин Федорович Куропатенко, доктор физико-математических наук, профессор.

главный научный сотрудник, Российский федеральный ядерный центр - Всероссийский

научно-исследовательский институт технической физики им. академика Е.И. Забабахина

(г. Снежинск, Челябинская обл., Российская Федерация), v.f.kuropatenko@rambler.ru.

Bulletin of the South Ural State University. Series "Mathematical Modelling, Programming & Computer Software",

2014, vol. 7, no. 1, pp. 62-75.

MSC 76.L, 74.S DOI: 10.14529/mmp140106

A Shock Capturing Method

V.F. Kuropatenko, Russian Federal Nuclear Center - Zababakhin Institute of Applied Physics, Snezhinsk, Russian Federation, v.f.kuropatenko@rambler.ru

Strong discontinuities, or shocks in continua are a result of external dynamic loads. On the shock surface the conservation laws take the form of nonlinear algebraic equations for jumps across the shock. Entropy jumps across a strong discontinuity, and just this jump differs shocks from waves where the quantities vary continuously. In the heterogeneous difference schemes, the shock is treated as a layer of a finite thickness comparable with the cell size. This property of finite-difference schemes was called distraction. Since the state behind a shock is related to the state before it by the Hugoniot, in the distraction region there must act a mechanism that increases entropy. The physical viscosity and heat conductivity in continuum mechanics equations do not make it unnecessary to introduce a shock surface and hence cannot make the distraction length comparable with a few cells of the difference mesh. The paper considers a number of finite difference schemes where energy dissipation in the distraction region is defined by equations which are valid on the shock surface.

Keywords: shock wave; differential method; distraction; energy dissipation;

conservation laws.

References

1. Kuropatenko V.F. Finite Difference Methods for Hydrodynamics Equations [O raznostnykh metodakh dlya uravneniy gidrodinamiki]. Trudy matematicheskogo instituía im. V.A. Steklova [Proceedings of the Steklov Institute of Mathematics], 1966, vol. 74, part 1, pp. 107-137.

2. Neumann J., Richtmayer R. A Method for the Numerical Calculation of Hydrodynamic Shocks. J. Appl. Phys., 1950, vol. 21, no. 3, pp. 232-237. DOI: 10.1063/1.1699639

3. Lax P.D. Weak Solution of Nonlinear Hyperbolic Equations and Their Numerical Computations. Comn. Pure and Appl. Math., 1954, vol. 7, pp. 159-193. DOI: 10.1002/cpa.3160070112

4. Godunov S.K. A Finite-Difference Method for Shock Calculation [Raznostnyy metod rascheta udarnykh voln]. Uspekhi Matematicheskikh Nauk [Russian Mathematical Surveys], 1957, no. 12, issue 1, pp. 176-177.

5. Kuropatenko V.F. A Shock Calculation Method. DAN SSSR, 1960, vol. 3, no. 4, pp. 771-772.

6. Rohzdestvensky B.l, Yanenko N.N. Sistemy kvazilineynykh uravneniy i ikh prilozheniya k gazovoy dinamike [Systems of Quasi-Linear Equations and Their Applications to Hydrodynamics]. Moscow, Nauka, 1968. 592 p.

7. Kuropatenko V.F., Makeyeva I.R. Discontinuity Distraction in Shock Calculation Methods [Issledovanie distraktsii razryvov v metodakh rascheta udarnykh voln]. Matematicheskoe modelirovanie [Mathematical Models and Computer Simulations], 2006, vol. 18, no. 3, pp. 120-128.

8. Stupochenko E.V., Losev S.A., Osipov A.I. Relaksatsionnye protsessy v udarnykh volnakh [Relaxation Processes in Shock Waves]. Moscow, Nauka, 1965. 484 p.

9. Kuropatenko V.F, Makeyeva I.R. Raznostnyy metod rascheta udarnykh voln s povyshennymi svoystvami monotonnosti [A Higher-Monotonicity Finite-Difference Shock Capture Method]. VNIITF Preprint, 1997, no. 120.

10. Kuropatenko V.F. Local Conservatism of Difference Schemes for Hydrodynamics Equations [Lokal’naya konservativnost’ raznostnykh skhem dlya uravneniy gazovoy dinamiki]. Zhurnal vychislitel’noy matematiki i matematicheskoy fiziki [Computational Mathematics and Mathematical Physics], 1985, vol. 25, no. 8, pp. 1176-1188.

11. Kuropatenko V.F. Ultimate Conservatism of Finite-Difference Conservation Laws [O

polnoy konservativnosti raznostnykh zakonov sokhraneniya]. Voprosy atomnoy nauki i tekhniki. Seriya: Chislennye metody resheniya zadach matematicheskoy fiziki [Atomic Science and Engineering. Series: Numerical Methods of Mathematical Physics], Moscow, 1982.

issue 3 (11), pp. 3-5.

12. Kuropatenko V.F. Entropy Accuracy in Finite Difference Schemes for Hydrodynamics Equations [O tochnosti vychisleniya entropii v raznostnykh skhemakh dlya uravneniy gazovoy dinamiki]. Chislennye metody mekhaniki sploshnoy sredy [Numerical Methods for Continuum Mechanics], Novosibirsk, 1978, vol. 9, no. 7, pp. 49-59.

13. Kuropatenko V.F. Divergence and Conservatism of Finite-Difference Schemes for Hydrodynamics Equations [Svyaz’ divergentnosti s konservativnost’yu raznostnykh skhem dlya uravneniy gazovoy dinamiki]. Voprosy atomnoy nauki i tekhniki. Seriya: Matematicheskoe modelirovanie fizicheskikh protsessov [Atomic Science and Engineering. Series: Mathematical Modeling of Physical Processes], 1990, issue 2, pp. 63-69.

14. Kuropatenko V.F. Shock Calculation Methods [Metody rascheta udarnykh voln]. Entsiklopediya nizkotemperaturnoy plazmy. Seriya B. Matematicheskoe modelirovanie v nizkotemperaturnoy plazme [Encyclopedia of Low-Temperature Plasma. Series B. Mathematical Modelling of Low-Temperature Plasma], part 2, vol. VII-I, 2008, pp. 496-506.

15. Kuropatenko V.F. Modeli mekhaniki sploshnoy sredy [Continuum mechanics models]. Chelyabinsk, Chelyabinsk State University, 2007, 302 p.

16. Kuropatenko V.F., Dorovskikh I.A., Makeyeva I.R. The Properties of Finite Difference Schemes and Simulation of Dynamic Processes [O vliyanii svoystv raznostnykh skhem na matematicheskoe modelirovanie dinamicheskikh protsessov]. Vychislitel’nye tekhnologii [Computating Technologies], 2006, vol. 11, part 2, pp. 9-11.

17. Kuropatenko V.F., Kovalenko G.V., Kuznetsova V.I., Mikhaylova G.I. Complex Programs VOLNA and Method for Transient Flows of Continua [Kompleks programm VOLNA i neodnorodnyy raznostnyy metod dlya rascheta neustanovivshikhsya dvizheniy szhimaemykh sploshnykh sred]. Voprosy atomnoy nauki i tekhniki. Seriya: Matematicheskoe modelirovanie fizicheskikh protsessov [Atomic Science and Engineering, Series: Mathematical Modeling of Physical Processes], 1989, issue 2, pp. 9-25.

Поступила в редакцию 13 декабря 2Ü13 г.