1612725465-eb5831e2a489b993bf8c6c33f84310db (Хакимзянов, Черный - Методы вычислений(часть 3)), страница 9

Консервативные схемы3.1. Ранее при построении разностных схем мы аппроксимировалидифференциальные операторы разностными и доказывали сходимостьразностных схем, предполагая, что дифференциальная задача имеет достаточно гладкое решение. Однако не все физические процессы описываются с помощью дифференцируемых функций. Например, на ударнойволне скорость газа, его плотность, давление и температура являютсяразрывными функциями. Соответствующие дифференциальные задачине имеют гладких решений. Для уравнения теплопроводности с переменным коэффициентом теплопроводностиut = [ν(x, t)ux ]x + f (x, t)(3.1)не существует гладкого решения, если коэффициент теплопроводностивсюду является гладкой функцией, но в одной точке x0 имеет разрывпервого рода, т.

е. νлев (x0 , t) 6= νпр (x0 , t). В этом случае классическаяпостановка задачи для уравнения (3.1) может быть заменена, например,на следующую неклассическую постановку: найти функцию u(x, t), которая во всех точках (x, t), за исключением точек (x0 , t), удовлетворяетуравнению (3.1), а в точках (x0 , t) — условиям сопряженияuлев (x0 , t) = uпр (x0 , t);(3.2)(νux )лев (x0 , t) = (νux )пр (x0 , t).(3.3)Первое из условий сопряжения следует из требования непрерывностирешения, второе — непрерывности потока тепла.В общем случае, когда дифференциальное уравнение не имеет смысла, для получения разностных уравнений будем пользоваться записьюфизических законов сохранения не в дифференциальной, а в интегральной форме.

Законы сохранения в интегральной форме имеют смысли для негладких функций, которые нельзя дифференцировать, но интегрировать можно.Вспомним, что в основе дифференциального уравнения теплопроводности лежит интегральный закон сохранения тепла:x+∆xZx+∆xZu(ξ, t + ∆t)dξ −xu(ξ, t)dξ =x55t+∆tZ=t+∆tZW (x, τ )dτ −tt+∆tx+∆xZZW (x + ∆x, τ )dτ +tf (ξ, τ )dξdτ, (3.4)txгде ∆t и ∆x произвольные числа, u(x, t) — температура, W (x, t) — потоктепла (количество тепла, протекшего в единицу времени через единичную площадку)∂uW (x, t) = −ν(x, t) (x, t).∂xИнтегральный закон (3.4) получается при подсчете баланса тепла напроизвольном участке стержня за произвольный промежуток времени и потому называется еще уравнением баланса тепла. Приведенныйзакон устанавливает, что изменение количества тепла в стержне на отрезке [x, x + ∆x] за время ∆t определяется:— разностью количества тепла, которое втекло и вытекло через сечения стержня x и x + ∆x за время ∆t;— количеством тепла, выделившимся на отрезке [x, x + ∆x] за время ∆t за счет распределенных на нем источников тепла с плотностьюf (x, t).Для получения разностного уравнения рассмотрим интегральное уравнениебалансатепланаэлементарнойячейкесеткиωjn = [xj−1/2 , xj+1/2 ] × [tn , tn+1 ], т.

е. на отрезке xj−1/2 ≤ x ≤ xj+1/2 запромежуток времени tn ≤ t ≤ tn+1xj+1/2xj+1/2ZZu(x, tn+1xj−1/2xj−1/2n+1tZ=u(x, tn )dx =)dx −n+1 xj+1/2tZZn+1tZW (xj−1/2 , t)dt −tnW (xj+1/2 , t)dt +tnf (x, t)dxdt.tn xj−1/2Аппроксимируем входящие в уравнение баланса интегралы приближенными формуламиxj+1/2Zu(x, tn )dx ∼ unj h,xj−1/2n+1tZn+1nW (xj±1/2 , t)dt ∼ (σWj±1/2+ (1 − σ)Wj±1/2)τ,tn56n+1 xj+1/2tZZf (x, t)dxdt ∼ τ hϕnj ,tn xj−1/2nnWj+1/2= −νj+1/2unj+1 − unj.h(3.5)nПри этом νj+1/2определяется формулойnνj+1/2=nνjn + νj+1,2νjn = ν(xj , tn ),(3.6)а для вычисления ϕnj можно использовать одну из формул, приведенных в подп.

1.3 (в зависимости от выбранного значения параметра σ).После такой замены интегральное уравнение баланса тепла превращается в дискретное уравнение баланса тепла для элементарной ячейки³´n+1nun+1h − unj h = σWj−1/2+ (1 − σ)Wj−1/2τ−j³´n+1n− σWj+1/2+ (1 − σ)Wj+1/2τ + ϕnj τ h,(3.7)которое после деления на τ h дает разностное уравнение"#n+1n+1un+1un+1− un+1− unj1 n+1 uj+1 − ujjjj−1n+1ν+=σ− νj−1/2τh j+1/2hh+(1 − σ)·¸unj+1 − unjunj − unj−11 nnνj+1/2− νj−1/2+ ϕnj ,hhhт. е.

схему с весамиun+1− unjj= σΛn+1 un+1+ (1 − σ)Λn unj + ϕnj ,jτгдеΛn unj·¸unj+1 − unjunj − unj−11 nn=ν− νj−1/2.h j+1/2hhЕсли ввести обозначенияnanj = νj−1/2,nanj+1 = νj+1/2,57(3.8)(3.9)то оператор Λn можно записать в виде второй разностной производнойс переменным коэффициентом:(an unx̄ )j+1 − (an unx̄ )j= (an unx̄ )x,j .hМетод получения разностных уравнений, основанный на аппроксимации интегральных законов сохранения для элементарных ячеек, называется интегро-интерполяционным методом. С помощью этого метода мы получили дискретное уравнение баланса тепла (3.7) толькодля элементарной ячейки.

Но исходный интегральный закон выполняется для произвольного отрезка по x и для произвольного промежуткавремени. Оказывается, что это свойство справедливо и на дискретномуровне, а именно из того, что дискретное уравнение баланса тепла (3.7)справедливо для каждой элементарной ячейки, будет следовать, чтои для области, составленной из некоторого количества элементарныхячеек, дискретное уравнение баланса тепла также будет выполняться.Покажем это в самом простом случае, в случае явной (σ = 0) схемы(3.8) и прямоугольной области Ω = [a, b] × [tp , tp+m ], где a = xk−1/2 ,b = xl+1/2 , tp+m = tp +mτ , k ≤ l, m ≥ 1, составленной из элементарныхячеек ωjn :l, p+m−1[Ω=ωjn .Λn unj =j=k, n=pДля такой области интегральный закон сохранения тепла (3.4) записывается так:ZbZbp+mu(x, t)dx − u(x, tp )dx =ap+mtZ=W (a, t)dt −tpap+mtZp+mtZZbW (b, t)dt +tpf (x, t)dxdt.tp(3.10)aВыпишем дискретные уравнения баланса тепла (3.7) для элементарных ячеек при всех значениях j = k, .

. . , l:³´nnun+1h − unk h = Wk−1/2− Wk+1/2τ + ϕnk τ h,k³´nnnnun+1h−uh=W−Wk+1k+1/2k+3/2 τ + ϕk+1 τ h,k+1. . . . . . . . . . . . . . .³. . . . . . . . . . . . . . . . . .´. . . . . . . . .un+1hlnn− unl h = Wl−1/2− Wl+1/2τ + ϕnl τ h58и просуммируем их почленно сначала по индексу j:lXun+1h−jj=klXnnunj h = Wk−1/2τ − Wl+1/2τ +τj=klXϕnj h,j=kа затем по n = p, . . . , p + m − 1. В результате получим дискретныйаналогllXXup+mh−upj h =jj=k=p+m−1Xn=pnWk−1/2τ−j=kp+m−1Xnτ+Wl+1/2n=p(3.11)p+m−1lX Xn=pϕnj hτj=kинтегрального закона сохранения тепла (3.10) во всей области Ω.Важно отметить, что сокращение слагаемых при суммировании произошло благодаря тому, что для соседних элементарных ячеек поток наразделяющей их общей границе аппроксимировался по одной и той жеформуле (3.5).Определение. Разностные схемы, выражающие дискретные законы сохранения в элементарных ячейках и в составных областях,полученных объединением элементарных ячеек, называются консервативными.Отметим еще раз, что для консервативных схем дискретный законсохранения в составной области является алгебраическим следствиемразностных уравнений.3.2.

Введем теперь понятие дивергентных схем. Так называютсясхемы, полученные в результате конечно-разностной аппроксимациидифференциальных уравнений, записанных в дивергентной форме. Так,уравнение (3.1) имеет дивергентную форму, поскольку дифференциальный оператор этого уравнения имеет дивергентную формуdiv a ≡ ut + wx = f (x, t),где w = −νux — тепловой поток, a(x, t) = (u(x, t), w(x, t)). Сравниваяявную дивергентную схему·¸un+1− unjunj+1 − unjunj − unj−11 njn=νj+1/2− νj−1/2+ ϕnj ,τhhh59(3.12)где ϕnj = f (xj , tn ), с явной (σ = 0) консервативной схемой (3.8), видим,что они совпадают. Этот факт имеет место и в общем случае в томсмысле, что применяя интегро-интерполяционный метод к интегральным законам сохранения, можно построить разностные схемы, совпадающие с дивергентными схемами, аппроксимирующими дифференциальные уравнения, записанные в дивергентной форме.

Поэтому можносказать, что консервативные схемы являются дивергентными.Если уравнение (3.1) записать в недивергентной формеut = νuxx + νx ux + f,(3.13)то аппроксимирующая его разностная схемаnnun+1− unjνj+1− νj−1unj+1 − unj−1j= νjn unx̄x,j +·+ ϕnjτ2h2h(3.14)имеет название недивергентной схемы. Недивергентная схема, как правило, не является консервативной. Покажем это на примере недивергентной схемы (3.14), положив для простоты, что функция ν не зависитот времени ( ν = ν(x)) и f ≡ 0.Лемма 3.1. Для недивергентной схемы (3.14) дискретный аналогзакона сохранения (3.11) не выполняется.Д о к а з а т е л ь с т в о. Перепишем недивергентную схему (3.14)в виде¡¢1ut,j = νj unx̄x,j + (νx,j + νx̄,j ) unx,j + unx̄,j(3.15)4и выясним, чем она отличается от явной консервативной схемы (3.8)ut,j = (aunx̄ )x,j ,где aj = νj−1/2 , aj+1 = νj+1/2 .Для этого преобразуем правую часть консервативной схемы:µ¶1 νj + νj+1νj + νj−1(aux̄ )x,j =ux̄,j+1 −ux̄,j =h22µ¶1 νj + νj + hνx,jνj + νj − hνx̄,j=ux̄,j+1 −ux̄,j =h2211= νj ux̄x,j + νx,j ux,j + νx̄,j ux̄,j .2260(3.16)Следовательно,11νj ux̄x,j = (aux̄ )x,j − νx,j ux,j − νx̄,j ux̄,j22и недивергентная схема (3.15) может быть записана в следующем видеut,j = (aunx̄ )x,j −¢11¡1νx,j unx,j + νx,j unx̄,j + νx̄,j unx,j + νx̄,j unx̄,j =− νx,j unx,j − νx̄,j unx̄,j +224¢1¡n= (aux̄ )x,j +νx̄,j unx,j − νx,j unx,j + νx,j unx̄,j − νx̄,j unx̄,j =4hh2h= (aunx̄ )x,j − unx,j νx̄x,j + unx̄,j νx̄x,j = (aunx̄ )x,j − unx̄x,j νx̄x,j .444Таким образом, недивергентная схема преобразована к консервативной схеме с дополнительным членом в правой части:ut,j = (aunx̄ )x,j −h2 nuνx̄x,j .4 x̄x,j(3.17)Умножим полученные уравнения на hτ :³´h2nnun+1h − unj h = Wj−1/2− Wj+1/2τ − unx̄x,j νx̄x,j τ hj4и просуммируем почленно по всем элементарным ячейкам составнойобласти Ω.

В результате вместо дискретного закона сохранения (3.11)lXup+mhj−j=klXupj h=p+m−1XnWk−1/2τ−n=pj=kp+m−1XnWl+1/2τ(3.18)n=pполучим следующее равенство:lXj=kup+mhj−lXupj hj=k=p+m−1XnWk−1/2τn=p−p+m−1XnWl+1/2τ − R,n=pгдеR=p+m−1 lh2 X X nux̄x,j νx̄x,j hτ .4 n=pj=kТаким образом, действительно, недивергентная схема не является консервативной схемой.61Теперь понятно, в каких случаях важно, чтобы схема была консервативной, а в каких можно пренебречь этим свойством. Консервативнаясхема после умножения на hτ и суммирования по элементарным ячейкам области Ω = [a, b] × [tp , tp+m ], где a = xk−1/2 , b = xl+1/2 , привелак дискретному аналогу интегрального закона сохранения теплаZbp+mtZZbu(x, tp+m )dx −au(x, tp )dx =p+mtZW (a, t)dt −tpaW (b, t)dt,tpа для недивергентного уравнения мы получили дискретный аналог следующего интегрального уравнения:ZbZbu(x, tp+map+mtZ=ap+mtZW (a, t)dt −tpu(x, tp )dx =)dx −tph2W (b, t)dt −4p+mtZZbνxx uxx dxdt,tpaкоторое отличается от уравнения баланса тепла последним членом: дисбалансом.Если все функции гладкие, то отличие будет небольшим — в пределах погрешности аппроксимации.