VIII Агафонов и др. Дифференциальные уравнения (1081404), страница 42
Текст из файла (страница 42)
Вместе с тем известно [Ч1], что более точными являются квадратурные 4ормулы трапеции и прямоугольника. Но если применить, например, формулу трапеции непосредственно к (12.13), то получим, вообще говоря, нелинейное уравнение Ь х„ь1 = х„+ — (~„+ ~($„.ьм х„ь1)) (12.14) относительно искомого значения х„~.ь Чтобы избежать необходимости решать такое уравнение, вычисления проводят в два этапа. Сначала прогнозируют при помощи (12.11) ожидаемое значение х„+м а затем уже используют формулу трапеции длякоррекции этого значения. В итоге рабочая формула принимает вид Ь х„~.1 = х„+ — (~„+ ~(1„.ьм х„+ Ц„)). 2 (12.15) хи+1 = хи+ П~ +Из> х .ь1!з)Ь, Она характеризует метпод Эй 1ера — Хаши.
Аналогичная идея лежит в основе усовершенспьвованного мепьода ломаных. В этом методе сначала прогнозируют ожидаемое значение х„+1~э — — х(1„+1~э) = х„+ Т"„Ь/2 в промежуточной точке 1„+1~э — — Ф„+ Ь/2, а затем используют формулу прямоугольника 313 12,З.
Метод ломаных Эйлера так что рабочая формула имеет вид х„+1 = х„+ ~(1н+~12, х„+ У„Ь/2)Ь. (12.16) Пример. Применим метод ломаных Эйлера и его модификации для численного решения задачи Коши Ых — =21х, х(0) =1, ФЕ[0, 1], (12.17) имеющей точное аналитическое решение х(е) = е (12.18) В данном случае благодаря линейности ОДУ (12.17) удобно использовать и (12.14), принимающее вид Ь2 х„+1 = х„— а 1-Ь2( +1)' Результаты расчетов при Ь = 0,1 (Ф = 10) представлены в табл. 12.1. Таблица 19.1 Рабочая формула (12.14) (12.18) (12.11) (12.15) (12.16) Из таблицы видно, что модификации позволяют уменьшить погрешность вычислений метода ломаных Эйлера.
0,1 0,2 0,3 0,4 0,5 0,6 0,7 0,8 0,9 1,0 1,00000 1,02000 1,06080 1,12445 1,21440 1,33584 1,49615 1,70561 1,97850 2,33463 1,01010 1,04102 1,09468 1,17450 1,28577 1,43624 1,63700 1,90390 2,25958 2,73660 1,01000 1,04070 1,09399 1,17319 1,28347 1,43236 1,63059 1,89345 2,24260 2,70906 1,01000 1,04060 1,09367 1,17253 1,28228 1,43038 1,62749 1,88870 2,23546 2,69843 1,01005 1,04081 1,09417 1,17351 1,28403 1,43333 1,63232 1,89648 2,24791 2,71828 314 12. привлиженнык мктоды ркшкния оду 12.4. Метод Рунге — Кутты Построение метода численного решения задачи (12.7) Коши для обыкновенного ди44еренииального уравнения (ОДУ) первого порядка можно провести формальным путем. Рассмотрим один из частпичных отрезков [Фп, оп+1) (и = О, 1У вЂ” 1) разбиения отрезка [1в, Т) на Ф равных частей с шагом сетки Ь = (Т вЂ” зв)/№ Предположим, что известна ордината хп = = х(еп) точки искомой интегральной кривой х(1) для задачи (12.7) Коши в узле сетки З„=зв+пЬ, являющемся левым кон- цОМ ЧаетИЧНОГО ОтрЕЗКа [Фп, оп+1).
БудЕМ ИСКатЬ ПрИбЛИжЕННОЕ значение ординаты, соответствующей правому концу этого отРЕЗКа В УЗЛЕ оп+1, В ВИДЕ хп+1 = х(2п 11) = хп+ Ь~) 7;Й;, ~~~ УЬ = 1, (12.19) где Й1 1(~п~ хп)1 Й2 = 2 (Фп + оз Ь, хп + ~321 ЬЙ1); Йз = 1(зп+озЬ, хп+Ь(Р21Й1+РззЙ2)); (12.20) тп-1 Й~п = 1(Зп+отпЬ, хп+ Ь ~~' 17пОЙ1). Число тЕ14 и коэффициенты а;, Д1,7ь (1 =2, т, у =1,т — 1, Й = 1, т ) выбирают, руководствуясь различнь1ми соображениями, но одно из основных соображений связано с желанием повысить точность рабочей формулы (12.19). Эта формула характеризует метод, названный по имени немецких математиков К.Д.Т. Рунге (1856 — 1927) и В.
Кутты (1867-1944) мепзодом Рунге — Кутпьы. Пусть х(1) — точное решение задачи (12.7) Коши, а хп+1— приближенное значение ординаты интегральной кривой для 12А. Метод Рунге — Кутты ОДУ в (12.7), проходящей через точку ($в,хв), найденное из (12.19) в узле Ф„е1 при условии, что х„= х(1„), т.е. по точному значению ординаты этой кривой при 1„. Ясно, что в этом случае разность р(Ь) = х(1х+1) - хи+1 = =х(1„+Ь) — х(Ф„) — Ь~~~~ у;Ь1~, (12.21) 1=1 хх=хрх) называемая погрешностью метпода на шаге, будет функцией шага Ь.
Предположим, что функция 1(8, х) в (12.7) дифференцируема в+ 1 раз в области изменения своих аргументов. Тогда, согласно (12.20) и (12.19), столько же раз можно продифференцироватьпо Ь ифункцию <р(Ь). Представимеев окрестности точки Ь = 0 формулой Тейлора с остаточным членом в форме Лагранхеа ~р(Ь) = ~~) Ь + Ь'+1, 6 Е (О, 1). (12.22) ~") 0 ~'+1) 6Ь р к) (8+ 1)1 В этой записи принято, что при Ь = 0 ~р~о)(0) = <р(0) и О! = 1, а также обозначено ~р~")(0) = Н"<р/аЬх~, Задачу подбора коэффициентов а;, )УО, уь (1 =2,т, 1 = = 1, т — 1, Ь = 1, т), обеспечивающего повышение точности рабочей формулы (12.19) теперь можно сформулировать так: при заданном т Е М найти такое сочетание значений этих коэффициентов, чтобы при наибольшем возможном в е М в (12.19) было выполнено условие <р(0) = у'(0) = ...
= ~р(') (0) = О. В этом случае погрешность метода на шаге будет определять в (12.22) лишь остаточный член, пропорциональный сомножителю Ь'+1. Степень в + 1 этого сомножителя называют порядком точностпи метпода на шаге. Поскольку погрешность возникает на каждом шаге, то при 1т' вычислениях на 316 нь ИРиелиженные метОДы РешениЯ ОДУ отрезке [$е, Т) суммарная погрешность будет пропорциональна ЖЬ'+1=(Т вЂ” $е)Ь', т.е. порядок л точностпи метода ка огпреэке будет на единицу меньше, чем на шаге. Ясно, что при тп = 1 (12.19) переходит в (12.11), т.е. метод ломаных Эйлера можно рассматривать как частный случай метода Рунге — Кутты. В этом случае, согласно (12.19)— (12.21), имеем у1 = 1 и р(Ь) = х(1„+ Ь) - х($„) — ЬУ($„, х(1„)).
Отсюда при Ь = 0 получаем у(0) = 0 и, учитывая, что Их/й1,, =Па„,х(ф„)), 1~'(0) = — ~ = [ — ! 1 — ~(й„, х(й„))] = О, а также ~л(Ь) = Р = ~*~ 1=9"(~„+Ь). с Ьз йз!с=с„+и Поскольку в общем случае р" (0) = х" (Ф„) ф О, то в соответствии с (12.22) устанавливаем, что метод ломаных Эйлера имеет второй порядок точности на шаге и первый порядок точности на отрезке, причем погрешность этого метода на шаге будет р(Ь) = хи„+ 0 Ь) Ьз!2, Е б (О, 1). При т = 2 из (12.19) и (12.20) с учетом условия у1+ уз =1 следует х .~~ =х +Ь((1 — уг)й1+Щ); (12.23) й1 = ДВ„, х„), йз = Д1„+азЬ, х„+ЩЬй1).
В этом случае, согласно (12.21), будем иметь р(Ь) = х(~. +Ь) -х(~.) -(1- уз) Ь~(~., ж(1.))— ~ИЬИН(1 + азЬ х(Ф ) + Я1М(~ х(Ф ))) ° 317 12.4. Метод Рунге — Кутты Обозначив Хо = Х($о), $ = Фо + агЬ) 1о = ~($о) Хо) > Х = Хо + Щ1Ьуо, трижды продифференцируем функцию ~р(Ь): Р'(Ь) =х'( +Ь) — (1 — )У вЂ” УИ *)— 12Ьа266 Х) 72ЬР21ых)1 <ро(Ь) = х" (Ф„+ Ь) — 2 уг (аг Я1, х) + Я1 ~, 'Я х) Уо)— — угЬ(аг~Уй(1, х) + 2агД1Уе~~(Ф х) Хо+Щ1У (1 х)~ )~ ро'(Ь) = хо'(Ф„+ Ь)— — 3 12 (агг ~Я х) + 2аг дгд Ух (1, х) Уо + Щ У, е(1, х) Уо) + 0(Ь) где использованы обозначения д~, д~ о д~~ е дг~ я д9 дг ' * дх' " дгг ' " д1дх' ** дхг ' Учитывая, что для ОДУ в (12.7) х =Я+Я х =Я+2Д1+1 ~ +1 е, х' = 1"., <р(0) = х(Ф„) — х($„) = 0; р'(0) = у„-(1 -72) ~„-72~„=0; Чо( ) = (1 — 2Ъаг) Й1о, х )+(1 — 27гд21)1,'(1 хо) Ло; <р'"(Ь) = (1 — Зуга~~) Д (1„, х„) + 2(1 — 372 агА1) Д'.
(1„, х„) ~о + + (1 — 3 Уг)121) У,",(Йо, ххв) Х„'+ Д(1„, х„) хо(1„). Равенство уо(0) = 0 имеет место для произвольной дифференцируемой функции 1(1, х) при условиях 1 — 2 12аг = О и 1 — 2 угД1 = О. (12.24) подставим в выражения для функции ~р(Ь) и ее производных значение 6 = 0 и получим 318 12. ПРИБЛИЖЕННЫЕ МЕТОДЫ РЕШЕНИЯ ОДУ Выполнение условий (12.24) обеспечивает при использовании (12.23) третий порядок точности на шаге и второй — на отрезке. Если принять ух = 1/2, то из (12.24) найдем ах = фг~ = 1, что соответствует рабочей формуле (12.15) нетиода Эйлера— Коши. Если же взять уз = 1, то получим ах = )3х1 = 1/2, т.е. придем к рабочей формуле (12.16) усовершенсшвованного метода ломаных. Таким образом, эти методы имеют второй порядок точности на отрезке.
Задавая различные значения параметра уг, можно построить так называемое однопараметрическое семейство методов Рунге — Кутты, имеющих второй порядок точности на отрезке. Увеличить порядок точности при гн = 2 путем выполнения равенства фв(0) = 0 в общем случае не удается. Например, для ОДУ дх/~й=х при любых значениях ух, аг и рг1 имеем ф" (0) = 1. Поскольку главная часть погрешности этого семейства методов, согласно (12.22), равна ~рв'(0) йз/6, в некоторых конкретных случаях эту величину можно уменьшить, пользуясь полученным выражением для у'в(0).
Так, если для некоторого класса ОДУ мала величина Дх", то целесообразно подобрать коэффициенты так, чтобы помимо выполнения условий (12.24) в выражении для ~раз(0) обратились в нуль первые три слагаемых, т.е. были бы выполнены равенства 1 — 3 узах — — О, 1 — 3'угааВх1 = О, 1 — 3 ухрх1 — — О. Несложно убедиться, что всем этим условиям удовлетворяют значения ух = 3/4, ах = Дх1 = 2/3, что приводит к рабочей формуле х„+1 = х„+ — (/„+3/(Ф„+26/3, х„+26/„/3)). А 4 Аналогично рассмотренному случаю тп = 2 можно построить семейства методов Рунге — Кутты для тп > 3.
При тп = 3 в общем случае в < 3, т.е. удается обеспечить третий порядок 12.4. Метод Рунге — Кутты точности на отрезке, причем одной из наиболее употребитель- ных рабочих формул является й хо+1 = хп+ — (й1 + 4йг+ йз); 6 й1 = /п~ йг = /(Зп+ й/2 хо + й й1/2) ~ йз = /(г„+ Ь, хп — й й1 + 25 йг) Из семейства методов четвертого порядка точности на отрезке, получаемых при пг = 4, наиболее часто используют вариант, приводящий к рабочей формуле Ь хо+1 — — хп+ — (й1 + 2йг+ 2йз+ й4); 6 й1=/ йг=/(З +й/2, х„+'пй1/2), йз=/(Зп+П/21 Хп+Пйг/2), й4=У(еп+й~ Хп+ййз). (12.25) Отметим, что если правая часть ОДУ в (12.7) не зависит от х, то обе последние рабочие формулы переходят в формулу Симпсона хо+1 = хп + — (.~о + 4/(еп + й/2) + /(зп + й)) й для отрезка длиной и, имеющую пятый порядок точности на шаге и четвертый — на отрезке длиной Т вЂ” 1о = Ф Й.