Диссертация (Консервативные дискретизации одномерных квазигазо- и квазигидродинамических систем уравнений), страница 7

Обратим такжевнимание на то, что множители Dhρ и Dhθ обращаются в 0 приp ε1s−+adθρ11ρ ρθθ:=− −+,(2.55)1[ρ]0ρaθ θ dρ p pε11+−aρ dθ(2.56)[ε]3 := θ− θ+ aρ dθ s −θ[ρ]0ρθсоответственно. Поскольку по формулам (2.46) и равенству Максвелла (2.10)имеемεpεε pρ2 ε ρs−=− 2 , s−= ,=− 2 ,(2.57)θ ρρθθ θ θ2θ θθто с помощью формул (2.40) для построенных дискретизаций можно записатьинтегральные представленияD E pp11 aθ−2ρ ρ θρ− ρ+ θ11=+,(2.58)1[ρ]0ρa hp iθ[ε]3 = θ− θ+ aρDεEθ2−θρ ρ θ D E 11ερ−aρ ρ 2 2.[ρ]0ρθ θ43(2.59)Отметим, что первые слагаемые в формулах (2.58) и (2.59) являются главными, а остальные слагаемые имеют порядок аппроксимации O(h2max ) на гладких решениях.

Из последних формул следует, что выражения (2.55) и (2.56)действительно дают аппроксимации плотности и внутренней энергии, хотя ивесьма нетривиального вида.Лемма 2.1. Если pρ > 0 и(3.19) положительна.136 ρ+ /ρ− 6 3, то правая часть формулыДоказательство. Приведем правую часть формулы (2.58) к общему знаменателю (который положителен) и обозначим числитель через через ∆.

Достаточно считать, что ρ− 6= ρ+ (иначе результат очевиден). Сделаем заменупеременной r = [ρ](α) . С помощью интегрирования по частям получим p 1 1 hpiaθ−=ρ2 ρ θρ− ρ+ θZ ρ+ 11p(ρ− , θ− ) p(ρ+ , θ+ )p(r, θ− ) p(r, θ+ ) dr=+−+=2(ρ+ − ρ− ) ρ−θ−θ+r2 2ρ− ρ+θ−θ+Z ρ+ 11 pr (r, θ− )11 pr (r, θ+ )1=−+−dr.2(ρ+ − ρ− ) ρ−r ρ+θ−r ρ−θ+Поэтому верно интегральное представлениеZ ρ+ 111 pr (r, θ− )1∆=+ −+2(ρ+ − ρ− ) ρ−ρr ρ+θ− 11 pr (r, θ+ )1+ −dr.+ρr ρ−θ+В нем при 13 6 ρ+ /ρ− 6 3 оба множителя в круглых скобках неотрицательны.

Как следствие ∆ > 0 и правая часть формулы (2.58) (значит, и (2.55))положительна.Для уравнений состояния совершенного политропного газа (2.14) построенные аппроксимации резко упрощаются и принимают вид, найденный в [31][ρ]0 =1, [ε]3 = ε− ε+ ln(ε− ; ε+ ),ln(ρ− ; ρ+ )где ln(α; β) — разделенная разность для логарифмической функцииln(α; β) =ln β − ln α1при α 6= β, ln(α; α) = , α > 0, β > 0.β−αα44(2.60)Действительно, в этом случае 1 hpi[ρ]ε1p=R= R ln(ρ− ; ρ+ ),aθ dρ = R,=R, −aθ dρ s −θρ− ρ+ θρ− ρ+ρθεpaρ dθ s −= θ− θ+ c2V ln(ε− ; ε+ ), dθ = 0.θθ2.4. Преобразуем теперь сумму основных первых трех слагаемых в правойчасти равенства (2.48).

В силу уравнения баланса внутренней энергии (2.44)и формул (2.28) и (2.30) получим для нее промежуточное представление11{∂t (ρε) + δ ∗ (j[ε]3 ) + pδ ∗ ([u] − w)} = (−δ ∗ q + [Πδu + wδp + Q]∗ ) =θθ ∗ 111(2.61)= qδ + (Πδu + wδp + Q)+ B2h− δ∗ qθθθс сеточным дисбалансным слагаемым1B2h := 0.25h2+ (Πδu + wδp + Q)δ .θ(2.62)Из определений −q и Π, см. (2.38) и (2.37), с учетом формул (2.49) выводим 1κ(δθ)2 4A12 1=+ µ(δu)+,(2.63)qδ + (Πδu + wδp + Q)θθθ− θ+3θθ− θ+где введено обозначение[p]12A := τ [ρ][u] δε − 2 δρ − [u]Q δθ+[ρ]+ [ρ][u]ŵ + τ ([u]δp + [ρ][Cs2 ]4 δu − ([γQ ]5 − 1)Q) δu + wδp + Q [θ].Разобьем A на слагаемое без множителя Q и слагаемые с ним:A = Ã − τ {[u]δθ + ([γQ ]5 − 1)[θ]δu} Q + [θ]Q.(2.64)Предварительно заметим, вспомнив определения w и ŵ, см.

(2.36), чтоτ[ρ] 2[u][ρ][u]ŵδu + wδp = ŵ[ρ][u]δu + ŵ +[u]δ(ρu) δp =ŵ + τδ(ρu) · δp.[ρ]τ[ρ]Отсюдаà =[ρ][θ] 2ŵ + τ Â,τ45(2.65)где введено обозначение[u][p]1222 := [θ]δ(ρu) · δp + [u]δu · δp + [ρ][Cs ]4 (δu) + [ρ][u] δε − 2 δρ δθ.[ρ][ρ]Подставим в последнее выражение формулыδp = (aθ dρ p)δρ + (aρ dθ p)δθ, δε = (aθ dρ ε)δρ + (aρ dθ ε)δθ,см.

(2.39). Существенно, что в силу условий термодинамической устойчивости(2.11) и интегральных представлений (2.40) имеемaθ dρ p > 0, aρ dθ ε > 0.Перегруппировав слагаемые и учтя формулу δ(ρu) = δρ · [u] + [ρ]δu, получимaθ d ρ p [u]δ(ρu) · δρ + [ρ][u]δu · δρ + ([ρ]δu)2 +[ρ]2[θ](adp)ρθ(δu)2 ++[ρ][θ] [Cs2 ]4 − aθ dρ p −2[ρ] aρ dθ ε[θ]2 (aρ dθ p)2[u]+(δu)2 + [θ] (aρ dθ p) [ρ]δu · δθ + (aρ dθ p)[u]δu · δθ +[ρ]aρ dθ ε[ρ][u]2[p]1222δρ · δθ.+[ρ][u] (aρ dθ ε)(δθ) + [θ]aρ dθ p + [ρ][u] aθ dρ ε − 2[ρ][ρ]Â = [θ]В этом выражении выделим квадраты сумм и дисбалансные слагаемые2aθ dρ p[θ]aρ dθ p2Â = [θ]{δ(ρu)} + [ρ]aρ dθ εδu + [u]δθ +[ρ][ρ]aρ dθ ε+[θ]Dhu (δu)2 + [θ]Dhρ θ δρ · δθ(2.66)с множителями2[θ](adp)ρθ,Dhu := [ρ] [Cs2 ]4 − aθ dρ p −[ρ]2 aρ dθ ε[u]2[ρ]2[p]1Dhρ θ :=.aρ dθ p +aθ dρ ε − 2[ρ][θ][ρ](2.67)(2.68)Обратим внимание на то, что множители Dhu и Dhρ θ равны 0 при[Cs2 ]4 := aθ dρ p +[θ](aρ dθ p)2> 0, [p]1 := [θ]aρ dθ p + [ρ]2 aθ dρ ε2[ρ] aρ dθ ε46(2.69)соответственно.

Первое из этих выражений специальным образом аппроксимирует квадрат скорости звука именно в виде (2.12), а второе — равенствоМаксвелла (2.10) (а отнюдь не само уравнение состояния p = p(ρ, θ)). Подобно [26], отметим, что дисбалансное слагаемое с множителем Dhu будет такженеотрицательным, если величина [Cs2 ]4 будет больше или равна указанной в(2.69) величины (а не только равна ей).С помощью формул (2.40) и равенства Максвелла можно записать интегральное представление для второй формулы (2.69) !pθpθ[θ]2[p]1 = [θ]aρ hpθ iθ + [ρ]2 aθ hερ iρ = [ρ]2 aθahpi−a+[ρ].ρθθθρ2 ρ[ρ]2ρ2 ρДля уравнений состояния совершенного политропного газа (2.14) построенные аппроксимации (2.69) резко упрощаются и принимают вид из [31][ρ][Cs2 ]4 = γ[p]1 , [p]1 = R[ρ][θ].(2.70)Напомним, что вместо дискретизаций [ρ]0 и [ε]3 вида (2.60) и указанной [p]1в стандартном (описанном в предыдущем разделе) варианте используютсяпросто [ρ], [ε], [p].Кроме того, вспомнив слагаемые с множителем Q, см.

(2.64), запишем2[θ]aρ dθ pδu + [u]δθ − τ {[u]δθ + ([γQ ]5 − 1)[θ]δu} Q + [θ]Q =τ [ρ](aρ dθ ε)[ρ]aρ dθ ε2[θ]aρ dθ pQτQ= τ [ρ](aρ dθ ε)δu + [u]δθ −+ [θ]Q 1 −−[ρ]aρ dθ ε2[ρ]aρ dθ ε4[ρ][θ]aρ dθ ε−τ [θ]DhQ (δu)Q,(2.71)где последнее слагаемое является дисбалансным с множителемDhQ := [γQ ]5 − 1 −aρ dθ p.[ρ]aρ dθ ε(2.72)Этот множитель равен 0 при[γQ ]5 := 1 +aρ dθ p,[ρ]aρ dθ ε(2.73)что специальным образом аппроксимирует вторую из формул (2.12) (а дляуравнений состояния совершенного политропного газа дает просто [γQ ]5 = γ).В итоге, сведя воедино формулы (2.48), (2.51), (2.61), (2.63)–(2.66), (2.71),придем к следующему результату.47Теорема 2.1. Для пространственной дискретизации (2.31)–(2.38) КГДсистемы уравнений справедливо следующее уравнение баланса энтропии 1∗∗∂t (ρs) + δ (j[s]) = δ −q+ B1h − B2h +θκ(δθ)2[ρ]4+ ([u] − w)(Dhρ δρ + Dhθ δθ) ++µ(δu)2 + ŵ2 +θ− θ+3τ2τ aθ dρ pτ [ρ]aρ dθ ε [θ]aρ dθ pQ2+{δ(ρu)} +δu + [u]δθ −+[ρ][θ][ρ]aρ dθ ε2[ρ]aρ dθ ε ∗τQ1+τ Dhu (δu)2 + τ Dhρ θ δρ · δθ + Q 1 −− τ DhQ (δu)Q,4[ρ][θ]aρ dθ εθ(2.74)где дивергентные дисбалансные слагаемые B1h и B2h задаются формулами(2.52), (2.62), а множители пяти недивергентных дисбалансных слагаемыхDhρ , Dhθ , Dhu , Dhρ θ , DhQ — формулами (2.53), (2.54), (2.67), (2.68), (2.72).В этом уравнении пять стоящих подряд квадратичных слагаемых производства энтропии всегда неотрицательны, а слагаемое, содержащее величину Q дважды, неотрицательно при выполнении условия τ Q 6 4[ρ][θ]aρ dθ ε.Недивергентные дисбалансные слагаемые с множителями Dhρ , Dhθ , Dhu ,Dhρ θ , DhQ обращаются в 0 при выборе дискретизаций [ρ]0 , [ε]3 , [Cs2 ]4 , [p]1 ,[γQ ]5 согласно формулам (2.55), (2.56), (2.69), (2.73) соответственно.Выведенное дискретное по пространству уравнение баланса энтропии (2.74)является аналогом дифференциального уравнения (2.15).Как и в дифференцальном случае, свойство неотрицательности слагаемых производства энтропии сохраняет силу и при µ > 0, κ > 0, τ > 0;для выполнения последнего слагаемое ([ρ]/τ )ŵ2 следует переписать в виде(τ /[ρ])([ρ][u]δu + δp)2 .Отметим, что дисбалансные слагаемые B1h , B2h и с множителями Dhρ ,Dhθ являются величинами порядка аппроксимации O(h2max ), а с множителямиDhu , Dhρ θ , DhQ — даже O(τ h2max ), но это на гладких решениях; для случаяже негладких решений, который представляет основной интерес, их вкладможет быть существенным, и это проявляется в расчетах.При выводе уравнения баланса энтропии (2.74) вместо формулы (2.39)можно было опираться на “более аналитичную” формулуZ 1δψ(ρ, θ) = hψρ iδρ + hψθ iδθ, где hϕ(ρ, θ)i :=ϕ([ρ](α) , [θ](α) ) dα.

(2.75)048Тогда получилось бы уравнение баланса энтропии типа (2.74), но с заменойвсех величин вида aθ dρ ψ, aρ dθ ψ на hψρ i, hψθ i соответственно, в том числев формулах (2.53), (2.54), (2.67), (2.68), (2.72) для величин Dhρ , Dhθ , Dhu ,Dhρ θ , DhQ (в первых двух из которых надо учесть формулы (2.57)). Однако соответствующие выражения для [ρ]0 , [ε]3 , [Cs2 ]4 , [p]1 , [γQ ]5 , обнуляющиенедивергентные дисбалансные слагаемые, редко вычислимы в явном виде.2.4Результаты численных экспериментовРяд приводимых ниже понятий ранее уже были введены в предыдущемразделе; повторим их для полноты изложения.Аналогично предыдущему разделу, выполненные численные эксперименты относятся к решению системы уравнений Эйлера невязкого нетеплопроводного газа.

Как хорошо известно, в отличие от КГД системы (2.1)–(2.7) онане содержит диссипативных слагаемых с коэффициентами τ , µ, κ. При численном решении эти слагаемые вводятся как искусственные регуляризаторыи для последних используются, как и в [59], формулы указанного в [15, 43]типаµcpαSµ = αS τ ρpρ , κ ==τ ρcV Cs2(2.76)αPαPс учетом равенства Cs2 = (cp /cV )pρ [37], где αS и αP — числа Шмидта иПрандтля (положительные постоянные).Далее, дискретизации τ , µ, κ по пространству берутся в видеhαSτ = αp, µ = αS τ [ρ][pρ ], κ =τ [ρ][cV ][Cs2 ]4 .2αP[Cs ]4Для дискретизации по времени используется простейший явный метод Эйлера (за одним указанным ниже исключением) с переменным шагом, выбираемым в соответствии с условием устойчивости Куранта по формуле∆t = β min∗ωh|[u]| +hp[Cs2 ]4.Здесь 0 < α < 1 и 0 < β 6 1.1 — числовые параметры, наилучшие значениякоторых подбираются отдельно в каждом расчете.Во всех тестах рассматриваются задачи Римана о распаде разрывов и берутся кусочно-постоянные начальные данные((ρL , uL , pR ), x 6 X2(ρ0 , u0 , p0 )(x) =.(ρR , uR , pR ), x > X249При этом X = 1 в тестах 1, 2 либо X = 10 в тестах 3.