Численные методы. Соснин (2005) (1160462), страница 16
Текст из файла (страница 16)
Методы Адамса2 — семейство m -шаговых линейных методов решения задачи Коши,в которых беретсяa0 = 1, a1 = −1, a2 = a3 = . . . = am = 0.Таким образом, общая формула для нахождения приближения yn выглядит так:yn − yn−1= b0 fn + b1 fn−1 + . . . + bm fn−m .τ(5.10)Посмотрим, к чему приведет требование m-го порядка аппроксимации для данного метода.
Ранеебыло показано, что достаточным условием будет выполнение системы равенств: mXbi = 1;i=0mXai = 0;i=0mXiai = −1;i=1mXij−1 (iai + jbi ) = 0,j = 2, k.i=1Второе и третье равенство следуют из определения методов Адамса. Подставив значения ai впоследнее условие, получим такую систему: mXbi = 1;i=0(5.11)mXj−1 ji bi = 1, j = 2, k.i=1Эти уравнения на коэффициенты bi — достаточные условия k-го порядка аппроксимации. Их kштук, а должны они определять m + 1 неизвестное.
Очевидно, чтобы система была разрешима, необходимо выполнение неравенства:k 6 m + 1.Таким образом, максимально возможный порядок аппроксимации не может превышать m + 1.Если мы потребуем, чтобы k = m + 1, то система (5.11) даст единственное решение — то есть схема,отвечающая максимально возможному порядку аппроксимации, одна.Взглянув на формулу (5.10), легко заметить, что при b0 = 0 она становится явной относительноyn . Очевидно, что в этом случае система (5.11) содержит меньше неизвестных, и разрешима она будетуже при k 6 m.Если же b0 6= 0, то схема является неявной, и для нахождения yn приходится решать такое, вобщем случае нелинейное, уравнение:yn − τ b0 f (tn , yn ) = yn−1 + τmXi=1В общем случае оно решается методом Ньютона.2 Впервыепредложены Дж. Адамсом (1855).bm fn−i .82Глава 5. ЧИСЛЕННЫЕ МЕТОДЫ РЕШЕНИЯ КРАЕВЫХ ЗАДАЧПример 5.1.
Положим b0 = 0, m = 1. Это означает, что схема для вычисления yn будет явная,поэтому максимально возможный порядок аппроксимации будет равен m, то есть 1. Пусть будет так.От системы (5.11) останется только первое уравнение — на b1 , и из него просто получается, чтоb1 = 1. Общая схема метода будет такова: yn − yn−1 = f (tn−1 , yn−1 ), n = 1, 2, . . .τy0 = u0 — начальное условие.Пример 5.2. Здесь возьмем b0 = 0, m = 3. Максимально возможный порядок аппроксимациибудет равен 3 — его и потребуем.
При этом система (5.11) для нахождения bi будет записана так:23;b1 =12 b1 + b2 + b3 = 1;42(b1 + 2b2 + 3b3 ) = 1; ⇐⇒b2 = − ;33(b1 + 4b2 + 9b3 ) = 1. b3 = 5 .12Общая схема такова:yn − yn−123fn−1 − 16fn−2 + 5fn−3=, n > 3.τ12Значение y0 , как и раньше, берется равным u0 , а y1 и y2 обычно ищутся методами Рунге-Кутта.В данном случае они должны быть не менее, чем третьего порядка точности (так как k у нас равнотрем).Пример 5.3. Рассмотрим пример неявного метода Адамса при m = 1. Тогда можно взять k,равное 2 — максимально возможному порядку аппроксимации.
Точно так же решаем систему:b0 + b1 = 1;2b1 = 1.1 b0 = ;2⇐⇒ b = 1.12и получаем такую схему: yn − yn−1 = fn + fn−1 ,τ2y0 = u 0 .n>1В данном случае yn приходится находить из, вообще говоря, нелинейного уравнения. По-другомуего можно записать так:ττyn − f (tn , yn ) = yn−1 + fn−1 .22Обычно его решают методом Ньютона, где в качестве начального приближения берут yn−1 .Пример 5.4. Последним примером на метод Адамса будет неявный метод, с m = 2 и k = 3.Получаем для нахождения bi такую систему:5 b0 = 12 ; b0 + b1 + b2 = 1;22(b1 + 2b2 ) = 1; ⇐⇒b1 = ;33(b1 + 4b2 ) = 1. b2 = − 1 .125.4.
МЕТОДЫ АДАМСА И ГИРА83Отсюда получаем общую схему:yn − yn−15fn + 8fn−1 − fn−2=, n>2τ12y = u0 ; 0y1 ищется методом Рунге-Кутта 3-го порядка точности.yn ищется из неявной формулы по методу Ньютона.Методы ГираОпределение. Методами Гира называется семейство линейных m -шаговых методов, в которыхзаранее определяется b0 = 1, b1 = b2 = . . . = bm = 0. Общая расчетная формула будет такова:a0 yn + a1 yn−1 + .
. . + am yn−m= fn .τЗаметим, что все эти методы являются неявными, так как yn приходится находить из нелинейного,вообще говоря, уравнения:mXa0 yn − τ f (tn , yn ) = −ai yn−i .i=1В системе достаточных условий для k-го порядка аппроксимации mXbi = 1;i=0mXai = 0;i=0mXiai = −1;i=1mXij−1 (iai + jbi ) = 0,j = 2, k.i=1первое равенство будет выполняться всегда. Подставив в последнее условие значения из определенияметодов Гира, получим такую систему: mXai = 0;i=0m Xiai = −1;i=1mXij ai = 0, j = 2, k.(5.12)i=1Она содержит m + 1 неизвестное и состоит из k + 1 уравнения.
Таким образом, чтобы обеспечитьразрешимость системы и требование k-го порядка аппроксимации, приходится ограничивать k числомm. При k = m порядок аппроксимации будет максимально возможным, а схема, им определяемая,тоже будет одна.84Глава 5. ЧИСЛЕННЫЕ МЕТОДЫ РЕШЕНИЯ КРАЕВЫХ ЗАДАЧПример 5.5. В первом простейшем примере берем m = 1, k = 1. Система для нахождения aiбудет такова:a0 + a1 = 0;a0 = 1;⇐⇒a1 = −1.a1 = −1.Отсюда получаем простейшую разностную схему:( y −ynn−1= f (tn , yn ),τy0 = u 0 .n > 1;Пример 5.6.
Взяв в этом примере m = 2 и k = 2, мы получим такую систему уравнений: a0 + a1 + a2 = 0;a1 + 2a2 = −1;a1 + 4a2 = 0.3a0 = ;2⇐⇒a=−2;1 a = 1.22Соответствующая схема будет выглядеть так:3yn − 4yn−1 + yn−2= τ f (tn , yn ),2n>2y = u0 ; 0y1 ищется методом Рунге-Кутта 2-го порядка точности.То есть, для определения yn мы имеем такое нелинейное уравнение:13yn − τ f (tn yn ) = 2yn−1 − yn−2 .22Пример 5.7. В этом примере мы возьмем m = 3, k = 3. Тогда последнее условие в системе (5.12)разобьется на два уравнения, и мы получим такую систему:11a0 =;6a+a+a+a=0;0123 a1 = −3;a1 + 2a2 + 3a3 = −1;3⇐⇒a1 + 4a2 + 9a3 = 0;a2 = ;2a1 + 8a2 + 27a3 = 0. a3 = − 1 .3Этому будет соответствовать такие расчетные формулы:11yn − 18yn−1 + 9yn−2 − 2yn−3= τ f (tn , yn ),6n>2y = u0 ; 0y1 , y2 ищутся методом Рунге-Кутта 3-го порядка точности.Замечание. На практике используются методы Гира вплоть до десятого порядка аппроксимации(и соответствующей точности).
Это связано с тем, что эти методы обладают свойством вычислительнойустойчивости, о котором мы поговорим в следующем разделе.5.5. УСТОЙЧИВОСТЬ ЧИСЛЕННЫХ МЕТОДОВ РЕШЕНИЯ ЗАДАЧИ КОШИ5.585Устойчивость численных методов решения задачи КошиМы будем рассматривать численные методы для поиска функции, являющейся решением такой задачиКоши:du= f (t, u), t > 0;dt u t=0 = u0 .Для обоснования дальнейших действий сначала проведем теоретические рассуждения. Будем считать, что функция u является решением задачи Коши:dudtut=0= f (t, u),t > 0;(5.13)= u0 .Пусть функция U — решение аналогичной задачи, но с «возмущенными» начальными данными:dUdt Ut=0= f (t, U ),t > 0;(5.14)= u0 + u0 .Представим ее следующим образом: U = u + u, где u — функция-погрешность. Исследуем u,подставив представление для U в (5.14) и разложив f (t, u) в ряд Тейлора по второй переменной:(udu du+dt dt+ ut=0t=0= f (t, u + u) = f (t, u) + fu0 (t, u)u + Rf ,t > 0;= u0 + u0 .Используя (5.13), получим:(du dtut=0= fu0 (t, u)u + Rf ,t > 0;(5.15)= u0 .Из курса «Дифференциальные уравнения» известно, что если fu0 (t, u) < 0, то функция u(t) монотонно стремится к нулю на бесконечности.
В этом случае говорят, что задача (5.13) устойчива поначальным данным.Теперь вернемся к численным методам. Если их применять к возмущенной задаче (при условииfu0 < 0 ), то естественно требовать, чтобы функция u получалась убывающей (или, что тоже самое,функция U приближалась к u, то есть погрешность не накапливалась). После проведенных выкладокпонятно, что это то же самое, что требовать получение убывающей функции при решении системы(5.15).
Поэтому методы обычно тестируют на такой модельной задаче:(du dtut=0= λu,t > 0;λ = const < 0;= u0 ,справедливо ожидая, что если при ее решении будет получаться убывающая на бесконечности функция, то метод, примененный к системе (5.14), даст функцию, не сильно отличающуюся от u. Такиеметоды называются устойчивыми.Примечание.
Мы требуем убывания от функции, если u0 , задающее начальное условие, положительно. Если же оно отрицательно, то логично требовать возрастание, или, точнее, стремлениек нулю слева.86Глава 5. ЧИСЛЕННЫЕ МЕТОДЫ РЕШЕНИЯ КРАЕВЫХ ЗАДАЧБудем называть метод условно устойчивым, если он устойчив при некоторых значениях своихпараметров.Теперь рассмотрим несколько связанных с этой модельной задачей примеров.Пример 5.8. Явный одношаговый метод Адамса первого порядка аппроксимации. Разностная схема имеет вид:yn+1 − yn= f (tn , yn ).τПравая часть, согласно дифференциальному уравнению, выражается как f (tn , yn ) = λyn , тогдасхема будет иметь следующий вид:yn+1 = (1 + τ λ)yn .Как уже говорилось, метод называется устойчивым, если сеточная функция yn не возрастает поn.
В нашем примере |yn+1 | = |1 + τ λ| · |yn | не будет возрастать, если |1 + τ λ| 6 1, то есть при−2 6 τ λ 6 0.Правая часть неравенства выполняется всегда, так как τ > 0, а λ < 0. Из этого следует, что явныйодношаговый метод Адамса удовлетворяет условию устойчивости лишь приτ 6−22=.λ|λ|Этот метод условно устойчивый, и это не очень хорошо, так как, чтобы не набирать погрешность,надо учитывать ограничение на τ, где величина λ зависит от поведения функции f (t, u) (λ — разностный аналог производной функции f (t, u)).