Исследование А-устойчивости

Исследовать на А-устойчивость одношаговый метод:y m +1 − y mf (x m +1 , y m +1 ) + 3 f ( x m , y m )=.h4ËМодельное уравнение (получающееся при линеаризации):y m +1 − y my + 3 ym= λ m +1,h4или, полагая hλ = z ,(4 − z ) y m+1 − (4 + 3z ) y m = 0 .Ему соответствует характеристическое уравнение (4 − z )q − (4 + 3 z ) = 0 .4 + 3z.Решение характеристического уравнения q =4− zq < 1 , т.е. 4 + 3 z < 4 − z и можно предполагать, что область устойчивости лежит в левойполуплоскости комплексной плоскости z = x + iy .Границей области устойчивости будет множество таких точек z = 4iϕq = 1 , т.е. q = e : z ∗ = 4q −1, для которыхq+3(cosϕ − 1 + i sin ϕ )(cosϕ + 3 − i sin ϕ )e iϕ − 1cos ϕ − 1 + i sin ϕ=4=4iϕcos ϕ + 3 + i sin ϕe +3(cosϕ + 3)2 + sin 2 ϕ2(cos ϕ − 1) + i 4 sin ϕ.(cosϕ + 3)2 + sin 2 ϕТ.к. Re z ∗ ≤ 0 , то вся левая полуплоскость не может являться областью устойчивостисхемы, поэтому метод не А-устойчив.7Замечание.

При z = 1 q = > 1 - схема неустойчива. Область устойчивости внутри3контура.=4Исследовать на А-устойчивость схему четвертого порядка:25 y m − 48 y m −1 + 36 y m − 2 − 16 y m −3 + 3 y m −4= f (xm , y m ) .12hМетод Гира, чисто неявный.ËМодельное уравнение (получающееся при линеаризации):25 y m − 48 y m −1 + 36 y m − 2 − 16 y m −3 + 3 y m − 4 = 12hλy m ,или, полагая hλ = z ,25 y m − 48 y m −1 + 36 y m − 2 − 16 y m −3 + 3 y m − 4 = 12 zy m .Ему соответствует характеристическое уравнение(25 − 12 z )q 4 − 48q 3 + 36q 2 − 16q + 3 = 0 .(1)(2)ŒЧтобыопределить область устойчивости, нам нужно найти множество точек Gкомплексной плоскости z = u + iv = h Re λ + ih Im λ , для которых все корни (2) непревосходят по модулю единицу.Соответственно, границей области G Является множество таких точек z = λh ∈ C ,для которых q < 1 .Для этого выразим параметр z через переменную q :25 − 48q −1 + 36q −2 − 16q −3 + 3q −4.12Заметим, что если q = 1 , то q = e −iϕ , и (3) можно записать в виде:z=(3)25 − 48e iϕ + 36e i 2ϕ − 16e i 3ϕ + 3e i 4ϕz=.(4)12При изменении аргумента ϕ от 0 до 2π точка z описывает замкнутую кривую Γ ,симметричную относительно действительной оси (Im z = Im z ) .Заметим, что:q = 2 соответствует области неустойчивости.

При этом25 − 48 ⋅ 2 −1 + 36 ⋅ 2 −2 − 16 ⋅ 2 −3 + 3 ⋅ 2 −4131z (2) ==.121921соответствует области устойчивости. При этом212347⎛ 1 ⎞ 25 − 48 ⋅ 2 + 36 ⋅ 2 − 16 ⋅ 2 + 3 ⋅ 2z⎜ ⎟ == − . (Т.е. надежда на А1212⎝2⎠устойчивость есть.)25 − 48 + 36 − 16 + 3Re q ∈ [− 1; 1] z (1) == 01,1225 + 48 + 36 + 16 + 3 32,т.е. область устойчивостиz (− 1) ==123представляет собой внешность кривой Γ .q=1 способ (графический).Воспользоваться соответствующей программой и нарисовать кривую Γ .2 способ (аналитический).Перепишем (4) в виде (выделив вещественную и мнимую часть параметра z ):1Т.о. граница области устойчивости проходит, через начало координат.

Этот результат в данной задачеполучен случайно – вообще говоря, надо проводить соответствующее исследование.25 − 48 cos ϕ + 36 cos 2ϕ − 16 cos 3ϕ + 3 cos 4ϕ+12− 48 sin ϕ + 36 sin 2ϕ − 16 sin 3ϕ + 3 sin 4ϕ+i12z=(5)Введем обозначение: x = cos ϕ . Тогда sin ϕ = ± 1 − x 2 .Тогда,учитывая что cos 2ϕ = 2 cos 2 ϕ − 1 = 2 x 2 − 1 ,cos 3ϕ = 4 cos 3 ϕ − 3 cos ϕ = 4 x 3 − 3 x ,cos 4ϕ = 8 cos 4 ϕ − 8 cos 2 ϕ + 1 = 8 x 4 − 8 x 2 + 1 ,sin 2ϕ = 2 sin ϕ cos ϕ = ± 2 x 1 − x 2 ,[]()sin 3ϕ = 3 sin ϕ − 4 sin 3 ϕ = sin ϕ 4 cos 2 ϕ − 1 = ± 4 x 2 − 1 1 − x 2 ,[]()sin 4ϕ = sin ϕ 8 cos 3 ϕ − 4 cos ϕ = ± 8 x 3 − 4 x 1 − x 2 ,получим:()() ()25 − 48 x + 36 2 x 2 − 1 − 16 4 x 3 − 3 x + 3 8 x 4 − 8 x 2 + 1+12.− 48 + 36 ⋅ 2 x − 16 4 x 2 − 1 + 3 8 x 3 − 4 x2±i1− x12Раскрывая скобки, находим:− 8 + 48 x 2 − 64 x 3 + 24 x 4− 32 + 60 x − 64 x 2 + 24 x 3z=±i1− x21212илиz=() ()2i 1− x21 − 6 x 2 + 8 x 3 − 3x 4 ±− 8 + 15 x − 16 x 2 + 6 x 3 .33Заметим, что x ∈ [− 1; 1] .Заметим также, что если Re z ≥ 0 , то метод А-устойчив.z=−()(()(6))⎛ 2⎞Найдем min Re z , т.е.

min ⎜ − 1 − 6 x 2 + 8 x 3 − 3x 4 ⎟ .x∈[−1; 1]⎝ 3⎠2Итак:f (x ) = − 1 − 6 x 2 + 8 x 3 − 3x 4 ,3′⎛ 2⎞⎛ 2⎞f ′( x ) = ⎜ − 1 − 6 x 2 + 8 x 3 − 3 x 4 ⎟ = ⎜ − − 12 x + 24 x 2 − 12 x 3 ⎟ =⎝ 3⎠⎝ 3⎠22= 8 x 1 − 2 x + x = 8 x(1 − x ) .()((f (− 1) = −)()2(1 − 6 − 8 − 3) = 32 ,332f (0) = − ,32f (1) = − (1 − 6 + 8 − 3) = 0 .3min f ( x ) = min{ f (− 1), f (0 ), f (1)} = f (0) = −x∈[−1; 1]метод не А-устойчив.2<03)Исследуем метод на А(α)-устойчивость.2Согласно (6) u = Re z = − 1 − 6 x 2 + 8 x 3 − 3x 4 ,3()1− x2− 8 + 15 x − 16 x 2 + 6 x 3 2.3С одной стороны, уравнение касательной, проходящей через начало координат v = u tg α ,v′с другой стороны, угол наклона касательной определяется соотношением: tg α = vu′ = x .u ′xОтсюда получаем условие касания прямой, проходящей через начало координат, сграницей области устойчивости Γ :v v ′x.(7)=u u ′x(v = Im z =)u ′x = 8 x(1 − x ) ,2v ′x =−x(− 8 + 15x − 16 x2+ 6x3 ) +1− x2(15 − 32 x + 18 x 2 ) =33 1− x21=x (8 − 15 x + 16 x 2 − 6 x 3 ) + (1 − x 2 )(15 − 32 x + 18 x 2 ) =23 1− x18 x − 15 x 2 + 16 x 3 − 6 x 4 + 15 − 32 x + 18 x 2 − 15 x 2 + 32 x 3 − 18 x 4 ==23 1− x5 − 8 x − 4 x 2 + 16 x 3 − 8 x 41.15 − 24 x − 12 x 2 + 48 x 3 − 24 x 4 = ==3 1− x21− x2[][][]Т.о.,1 − x 2 (− 8 + 15 x − 16 x 2 + 6 x 3 ) 5 − 8 x − 4 x 2 + 16 x 3 − 8 x 4.(8)=22⎛ 2⎞234()x−x−x8113⎜ − ⎟(1 − 6 x + 8 x − 3x )⎝ 3⎠Заметим, что (1 − 6 x 2 + 8 x 3 − 3 x 4 ) = 0 , т.е.

(1 − 6 x 2 + 8 x 3 − 3x 4 ) = (1 − x )(a + bx + cx 2 + dx 3 ) ,x =1⎧a = 1,⎪b − a = 0,⎪⎪где ⎨c − b = −6, Т.е. 1 − 6 x 2 + 8 x 3 − 3x 4 = (1 − x ) 1 + x − 5 x 2 + 3x 3 .⎪d − c = 8,⎪⎪⎩− d = −3.((Но 1 + x − 5 x 2 + 3 x 3)x =1)⎧a = 1,⎪b − a = 1,⎪232= 0 и снова 1 + x − 5 x + 3x = (1 − x )(a + bx + cx ), где ⎨⎪c − b = −5,⎪⎩c = −3.Т.е. (1 + x − 5 x 2 + 3 x 3 ) = (1 − x )(1 + 2 x − 3 x 2 ) , 1 + 2 x − 3x 2 = (1 + 3x )(1 − x ) , откудаокончательно 1 − 6 x 2 + 8 x 3 − 3x 4 = (1 + 3x )(1 − x ) , и (8) принимает вид:3−(1 − x )(− 8 + 15x − 16 x22(1 − x )(1 + 3 x )2+ 6x 3) = 5 − 8x − 4 x+ 16 x 3 − 8 x 48x2или2Рассматриваем нижнюю полуплоскость, т.к.

криваяΓ , симметрична относительно действительной оси.− 4 x(1 + x )(− 8 + 15 x − 16 x 2 + 6 x 3 ) = (5 − 8 x − 4 x 2 + 16 x 3 − 8 x 4 )(1 + 3 x ) .Раскрываем скобки:32 x − 60 x 2 + 64 x 3 − 24 x 4 + 32 x 2 − 60 x 3 + 64 x 4 − 24 x 5 =,5 − 8 x − 4 x 2 + 16 x 3 − 8 x 4 + 15 x − 24 x 2 − 12 x 3 + 48 x 4 − 24 x 5или32 x − 28 x 2 + 4 x 3 + 40 x 4 = 5 + 7 x − 28 x 2 + 4 x 3 + 40 x 4 ,т.е.125 x = 5 или x = .53316 ⋅ 64⎛ 2⎞8⎛4⎞⎛ 2 ⎞⎛ 3 ⎞⎛ 1 ⎞u x = 1 = ⎜ − ⎟⎜1 + ⎟⎜1 − ⎟ = ⎜ − ⎟ ⎜ ⎟ = −,625 ⋅ 35⎝ 3⎠5⎝5⎠⎝ 3 ⎠⎝ 5 ⎠⎝ 5 ⎠23⎛⎞⎜ − 8 + 15 1 − 16⎛⎜ 1 ⎞⎟ + 6⎛⎜ 1 ⎞⎟ ⎟24(− 5 ⋅ 125 − 16 ⋅ 5 + 6)⎜⎟555⎝⎠⎝⎠⎝⎠5v x= 1 ===3 ⋅ 125352 6 (− 5 ⋅ 141 + 6)2 6 ⋅ 699.== −3 ⋅ 6253 ⋅ 6256 ⋅ 699v2 6 ⋅ 699tg α ==≈ 3.344 ≈ 73.35 o=u x= 116 ⋅ 64512⎛1⎞1− ⎜ ⎟⎝5⎠25метод А(α)-устойчив, α ≈ 73.35 o.