1612725465-542b3179b36a4849700e0b2ecf7da111 (828844), страница 15
Текст из файла (страница 15)
Поэтому необходимое условиеустойчивости записывается в виде неравенства d ≤ 0, или 0 ≤ q ≤ 2,или2r sin2 φ20≤≤ 2,1 + 4σr sin2 φ298которое должно выполняться для любого φ ∈ R. Левое неравенствовыполняется при1σ>− ,4rа правое – при11σ≥ − .4 4rПоэтому необходимое условие устойчивости запишется следующим образом:1h2σ ≥ − 2 2.(6.14)4 4a τВ частности, если σ ≥ 1/4, то необходимое условие устойчивости выполняется при любых шагах τ и h, а для явной схемы получаем следующеенеобходимое условие устойчивости:aτ≤ 1.h(6.15)6.3.
Обобщенное решение. Классическое (дважды непрерывнодифференцируемое или гладкое) решение уравнения (6.1) существуетне всегда. В таком случае можно искать обобщенное решение [19] этого уравнения – функцию u(x, t), имеющую в области D ограниченныекусочно-непрерывные производные ∂u/∂x, ∂u/∂t и удовлетворяющуюинтегральному соотношениюI∫∫∂u∂udx + a2 dt = −f (x, t)dxdt(6.16)∂t∂xCDдля произвольной ограниченной области D ⊆ D с границей C = ∂D.Если u(x, t) является гладким решением уравнения (6.1), то, интегрируя уравнение (6.1) по области D∫∫ (D2∂2u2∂ u−a∂t2∂x2)∫∫dxdt =f (x, t)dxdt(6.17)Dи применяя формулу Грина для интеграла в левой части этого равенства, получаем соотношение (6.16).
Обратное верно не всегда: еслиu(x, t) – обобщенное решение уравнения (6.1), то оно может и не являться классическим решением.99Известно [25], что линиями разрыва производных ut , ux обобщенногорешения могут быть только характеристики уравнения (6.1)x − at = const,(6.18)x + at = const.При этом величины разрывов удовлетворяют равенствам(ut + aux )лев = (ut + aux )пр ,(ut − aux )лев = (ut − aux )пр ,(6.19)где величины с нижними индексами «лев» и «пр» означают пределыпроизводных при приближении к соответствующей характеристике слева и справа.При построении разностных схем для вычисления обобщенного решения возьмем в качестве кривой C в интегральном соотношении (6.16)контур прямоугольника ωjn = [xj−1/2 , xj+1/2 ] × [tn−1/2 , tn+1/2 ], изображенный штриховой линией на рис. 24.
Тогда (6.16) примет вид:∫∫xj+1/2xj+1/2∂u(x, tn+1/2 )dx −∂txj−1/2t∫xj−1/2n+1/2+∂u(x, tn−1/2 )dx+∂tt∫(6.20)n+1/2a2∂u(xj+1/2 , t)dt −∂xtn−1/2a2∂u(xj−1/2 , t)dt =∂x∫∫f (x, t)dxdt.ωjntn−1/2Применение различных квадратурных формул для вычисления интегралов приводит к различным разностным схемам. В частности, таким способом могут быть получены разностные схемы сквозного счета.n+1nn-1j-1jj+1Рис.
24. Контур интегрирования (штриховая линия), используемыйпри построении консервативной схемы для вычисления обобщенногорешения100При вычислении обобщенного решения по этим схемам линии разрывов(6.18) в процессе расчета не выделяются. Конкретные виды схем будутвыведены в следующем пункте для более общего уравнения колебанийнеоднородной струны.6.4. Задача для уравнения колебаний неоднородной струны.Для неоднородной струны уравнение колебаний принимает следующийвид:()∂u∂2u∂=k(x,t)+ f (x, t), x ∈ (0, l); t ∈ (0, T ],(6.21)∂t2∂x∂xгде k(x, t) ≥ k0 > 0. В области D требуется найти решение этого уравнения, удовлетворяющее начальным условиям (6.2) и краевым условиямтипа (6.3), (6.4) или (6.5), либо комбинированным краевым условиям:условиям разного рода на разных концах отрезка [0, l].
Для уравненияс переменным коэффициентом k(x, t) несколько изменится вид краевыхусловий второго и третьего рода:∂u∂u(0, t) = µ0 (t); k(l, t) (l, t) = µl (t);(6.22)∂x∂x∂u∂uk+ γ0 u + γl u = µl (t).= µ0 (t); k(6.23)∂x∂xx=0x=lДля конструирования разностных схем будем использовать интегро-интерполяционный метод (см. аналогичный метод для параболических уравнений в § 3.3 и для нелинейного уравнения переноса в § 4).Записывая уравнение (6.21) в виде интегрального соотношения, аналогичного (6.16)I∫∫∂u∂udx + k(x, t) dt = −f (x, t)dxdt,(6.24)∂t∂xk(0, t)CDи выбирая в качестве кривой интегрирования замкнутый контур, показанный на рис. 24 штриховой линией, получаем соотношение, аналогичное (6.20) и отличающееся от него только тем, что вместо постоянной a2теперь берется функция k(x, t).
Используя, например, аппроксимацииинтеграловtn+1/2∫k(xj+1/2 , t)unj+1 − unj∂u(xj+1/2 , t)dt ≈ k(xj+1/2 , tn )·τ∂xhtn−1/2101и учитывая равенстваnk(xj+1/2 , tn ) = 0, 5(kjn + kj+1) + O(h2 ),(6.25)приходим к явной схемеun+1− 2unj + un−1jj= Λunj + fjnτ2со следующим оператором Λ:()unj+1 − unjunj − unj−11nΛuj =kj+1/2− kj−1/2,hhh(6.26)(6.27)гдеnnkjn + kj−1kjn + kj+1, kj−1/2 =.22Для той же схемы можно предложить другую форму оператора Λ.Будем вычислять k(xj+1/2 , tn ) по формуле, более точной, чем (6.25):kj+1/2 =( )2 2h∂ k(xj , tn ) + O(h3 ) =2∂x2(6.28)nnnn2 knk−k−2k+khhj+1j−1j+1jj−1= kjn +++ O(h3 ).22h8h2И вновь имеем явную схему (6.26), но теперь уже с оператором()unj+1 − unjunj − unj−11(1)(1)Λ(1) unj =k̄j+1/2− k̃j−1/2,(6.29)hhhh ∂k1k(xj+1/2 , t ) = k(xj , t ) +(xj , tn ) +2 ∂x2nnкоэффициенты которого задаются формулами:(1)k̄j+1/2 =nn−kj−1+ 6kjn + 3kj+1,8(1)k̃j−1/2 =nn−kj+1+ 6kjn + 3kj−1.8Еще один оператор Λ в явной схеме (6.26) может быть получен, если величины k(xj+1/2 , tn ) вычислять как среднее арифметическое двухзначений k, первое из которых находится интерполяцией в точку xj+1/2nвеличин kjn и kj+1, а второе – экстраполяцией в ту же точку величинnkj−1и kjn :()unj+1 − unjunj − unj−11(2)(2)Λ(2) unj =k̄j+1/2− k̃j−1/2,(6.30)hhh102где(2)k̄j+1/2 =nn+ 4kjn + kj+1−kj−1,4(2)k̃j−1/2 =nn−kj+1+ 4kjn + kj−1.4Все три схемы имеют второй порядок аппроксимации, однако явнаясхема (6.26) с оператором Λ из (6.27) является консервативной, а с операторами (6.29) и (6.30) – не является таковой.
Напомним, что для консервативной схемы, в отличие от неконсервативной, дискретный аналогинтегрального соотношения (6.24) выполняется не только для элементарной ячейки сетки ωjn , но и для любой области, полученной объединением произвольного числа элементарных ячеек. Консервативные схемы позволяют рассчитывать решения, когда коэффициент k являетсяразрывной функцией, в то время как неконсервативные для этого негодятся (см.
§ 11, задание 15 лабораторной работы 6).6.5. Случай линейной зависимости правой части от решения. Рассмотрим теперь более общий случай, когда правая часть уравнения (6.21) может линейно зависеть от решения, т. е.f (x, t, u) = a(x, t)u + b(x, t).(6.31)В этом случае для численного решения можно использовать те жесхемы, что и для функции f = f (x, t), причем на этапе реализацииизменения коснутся лишь неявных схем. Ради простоты, рассмотримздесь полностью неявную схему:()n+1un+1− 2unj + ujn−1un+1un+1− un+11j−1jj+1 − ujjn+1n+1=kj+1/2− kj−1/2+τ2hhh+an+1un+1+ bn+1,jjjj = 1, .
. . , N − 1.(6.32)Переписав уравнения (6.32) в виде трехточечныхn+1Aj un+1+ Bj un+1j−1 − Cj ujj+1 = Dj ,j = 1, . . . , N − 1,(6.33)гдеAj =τ 2 n+1k,h2 j−1/2Bj =τ 2 n+1k,h2 j+1/2Cj = 1 + Aj + Bj − τ 2 an+1,jDj = −τ 2 bn+1− 2unj + un−1,jj103и добавив к ним разностные краевые условияun+1= ξ0 un+1+ η0 ,01n+1un+1+ ηN −1 ,N −1 = ξN −1 uN(6.34)решаем полученную систему разностных уравнений (6.33), (6.34) методом прогонки [2]. Для корректности и устойчивости метода прогонкидостаточно потребовать (см. задачу 2.2.4), чтобы выполнялось неравенство(6.35)τ 2 · max a(x, t) ≤ 1.DОтметим, что неравенство (6.35) накладывает ограничение на шагпо времени только в случае a(x, t) > 0. Для этого случая член a(x, t)uв полностью неявной схеме (6.32) лучше брать с n-го слоя по времени,т. е.
заменить an+1un+1на anj unj . Тогда ограничения вида (6.35) на шаг τjjне возникает.6.6. Схемы для нелинейного уравнения колебаний струны.Рассмотрим наиболее общий случай, когда уравнение колебаний струны является нелинейным: коэффициент k и функция правой части fв уравнении (6.21) могут зависеть не только от x и t, но и от искомогорешения u, т. е. f = f (x, t, u), k = k(x, t, u).Для решения нелинейных задач используются консервативные схемы, например, схема с весами с оператором Λ, определенным по формуле (6.27), при этомnkj±1/2=k(xj , tn , unj ) + k(xj±1 , tn , unj±1 ),2fjn = f (xj , tn , unj ).Для явной схемы никаких трудностей в реализации не возникает, таккак оператор Λ действует на известную функцию un .
Таким образом,для вычисления un+1имеем явное выражениеjun+1= 2unj − un−1+ τ 2 (Λunj + fjn ).jj(6.36)При реализации неявных схем следует организовать итерации понелинейности. Например, в случае полностью неявной схемы получимуравнение (6.33) с коэффициентамиAj =τ 2 n+1k,h2 j−1/2Bj =τ 2 n+1k,h2 j+1/2Cj = 1 + Aj + Bj ,(6.37)104Dj = −τ 2 fjn+1 − 2unj + un−1.jНепосредственное применение метода прогонки для решения системы(6.33), (6.34) невозможно, так как уравнение (6.33) является теперьнелинейным, поэтому вместо (6.33) методом прогонки будем решать линейное уравнениеAνj un+1,ν+1− Cjν un+1,ν+1+ Bjν un+1,ν+1= Djνj−1jj−1(6.38)с коэффициентами, вычисленными по предыдущей итерации с номером ν:Aνj =τ 2 n+1,νk,h2 j−1/2Bjν =τ 2 n+1,νk,h2 j+1/2Djν = −τ 2 fjn+1,ν − 2unj + un−1,jn+1,νkj±1/2=Cjν = 1 + Aνj + Bjν ,fjn+1,ν = f (xj , tn+1 , un+1,ν),jk(xj , tn+1 , un+1,ν) + k(xj±1 , tn+1 , un+1,νjj±1 ).2В качестве начального приближения для итерационного процесса можно взять решение с n-го слоя по времени: un+1,0= unj .jПри использовании итераций по нелинейности следует иметь в виду,что коэффициенты краевых условий (6.34) могут зависеть от решения.В таком случае они тоже вычисляются по предыдущей итерации un+1,ν.jЗАДАЧИ6.1.