Амосов А.А., Дубинский В.А., Копченова Н.В. Вычислительные методы для инженеров (1994) (1095853), страница 79
Текст из файла (страница 79)
Епн 0 Рис. 14.10 Рис. 14.0 ~ 14.5. Модификации метода Эйлера второго порядка точности Медленная сходимость метода Эйлера (его погрешность убывает пропорционально лишь первой степени Ь) является серьезным препятствием для использования его на практике. Из рис. 14.7 видно, что уже один шаг по касательной к интегральной кривой приводит к значительной величине локальной погрешности 1„. Можно ли так подправить расчетную формулу метода, чтобы существенно уменьшить величину 1„7 Пусть у ф — решение дифференциального уравнения (14.55) '(Ю = У(1 у 0), удовлетворяющее условию у (1„) = у„. Далее, пусть Й„= = СД7а (14.56) Ь вЂ” угловой коэффициент секущей, проходящей через точки (1„, у (1„)) и (1,,~, у (1„~)) графика функции у (1) (рис.
14.10). 435 Значение Ьо„т, как правило, бывает очень трудно определить заранее. Однако если очень высокая точность не нужна, то необходимый для ее достижения шаг Ь обычно бывает много больше Ь„„,. Если же требуется высокая точность решения, то достичь ее с помощью метода Эйлера нельзя, даже если пойти на значительные затраты машинного времени (неизбежные при расчете с малым значением шага Ь). Необходимо иметь в своем распоряжении методы, имеющие более высокий порядок точности и позволяющие вести расчет со. сравнительно крупным шагом Ь.
Ясно, что "метод", состоящий в вычислении по формуле Уп+] — Уп + ~~п~ (14.57) дим к равенству у( ) - у( ) = l У(1, у(~)М (14.58) Из равенств (14.56) и (14.58) следует, что (14.59); ~п+1 ~ =-„ / П~, у(1))4~ Заметим теперь, что применение для приближенного вычисления интеграла, стоящего в правой части выражения (14.59), формулы левых прямоугольников / ~(1, у (1))Й1 в Ц'(1„, у (1„-)) немедленно приводит от (14.5Т) к методу Эйлера (14.46). Известно (см. гл. 13), что больший порядок точности имеет формула трапеций ~+~ й / й~ у(~)) 6~ У(~ у(1))+ Х(~ у( ))) Непосредственное ее применение к вычислению Фп приводит к правилу трапеций: Ь Уп+1 — Уп+ 2 У(~п) Уп) + 1(~п+ь Упм)) (14.60) (ср.
с (14.25)). Этот метод имеет второй порядок точности, но является неявным. Поэтому его реализация связана с необходимостью решения относительно уп+1 нелинейного уравнения (14.60). 436 имеет нулевую локальную погрешность. Для того чтобы воспользоваться этой формулой, нужно лишь "научиться вычислять значение Фпп. Интегрируя обе части уравнения (14.55) по 1 от 1п до Ц„~ и исполь-' Мп+1 зуя формулу Ньютона — Лейбница ( у'(М)АМ = у (М„,1) — у (М„), прихо-' Ь уп ) = у + — У(~п, уп) + У(1.), уп+ 1У(1п, уп))), (14.61) который называют А(етодо а Эйлера — Кожи (или л(етодол( Хьюпа).
Геометрическая иллюстрация этого метода представлена на рис. 14.11. Вычисления разбивают на два этапа. На первом этапе (этапе прогноза) в соответствии с методом Эйлера уп+) — — у„+ Ь1'„, 1п (о) (1) (1) у„) вычисляют грубое приближение к значению у (1„,)), В точке (1„,), у„,)) определяют угловой коэффициент Йп = ~(1„,), у„,)). (2) (О) На втором этапе (этапе коррекции) вычисляют усредненное значение углового коэффициента 1„= (14~' + й„'~~)/2. Уточненное значение у„+( находят по формуле уп+) — у„+ И„, что соответствует шагу по прямой, проходящей через точку (1„, - уп) и имеющей у) ловой коэффициент, равный Й„.
Рис. 14.11 3 а м е ч а н и е. Метод Эйлера — Коши относится к классу мето- дов прогноза и коррекции (иначе говоря, л(етодоо типа предик- тор — корректор'). ) От англ. Во рге(1)с1 — "предсказывать", "прогнозировать"; йо соггесй— "исправлять", "корректировать". 437 Построим на основе правила трапеций явный метод. Для этого подставим в правую часть формулы (14.60) значение у„,), "предсказываемое" методом Эйлера. В результате получается метод Метод (14.61), который можно рассматривать как модификацию метода Эйлера, имеет второй порядок точности. Еще одну модифика-, цию второго порядка точности можно получить с помощью формулы;,.
(центральных) прямоугольников оп+1 l У(1~ У (1))(11~ Ь~(1 1а У (1мд)), 1п+)а = 1п+ Ь/2, тп Уп+1/2 Уп + Ьп ~ Ьп — 1(1п Уп) (1) (1) 2 (14.62) Уп+1 — Уп + ЬЬп Ьп †,1(1п+1/2 Уп+1 /2) Геометрическая иллюстрация этого метода приведена на рис. 14.12. Пример 14.14. Применим рас— смотренные в этом параграфе методы для численного решений' задачи (14.45) с шагом Ь = 0.1.
Расчетная формула метода Эйлера †Ко принимает вид Ь у 1 = уп+ — (21„уп + 2 + 21п+1(уп+ Ь21пуп)) = (1463) — уп[1 + Ь (тп + (и+1 + 1а (к+Й=~п+$ ~и+/ Рис. Ц,12 + 2Ьтп(п+1)1. Вычисления усовершенствованного метода Эйлера производим по формуле Уп+1 = Уп + 2Ь(тп+ — )(Уп + — 21пуп) = Уп(1+ Ь(22п+ Ь)(1+ ЬЦ))),(14.64) Ь. Ь Правило трапеций (14.60) приводит к уравнению Ь Уп+1 = Уп + — (21пУп + 21п+1уп.1), 2 438 если для приближенного вычисления значения У (1п,) 2) применить метод Эйлера. В результате получим расчетные формулы Усовер(иенст- ° вованно(о метода Эйлера которое в данном случае линейно и легко разрешается относительно у„~..
у„., = у.(1+ Ы„)(1 — Ьт„,,) (14.65) Результаты вычислений по формулам (14.63), (14.64) и (14.65) приведены в 2 табл. 14.3. Там же для сравнения представлены значения решения у (1) = е Нижняя строка таблицы содержит величину абсолютной погрешности Е (Ь). Т а б л и ц а 14.3 Метод Эйлера- Усовершенство- Правило Коши; Ь = 0.1 ванный метод трапеций; Эйлера; Ь = 0.1 Ь = 0.1 Точное решение Е(Ь) 0.93 10 2 2.0.10 т 1,9 10- Как видно из сравнения табл. 14.3 с табл.
14.2, проведенные в этом параграфе модифйкации метода Эйлера действительно привели к повышению точности, 3 14.6. Методы Рунге — Кутты Наиболее популярными среди классических явных одношаговых методов являются методы Рунге — Кутты'. Методы Эйлера, Эйлера— Коши и усовершенствованный метод Эйлера можно рассматривать как простейших представителей этого класса методов. 1. Вывод расчетных формул. Поясним вывод расчетных формул метода Рунге — Кутты. Пусть (как и в 8 14.5) у (1) — решение диффе- ' В.М.Кутта (1867 — 1 944) — немецкий математик. 439 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 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.01010 1. 04102 1.09468 1.17450 1. 28577 1.43624 1.63700 1.90390 Ю. 25958 2.73660 1.01005 1.04081 1. 09417 1.17351 1. 28403 1. 43333 1. 63232 1. 89648 2.24791 2.71828 ренциального уравнения у = »" (», у), удовлетворяющеЕ условию у (»„) = у„. Запишем равенство (14.58) в следующем виде: М у(» ° ) = у(» ) + Г У(», у(»))6». (14.66) Если бы входящий в зто равенство интеграл можно было вычислить точно, то получилась бы простая формула, позволяющая последовательно вычислить значения решения в узлах сетки.
Поскольку в действительности эта невозможно, попробуем получить приближенную формулу, заменив интеграл квадратурной суммой (см. гл. 13). Введем на отрезке [»„, »„,1]и вспомогательных узлов»„= »„+ а(Ь> (1) 4,~~ = »„+ а2Л, ..., »~~,"~ = »„+ о((>Л, где О = а( ~ а2 ( ... ( а>)> ( 1.3аме- у (»~.1) у (»„) + Ь Е ~(~(»„'", у (»„'" )). (Ф 1 (14.67) Однако воспользоваться равенством (14.67) для вычисления у (»„,1) нельзя, так как значения функции у в точках»()' для» = 2, 3, ..., ш (() неизвестны. Чтобы найти эти значения, запишем равенства (() и у (»„'") = у (»„) + ~ У(», у (»))(1» (» = 2, 3, ..., т), »и (14.68) аналогичные равенству (14.66). Заменяя для каждого» входящий в формулу (14.68) интеграл соответствующей ему квадратурной форму- лой с узлами»„', »(,2), ..., 4', придем к приближенным равенствам у (»>) ) у (»а) + М21У(»й > у (»>> ))> у(»'") 'у(».) + Ь(Р Х(»'", у(»'")) + Р22Х(»'", у(»"'))), у (»'") ' у (» ) + Ь ~ Р4М' у Ф')) 1=1 440 тим, что»(,1) = »„, 4~~ 4»„,1.
Заменяя входящий в равенство (14.66) интеграл квадратурной суммой с узлами»„' ', ..., »„', получаем приб- лиженное равенство позволяющим последовательно вычислять приближения к значениям „(~йг) ) „(2йп)) ) Обозначим теперь через уп вспомогательные величины, имеющие ( () смысл приближений к значениям у (~йо); пусть Йй" = ~(1й", у('))— приближение к значению углового коэффициента й в точке 2~". Тогда расчетные формулы примут вид Уп+1 = Уп+ ~(~п> ~~п = 1((п > Уп )> Уп = Уп (1) (1] (1) (1) ~Сй = я~2) у(2)) уйг) = уп+ /))921М> (г) )сй = У(1й > Уп )> Уй = Уп+ "()эаФй + (гаг)сй )> »>-1 )сй = У(~й > Уй )> Уп = Уп + () ~ 1»оМй )=1 Часто из этих формул исключают вспомогательные величины уп(') и записывают формулы так: Уп+1 = Уп + )))еп, >сп = г сап (() 1=1 ~й =,)' (1п> Уп)> = ~(2п+ ог"> Уп+ )))эг1"й )> Ма) = 1(~„+ аА уп+ ЬРа))~йм + Рагйй")), Заметим, что выведенные формулы задают явный одношаговый метод вида у„+1 — уп + ЬФ (1п, уп, Ь), где для вычисления значений »> функции Ф (1п, уп, Ь) = г с(ййо используются значения правой части (=1 ,) в т вспомогательных точках.
Поэтому этот метод называют яоныл( т — этапныл( .я(етодол( Рук(е — Кутты. Выбор конкретных значений параметров сь аь )э( осуществляется исходя из различных соображений. Естественно, что одним из основных является желание сделать порядок аппроксимации р максимально возможным. 441 2. Устойчивость и сходимость. Следующая теорема позволяет в дальнейшем называть методы Рунге — Кутты, имеющие р-й порядок аппроксимации, методами р-го порядка точности. Т е о р е м а 14.7.
Пусть правая часть дифференциального уравнения удовлетворяет условиго ~ 1" ~ ч 1. Тогда всякий явный т-этапный У иетод Рунге — Кутпгы устойчив на конечно.и отреэке. С л е д с т в и е. Пусть выполнено условие ~~' ~ ч Ь. Тогда если У явный т-этапный лгетод Рунге — Крпты илгеет р-й порядок аппроксилгации, то он сходится с р-.и порядколг точности. Справедливость следствия вытекает из теоремы 14.4. 3. Семейство явных двухзтапных методов Выведем расчетные формулы семейства явных двухзтапных методов Рунге — Кутты второго порядка точности. Запишем формулы явного двухзтапного метода Уп+1 — Уп + Ь (с!Ьп + с2Ьп )) (11 (21 Ьп" = 1(1п, Уп), Ьп" = 1(~п+ оЬ, Уп+ МЬп") в виде = с1~(~п, Уп) + с 0(~п + оЬ, Уп + (УЬУ(~п, Уп)). Параметрами этого метода являются величины сг, с2, о, (!1.