1612725465-542b3179b36a4849700e0b2ecf7da111 (828844), страница 14
Текст из файла (страница 14)
рисунки 19, а и 22, а), но хуже, чем при использовании схемыЛакса – Вендроффа на адаптивной сетке (ср. рисунки 21, а и 22, а).На рис. 23 изображено поведение численного решения, полученногопо схеме предиктор-корректор с параметром [29]0g̃j+1/2 − g̃j+1/2−sθj+1/2 = θ0,j+1/2g̃j+1/2θ0,j+1/2припри g̃j+1/2 ≤ g̃j+1/2−s ,g̃j+1/2 · g̃j+1/2−s ≥ 0, g̃j+1/2 > g̃j+1/2−s ,(5.26)g̃j+1/2 · g̃j+1/2−s ≥ 0,приg̃j+1/2 · g̃j+1/2−s < 0,где() ng̃j+1/2 = |āj+1/2 | 1 − Crj+1/2 vq,j+1/2,θ0,j+1/2 =1Crj+1/2Crj+1/2 =−1,æ|āj+1/2 |<1,nJj+1/2(5.27))( nnā = a−xt , s = sgn(āj+1/2 ), æ = τ /h, vq,j+1/2= vj+1− vjn /h.
Сравнивая рисунки 23, а и 20, а, видим, что численное решение, полученное наадаптивной сетке, также не имеет осцилляций, но положение разрывапередается лучше, чем на равномерной сетке.90u1021.0t8160.5420.001020300x010а20x30бРис. 21. а – графики точного (1) и численного (2) решений в момент времени t = 10; б – траектории узлов адаптивной сетки. СхемаЛакса – Вендроффа, N = 60u10t11.0860.54220.001020300xа01020x30бРис. 22.
а – графики точного (1) и численного (2) решений в моментвремени t = 10; б – траектории узлов адаптивной сетки. Противопоточная схема, N = 605.5. Противопоточная схема на адаптивной сетке. Мы уже отмечали, что схема предиктор-корректор (3.3), (3.4) с параметром θ = θ0 ,заданным по формуле (3.5), совпадает с противопоточной схемой (3.7)на равномерной сетке. Получим теперь противопоточную схему на адаптивной сетке, взяв в схеме (5.14), (5.15) параметр θj+1/2 = θ0,j+1/2 . Длятакого параметра θ уравнение (5.14) и формула (5.27) дают следующее91∗выражение для величин vj+1/2:)][()( n∗n− vjn .vj+1/2= 0.5 vjn + vj+1− sgn āj+1/2 vj+1∗Подставляя найденное значение vj+1/2в уравнение шага корректор (5.15),приходим к следующей схеме:()nn(Jv)n+1− (Jv)nj+ vjnvjn + vj+1vj−11j+āj+1/2− āj−1/2+τh22(5.28)|āj−1/2 | n|āj+1/2 | n+vq,j−1/2 −vq,j+1/2 = 0.22u10t11.0860.54220.001020300x0а1020бx30Рис.
23. а – графики точного (1) и численного (2) решений в момент времени t = 10; б – траектории узлов адаптивной сетки. Схемапредиктор-корректор (5.14), (5.15), (5.26), N = 60Такая форма записи не совпадает с привычным для нас видом противопоточных схем. Поэтому преобразуем разностное уравнение (5.28)к другому виду, используя приведенное ниже утверждение.Лемма 5.1. Если для вычисления скорости узлов подвижной сетки и якобиана отображения используются формулы (5.16), (5.17), тосправедливо тождествоJjn+1 − Jjn− (xt )q,j = 0.τ(5.29)Д о к а з а т е л ь с т в о. С учетом формул (5.16), (5.17) получаем()n+1n+1nnJj+1/2+ Jj−1/2− Jj+1/2+ Jj−1/2Jjn+1 − Jjn==τ2τ92( n)n+1nnnxn+1+ xn+1− xn+1j+1 − xjjj−1 − xj+1 − xj + xj − xj−1==2τ hxt,j+1/2 − xt,j−1/2xt,j+1 − xt,j + xt,j − xt,j−1=== (xt )q,j ,2hhчто и требовалось доказать.Используя формулу (5.29) получаем, чтоJjn = Jjn+1 − τxt,j+1/2 − xt,j−1/2,hпоэтому(Jv)n+1− (Jv)njvjn+1 − vjnxt,j+1/2 − xt,j−1/2 nj= Jjn+1+vj .ττhСледовательно, схему (5.28) можно переписать так:[ nnnvjn+1 − vjnvj + vj+1vjn + vj+11Jjn+1 ·+a− xt,j+1/2+ xt,j+1/2 vjn −τh22( n)]nvj−1 + vjnvj−1+ vjnn− a− xt,j−1/2+ xt,j−1/2 vj+(5.30)22|āj−1/2 | n|āj+1/2 | nvq,j−1/2 −vq,j+1/2 = 0.22В силу очевидных равенств[ n]nn)vj + vj+1vj−1+ vjna( n1na−a=vq,j+1/2 + vq,j−1/2,h222[]nvjn + vj+1xt,j+1/2 n1nxt,j+1/2− xt,j+1/2 vj =vq,j+1/2 ,h22[]n+ vjnvj−1xt,j−1/2 n1xt,j−1/2− xt,j−1/2 vjn = −vq,j−1/2h22+и с учетом обозначения a = a − xt , схема (5.28) может быть записанав следующем окончательном виде:vjn+1 − vjnāj−1/2 + āj−1/2 n+vq,j−1/2 +τ2Jjn+1(5.31)āj+1/2 − āj+1/2 n+vq,j+1/2 = 0.2Jjn+193Видим, что противопоточная схема (5.31) на подвижной сетке имееттакой же вид, как на равномерной сетке противопоточная схема (3.7)для уравнения переноса с постоянным коэффициентом или противопоточная схема (4.40) для нелинейного скалярного уравнения (4.1).
Аналогичными будут и свойства этих схем. В частности, как и на равномерной сетке, противопоточная схема на подвижной сетке сохраняетмонотонность численного решения.Теорема 5.1. Выполнение при всех j условийæ āj+1/2 {}(5.32)n+1 ≤ 1,min Jjn+1 , Jj+1достаточно для того, чтобы противопоточная схема (5.31) сохраняламонотонность численного решения.Д о к а з а т е л ь с т в о.
Противопоточная схема (5.31) может бытьзаписана в виде схемы (2.22) из условия теоремы 2.3:vjn+1 − vjn−+nn+ Cj−1/2vq,j−1/2− Cj+1/2vq,j+1/2= 0,τгде−Cj−1/2=āj−1/2 + āj−1/22Jjn+1+Cj+1/2,=āj+1/2 − āj+1/22Jjn+1.Оба коэффициента неотрицательны, а в силу условия (5.32) и очевидных неравенств1Jjn+1≤1{ n+1 n+1 } ,min Jj , Jj+1получаем, что−Cj+1/2++Cj+1/21n+1Jj+1≤1{ n+1 n+1 } ,min Jj , Jj+1āj+1/2 + āj+1/2 + āj+1/2 − āj+1/21{ n+1 n+1 }≤≤ .æ2 min Jj , Jj+1Таким образом, условия (2.23) теоремы 2.3 выполнены, поэтому противопоточная схема (5.31) обладает свойством сохранения монотонностичисленного решения.Замечание. Схема предиктор-корректор (5.14), (5.15) с переменным параметром (5.26) также сохраняет монотонность численного решения при использовании подвижных сеток [29].94§ 6.
Разностные схемы для уравненияколебаний струныВ предыдущих параграфах мы познакомились с конечно-разностными схемами, аппроксимирующими гиперболические уравнения первого порядка – линейные и нелинейные уравнения переноса. В настоящем параграфе мы рассмотрим схемы для гиперболического уравнениявторого порядка, взяв в качестве примера уравнение колебаний струны.6.1. Математическая постановка начально-краевой задачи длягиперболического уравнения колебаний однородной струны заключается в следующем: требуется найти функцию u(x, t), определенную в замкнутой области D = [0, l] × [0, T ], являющуюся решением уравнения2∂2u2∂ u=a+ f (x, t), x ∈ (0, l),∂t2∂x2и удовлетворяющую заданным начальнымu(x, 0) = u0 (x),ut (x, 0) = v0 (x),t ∈ (0, T ](6.1)x ∈ [0, l](6.2)и краевымu(0, t) = µ0 (t),u(l, t) = µl (t),t ∈ [0, T ](6.3)условиям. Здесь a = const > 0 – скорость распространения колебаний,f – внешняя нагрузка на струну, концы струны движутся по заданнымзаконам (6.3).
В момент времени t = 0 заданы начальное отклонениеструны u0 (x) и скорости ее точек v0 (x).Краевые условия (6.3) называются условиями первого рода. Краевыеусловия второго рода заключаются в задании на концах отрезка [0, l]значений производных от решения∂u(0, t) = µ0 (t),∂x∂u(l, t) = µl (t).∂xУсловия третьего рода выглядят следующим образом:∂u∂u+ γ0 u + γl u = µl (t).= µ0 (t),∂x∂xx=0x=l(6.4)(6.5)Возможно также задание на концах отрезка [0, l] условий разного рода.Например, при x = 0 может быть задано условие первого рода, а приx = l – второго.95Физическая интерпретация краевых условий (6.3)–(6.5) дана в работе [25].
Там же доказана корректность рассмотренных постановок задач. В частности, показано, что для третьей краевой задачи необходимовыполнение условийγ0 ≤ 0, γl ≥ 0.(6.6)Решение задачи Коши (уравнение (6.1) и начальные условия (6.2)рассматриваются в бесконечной области −∞ < x < ∞, граничные условия не задаются) выписывается явно [25]:1u0 (x + at) + u0 (x − at)+u(x, t) =22aa+2∫t (x+at∫v0 (ζ)dζ+x−atx+a(t−η)∫)(6.7)f (ζ, η)dζ dη.0x−a(t−η)6.2. Схема с весами. Введем в области D равномерную прямоугольную сетку (xj , tn ) (j = 0, .
. . , N , n = 0, . . . , M ) с шагами h = l/N ,τ = T /M и рассмотрим на ней трехслойную схему с весами для решенияпервой начально-краевой задачи (6.1)–(6.3):[]unt̄t,j = a2 Λ σun+1+ (1 − 2σ) unj + σun−1+ φnj ,jjn = 1, . . . , M − 1, j = 1, . . . , N − 1,(6.8)nnu0 = µ0 (t ), unN = µl (tn ), n = 0, . . . , M,u0j = u0 (xj ),u0t,j = ηj ,j = 0, . . . , N,где σ – произвольный вещественный параметр (вес схемы),unt̄t,j =un+1− 2unj + ujn−1j;τ2Λunj = unx̄x,j =u0t,j =unj+1 − 2unj + unj−1;h2u1j − u0j.τКраевые условия и первое начальное условие аппроксимированыточно.
Правая часть f (x, t) аппроксимируется некоторой сеточной функцией φnj , которая может быть и не равной f (xj , tn ), и, наконец, второеначальное условие ut (x, 0) = v0 (x) будем аппроксимировать разностным96выражением u0t,j = ηj так, чтобы погрешность аппроксимации была величиной порядка O(τ 2 ). Для этого достаточно положитьηj = v0 (xj ) +)τ ( 2 u0 (xj+1 ) − 2u0 (xj ) + u0 (xj−1 )a+f(x,0).j2h2(6.9)При σ = 0 схема называется явной. При σ ̸= 0 получаем неявнуюсхему, в которой на (n + 1)-м временно́м слое связаны три неизвестныхn+1значения un+1и un+1j−1 , ujj+1 .Выясним, как влияет выбор веса σ на порядок аппроксимации дифференциального уравнения разностным.
Пусть n ≥ 1 и 0 < j < N .Тогда()nψjn = Lh (u)h − fh j ==u(xj , tn+1 ) − 2u(xj , tn ) + u(xj , tn−1 )(σ)− a2 Λuj − φnj ,τ2(σ)где u – достаточно гладкое решение задачи (6.1)–(6.3), uj(6.10)= u(σ) (xj ),u(σ) (x) = σu(x, tn+1 ) + (1 − 2σ)u(x, tn ) + σu(x, tn−1 ).Оценим порядок погрешности аппроксимации в узле (xj , tn ). Применяя формулу Тейлора, получаем( )u(xj , tn+1 ) − 2u(xj , tn ) + u(xj , tn−1 )τ2n=u(x,t)+utttt (xj , tn )+O τ 4 ;ttj2τ12( )u(σ) (xj ) = u(xj , tn ) + τ 2 σutt (xj , tn ) + O τ 4 ;(σ)Λujh2uxxxx (xj , tn )+12()τ 2 σh2+uttxxxx (xj , tn ) + O τ 4 + h4 .12= uxx (xj , tn ) + τ 2 σuttxx (xj , tn ) +Следовательно,ψjn = utt +()()τ2h2utttt − a2 uxx + στ 2 uttxx + uxxxx − φnj + O τ 4 + h4 .1212Учитывая, что для гладкого решения задачи (6.1)–(6.3) выполняются равенстваutt (x, t) = a2 uxx (x, t) + f (x, t);97uttxx (x, t) = a2 uxxxx (x, t) + fxx (x, t);utttt (x, t) = a4 uxxxx (x, t) + a2 fxx (x, t) + ftt (x, t),перепишем выражение (6.10) для погрешности аппроксимации в следующем виде:[]1h2ψjn = −τ 2 a4 σ −+uxxxx + f − φnj −12 12a2 τ 2(6.11)()( 4)1τ22 24−a τ σ −fxx + ftt + O τ + h .1212В этом равенстве функции u и f , а также их производные вычисляютсяв одной и той же точке (xj , tn ).Из выражения (6.11) следует, что если φnj = f (xj , tn ), то ψjn =( 2)O τ + h2 при любом σ.
Для схемы с весомσ ≡ σ∗ =1h2−12 12a2 τ 2(6.12)и правой частьюφnj = f (xj , tn ) +h2τ2fxx (xj , tn ) + ftt (xj , tn )1212(6.13)()получается схема (6.8) повышенного порядка аппроксимации O τ 4 + h4 .Для исследования устойчивости схемы с весами применим спектральный признак. Множитель перехода этой схемы удовлетворяет уравнениюλ2 − 2(1 − q)λ + 1 = 0,гдеq=2r sin2 φ21 + 4σr sin2φ2;r=a2 τ 2.h2Дискриминант квадратного уравнения равен d = q(q−2). Если d ≤ 0,то корни по модулю равны единице, а если d > 0, то один из вещественных корней по модулю больше единицы.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
