semopt13 (Практические занятия)

Занятие 13. ЧИСЛЕННЫЕ МЕТОДЫ РЕШЕНИЯ СИСТЕМОБЫКНОВЕННЫХ ДИФФЕРЕНЦИАЛЬНЫХ УРАВНЕНИЙПример 1. Приближенно решить задачу Кошиy ′ = −2 y − 3 x + 2 , y (0) = 0на отрезке [0; 1] с шагом h = 0,1 различными методами: Эйлера (явным и неявным),предсказания и коррекции второго порядка, Рунге–Кутты четвертого порядка иметодом трапеций.† Поскольку явный метод Эйлера и метод Рунге–Кутты четвертого порядкаотносятся к классу ограниченно устойчивых, то для них требуется определитьвеличины критического шага.

Сравнивая данное уравнение с тестовым примером,2заметим, что μ = −2 . Тогда для явного метода Эйлера hкр = − = 1 , а для методаμ2,78Рунге–Кутты hкр = −= 1,39 . Следовательно, интегрирование с шагом h = 0,1 < hкрμобеспечивает устойчивость этих методов.Явный метод Эйлера.

Из общей формулыyˆi +1 = yˆi + hi +1 f ( x i , yˆi ) ,i = 0, n − 1 , yˆ0 = y 0 ;получаем расчетную формулу явного метода Эйлера:yˆi +1 = yˆi + 0,1 ⋅ (−2 yˆi − 3 x i + 2) , ŷ 0 = 0 , i = 0,9 .Метод предсказания и коррекции второго порядка.Шаг «предиктор»:yˆi(+П1) = yˆi + hi +1 f ( x i , yˆi ) ;Шаг «корректор»:)ˆi +yˆ i +1 ≡ yˆ i(К+1 = yhi +1[ f ( x i , yˆ i ) + f ( x i + hi +1 , yˆ i(П)+1 )] .2Отсюда следуют расчетные формулы метода предсказания и коррекции:yˆi(+П)1 = yˆi + 0,1 ⋅ (−2 yˆi − 3 x i + 2) , ŷ 0 = 0 , i = 0,9 ;0,1yˆi +1 = yˆi +(−2 yˆi − 3x i + 2 − 2 yˆi(+П)1 − 3x i +1 + 2) ; ŷ 0 = 0 , i = 0,9 .2Метод Рунге–Кутты.

Из общей формулыyˆi +1 = yˆi +hi +16(K 1,i + 2K 2,i + 2K 3,i + K 4,i ),249yˆ0 = y 0 , i = 0, n − 1 ,гдеhh⎛⎞K 2,i = f ⎜⎜ x i + i +1 , yˆi + i +1 K 1,i ⎟⎟ ,22⎝⎠hh⎛⎞= f ⎜⎜ x i + i +1 , yˆi + i +1 K 2,i ⎟⎟ , K 4,i = f x i + hi +1 , yˆi + hi +1 ⋅ K 3,i ,22⎝⎠K 1,i = f i = f ( x i , yˆi ),K 3,i()следуют формулы метода Рунге–Кутты четвертого порядка:0,1yˆi +1 = yˆi +(K 1,i + 2K 2,i + 2K 3,i + K 4,i ) ,6K 1,i = −2 yˆi − 3 x i + 2 ,K 2,i = −2 ⋅ ( yˆi + 0,05 ⋅ K 1,i ) − 3 ⋅ ( x i + 0,05) + 2 ;K 3,i = −2 ⋅ (yˆi + 0,05 ⋅ K 2,i ) − 3 ⋅ (xi + 0,05) + 2 ; K 4,i = −2 ⋅ ( xi + 0,1 ⋅ K 3,i ) − 3 ⋅ (xi + 0,1) + 2 .Неявный метод Эйлера. Из общей формулыyˆi +1 = yˆi + hi +1 f ( x i +1 , yˆi +1 ) ≡ Φ( x i , x i +1 , yˆi +1 ) , i = 0, n − 1 ,получаем расчетную формулу неявного метода Эйлера:yˆi +1 = yˆi + 0,1 ⋅ (−2 yˆi +1 − 3 x i +1 + 2) ,yˆi +1 =откудаМетод трапеций.yˆi − 0,3x i +1 + 0,2.1,2Из общей формулыyˆi +1 = yˆi +hi +12[ f i + f (xi +1 , yˆi +1 )] ≡ Φ(xi , xi +1 , yˆi +1 ) , i = 0, n − 1 ,получаем расчетную формулу метода трапеций:откуда0,1yˆi +1 = yˆi +⋅ (−2 yˆi − 3x i + 2 − 2 yˆi +1 − 3x i +1 + 2) ,2yˆi +1330,9 ⋅ yˆi − ⋅ 0,1 ⋅ x i + 2 ⋅ 0,1 − ⋅ 0,1 ⋅ x i +122.=1,1Очевидно, в данном примере удалось получить явные формулы для нахожденияŷ i +1 неявным методом Эйлера и методом трапеций лишь в силу линейностирешаемого уравнения.

В общем случае применяются методы простых итераций илиНьютона.Точное решение рассматриваемой задачи Коши: y( x ) = 1,75 − 1,5x − 1,75 e −2 x .Результаты расчетов приведены в табл. 1, в последней строке которой указаныфактические погрешности.Анализ результатов показывает, что при решении данной задачи методпредсказания и коррекции точнее явного и неявного методов Эйлера, но уступаетметоду Рунге–Кутты и совсем немного методу трапеций (порядок погрешностиодинаков).„250Таблица 1ЯвныйметодЭйлераНеявныйметодЭйлераМетодРунге–КуттыМетодпредсказанияи коррекцииМетодтрапецийy (x )0,00,0000,0000,0000,0000,0000,0000,10,2000,1420,1670,1650,1680,1670,20,3300,2380,2770,2730,2790,2770,30,4040,2870,3400,3350,3420,3400,40,4330,3060,3640,3590,3660,3640,50,4270,2970,3560,3510,3580,3560,60,3910,2640,3230,3180,3250,3230,70,3330,2120,2680,2640,2700,2680,80,2560,1430,1970,1920,1990,1970,90,1650,0610,1110,1070,1120,1111,00,062-0,0330,0130,0090,0150,013max εi0,0680,0590,00001040,0050,002xiПример 2.

Найти приближенное решение задачи Кошиy′ = z − 1 ,y (0) = 1 ,z ′ = − y − 2z , z (0) = −1на отрезке [0; 1] с шагом h = 0,1 методами Эйлера (неявным и модифицированным),Адамса–Бэшфорта третьего порядка и трапеций.† Путем прямой подстановки в систему легко убедиться в том, что точноерешение задачи имеет вид y ( x) = −2 + 3 e − x + x e − x ; z ( x) = 1 − 2 e − x − x e − x .Выписываем формулы для нахождения приближенного решения указаннымиметодами (при этом применяется векторная форма записи).Для неявного метода Эйлераyˆi +1 = yˆi + hi +1 f ( x i +1 , yˆi +1 ) ≡ Φ( x i , x i +1 , yˆi +1 ) , i = 0, n − 1 ,имеем⎛ yˆi +1 ⎞ ⎛ yˆi ⎞⎛⎜⎟ ⎜⎟⎜⎜ zˆ ⎟ = ⎜ z ⎟ + 0,1 ⋅ ⎜ −⎝ i +1 ⎠ ⎝ i ⎠⎝zˆi +1 − 1 ⎞ ⎛ yˆi + 0,1 ⋅ zˆi +1 − 0,1 ⎞⎟=⎜⎟,yˆi +1 − 2zˆi +1 ⎟⎠ ⎜⎝ zˆi − 0,1 ⋅ yˆi +1 − 0,2 ⋅ zˆi +1 ⎟⎠yˆ0 = y (0) = 1,zˆ0 = z (0) = −1.Разрешая эту систему относительно ŷi +1 , ẑ i +1 , окончательно получаемyˆi +1 =0,1 ⋅ (zˆi − 0,1 ⋅ yˆi + 0,001)− 0,1 ;1,21251zˆi +1 =zˆi − 0,1 ⋅ yˆi + 0,01.1,21Соотношенияyˆi+12= yˆi +hi +1f ( x i , yˆi ) ,2i = 0, n − 1 ,h⎛⎞yˆi +1 = yˆi + hi +1 f ⎜⎜ x i + i +1 , yˆi + 1 ⎟⎟ ,22 ⎠⎝i = 0, n − 1 ,для модифицированного метода Эйлера принимают вид0,10,1yˆ i + 1 = yˆ i +( zˆ i − 1); zˆ i + 1 = z i +(− yˆ i − 2 zˆ i );2222()yˆ i +1 = yˆ i + 0,1 ⋅ zˆ i + 1 − 1 ;2()zˆ i +1 = zˆ i + 0,1 ⋅ − yˆ i + 1 − 2 zˆ i + 1 .22Для метода Адамса–Бэшфорта третьего порядка изhyˆi +1 = yˆi +[23 f i − 16 f i −1 + 5 f i − 2 ] , i = 2, n − 1 ,12находим0,1[23 ⋅ ( zˆ i − 1) − 16 ⋅ ( zˆ i −1 − 1) + 5 ⋅ ( zˆ i −2 − 1)] ,120,1[23 ⋅ (− yˆ i − 2 zˆ i ) − 16 ⋅ ( − yˆ i −1 − 2 zˆ i −1 ) + 5 ⋅ (− yˆ i −2 − 2 zˆ i −2 )] .= zˆ i +12yˆ i +1 = yˆ i +zˆ i +1Для определения «разгонных» точек ( x 0 , yˆ0 , zˆ0 ), ( x1 , yˆ1 , zˆ1 ), ( x 2 , yˆ2 , zˆ2 ) воспользуемсямодифицированным методом Эйлера.

Точка ( x 0 , yˆ0 , zˆ0 ) определяется начальнымиусловиями (0; 1; − 1) .Выпишем соотношения для метода трапецийyˆi +1 = yˆi +hi +12[ f i + f (xi +1 , yˆi +1 )] ≡ Φ(xi , xi +1 , yˆi +1 ) , i = 0, n − 1 ,и разрешим их относительно неизвестных:⎛ yˆ i +1 ⎞ ⎛ yˆ i ⎞ 0,1 ⎡⎛ zˆ i − 1 ⎞ ⎛ zˆ i +1 − 1 ⎞ ⎤⎢⎜⎜ ˆ ⎟=⎜ ˆ ⎟+⎟+⎜⎟⎥ =⎝ z i +1 ⎠ ⎝ z i ⎠ 2 ⎣⎝ − yˆ i − 2 zˆ i ⎠ ⎝ − yˆ i +1 − 2 zˆ i +1 ⎠ ⎦yˆ i + 0,05 zˆ i + 0,05 zˆ i +1 − 0,1⎛⎞=⎜⎟,⎝ −0,05 yˆ i − 0,05 yˆ i +1 + 0,9 zˆ i − 0,1zˆ i +1 ⎠откуда путем разрешения системы относительно yˆi +1 , zˆi +1 получаемyˆ i +1 = yˆ i + 0,05 zˆ i + 0, 05zˆ i +1 =−0,1 yˆ i + 0,8975 zˆ i + 0, 051,1025−0,1 yˆ i + 0,8975 zˆ i + 0,051,1025− 0,1 ;.Результаты проведенных расчетов даны в табл. 2 дляyˆi , y ( x ) и табл.

3 дляzˆi , z ( x ) соответственно. Из их анализа вытекает, что наиболее точно величину ŷiможно рассчитать методом Адамса–Бэшфорта, а величину ẑ i – методом трапеций. „252Таблица 2xiНеявный методЭйлераМодиф. методЭйлераМетод Адамса–БэшфортаМетодтрапецийТочноерешение0,01,000001,000001,000001,000001,000000,10,809920,805000,804990,804990,805000,20,629600,619970,619910,619910,619940,30,458850,444790,444690,444640,444700,40,297400,279240,279080,2790000,279090,50,145000,123090,122850,122730,122860,60,00131-0,02397-0,02428-0,02444-0,024280,7-0,13397-0,16224-0,16265-0,16284-0,162630,8-0,26119-0,29208-0,29257-0,29279-0,292550,9-0,38072-0,41384-0,41441-0,41465-0,414381,0-0,49287-0,52787-0,52853-0,52879-0,52848max ε i0,037204150,000672910,000059420,00033550Таблица 3xiНеявный методЭйлераМодиф.

методЭйлераМетод Адамса–БэшфортаМетодтрапецийТочноерешение0,0-1,00000-1,00000-1,00000-1,00000-1,000000,1-0,90083-0,90000-0,90023-0,90023-0,900160,2-0,80316-0,80095-0,80132-0,80132-0,801210,3-0,70753-0,70357-0,70402-0,70401-0,703880,4-0,61440-0,60844-0,60893-0,60890-0,608770,5-0,52408-0,51601-0,51650-0,51645-0,516320,6-0,43684-0,42663-0,42709-0,42702-0,426910,7-0,35287-0,34054-0,34095-0,34087-0,340780,8-0,27229-0,25794-0,25828-0,25818-0,258120,9-0,19518-0,17894-0,17920-0,17908-0179051,0-0,12158-0,10357-0,10378-0,10364-0,10364max ε i0,019581330,000325860,000121450,00003005-253.