1626435587-51311eae4652e8ad616b5bdef025cbb3 (844239), страница 19
Текст из файла (страница 19)
Каждаяпара шагов численного интегрирования дублируется при этом шагомдвойной величины; из сравнения полученных результатов делается вывод относительно погрешности вычислений и необходимости корректировки шага.Оценим замедление расчётов, обусловленное наличием дополнительных вычислений правой части уравнения (70). На каждом шагечисленного интегрирования по схеме Рунге — Кутты четвёртого порядка (85) правую часть (, ) необходимо вычислять четыре раза (одинраз в начале шага, два — в середине, и ещё один — на правом концеотрезка). Следовательно, на паре шагов правая часть (, ) должнабыть вычислена восемь раз. При пересчёте с двойным шагом необходимо вычислить (, ) ещё четыре раза.
Итого получаем двенадцать вычислений (, ). Заметим, однако, что начальная точка для двойногошага совпадает с начальной точкой для пары двойных шагов (рис. 8),что позволяет сократить число вычислений (, ) до одиннадцати накаждой паре шагов. Учитывая, что без пересчёта на сетке с шагом 2ℎмы должны были сделать 8 шагов, получаем коэффициент замедления11/8 = 1,375. Однако, «теряя» 40 % операций на дублирование расчётов на сетке с двойным шагом, при решении ряда задач можно добитьсясокращения общего времени вычислений на порядок величины и болеепри фиксированной точности за счёт увеличения шага сетки там, гдефункция (, ) мала, и решение меняется медленно.Рис.
8. Вычисление (, )при пересчёте с двойным шагомПосле выполнения двух шагов ℎ численного интегрирования по схеме (85) и одного шага 2ℎ мы получим пару значений сеточной функции , назовём их 1 и 2 (NB: речь идёт о значениях численного решенияв одном и том же узле + 2ℎ, полученных в результате вычислений насетках с одинарным и двойным шагом):( + 2ℎ) = 1 + 2 · ℎ5 + (ℎ6 ),56( + 2ℎ) = 2 + · (2ℎ) + (ℎ ),97(86)где — коэффициент пропорциональности в остаточном члене =(ℎ5 ) = ℎ5 . В первом случае (расчёт на сетке с шагом ℎ) суммарнаяпогрешность складывается из погрешностей на двух шагах ℎ, поэтому перед коэффициентом пропорциональности в остаточном членестоит множитель 2.Вычитая одно уравнение из другого, несложно получить оценку дляпогрешности вычислений:1 ≡ 2ℎ5 + (ℎ6 ) =2 − 1+ (ℎ6 ).15(87)Полученная разность является мерой погрешности численного интегрирования на одном шаге и может быть использована для принятиярешения о корректировке шага численного интегрирования.
Для этогонеобходимо сравнить фактическую погрешность (87) с заданным значением так, чтобы выполнялось неравенство1 < + || · ,где и — максимально допустимая на одном шаге абсолютная и относительная погрешность соответственно. В случае невыполнения указанного неравенства нужно уменьшить шаг ℎ и выполнить пересчётпоследнего шага.
Если 1 /( + ||) < 1, можно оценить безопасноеувеличение величины шага сетки на следующем шаге интегрирования,используя известную степенную зависимость (86) остаточного члена отшага сетки ℎ.Заметим также, что учитывая найденную погрешность в выражении (86), можно повысить порядок точности численной схемы (фактически, используя):метод Рунге — Ромберга( + 2ℎ) = 1 +2 − 1+ (ℎ6 ).15(88)Проводя вычисления с использованием (88), можно повысить порядокточности метода, однако при этом нельзя будет оценить погрешностьвычислений, не прибегая к дублированию расчётов ещё на одной сетке.6.6.
Многошаговые методыВ случае, если правая часть уравнения (70) имеет большое числонепрерывных производных, может быть целесообразно использоватьчисленные схемы более высокого порядка точности. Время, затраченное на написание и отладку программ, реализующих более сложные98методы высокого порядка точности, зачастую с лихвой компенсируетсяповышением скорости расчётов за счёт использования сетки с меньшимчислом узлов при сохранении заданного уровня точности вычислений.Ниже будет кратко рассмотрен класс относительно простых по формезаписи и в реализации многошаговых методов высокого порядка точности.Идея, лежащая в основе многошаговых методов, заключается в использовании в расчётах вместо (, ) в правой части уравнения (70)некоторой другой функции (), аппроксимирующей с заданным порядком точности: () = (, ()) + (ℎ ).В наиболее простом и часто встречающемся случае в качестве ()используют интерполяционный полином.
Очевидно, что для построения интерполяционного полинома степени ≥ 1 и, следовательно, длявыполнения шага численного интегрирования ОДУ уже недостаточнознания одного лишь значения ( ) — необходимо иметь значения ()в + 1 узле сетки (на + 1 предыдущем шаге численного интегрирования). Именно поэтому рассматриваемые методы называются.
Считая значения , −1 , . . . , − известными, запишеминтерполяционный полином в форме Ньютона (см. (54) на с. 75):мно-гошаговыми ()= ( ) +(89)+ ( − ) ( , −1 ) ++ ( − )( − −1 ) ( , −1 , −2 ) ++ ( − )( − −1 )( − −2 ) ( , −1 , −2 , −3 ) + . . .Для вычисления решения () в следующем узле сетки запишем дифференциальное уравнение (70) в интегральной форме∫︁+1∫︁+1 (, ()) ≈ ++1 = + () и подставим в него интерполяционный многочлен (89) третьей степени.Вычисляя возникающие при этом элементарные интегралы вида∫︁+1( − ) =ℎ2,2где ℎ = +1 − ,∫︁+1( − ) · ( − −1 ) =99ℎ3ℎ2+ ℎ−1 ,32∫︁+1(︃∏︁( − ) ==−2=ℎ4ℎ3ℎ2+ (2ℎ−1 + ℎ−2 ) + ℎ−1 (ℎ−1 + ℎ−2 ),432приходим в итоге к+1)︃формуле Адамса для переменного шага :1= + ℎ ( ) + ℎ2 ( , −1 ) +[︂]︂21 3 1 2+ ℎ + ℎ ℎ−1 · ( , −1 , −2 ) +(90)32[︂ 4]︂ℎ11+ + ℎ3 (2ℎ−1 + ℎ−2 ) + ℎ2 ℎ−1 (ℎ−1 + ℎ−2 ) ×432× ( , −1 , −2 , −3 ).Схема (90) имеет четвёртый порядок точности.
Если отбросить последнее слагаемое, получим формулу третьего порядка точности и т. д.Оставляя в (90) только одно слагаемое +1 = + ℎ ( ), получимсхему ломаных (73), имеющую первый порядок точности.Рассмотренный выше метод Адамса четвёртого порядка является, т.к. значение +1 может быть вычислено по предыдущим значениям , −1 , . . . с использованием выражения (90) за фиксированное число операций. Явные методы Адамса также называют в литературе методами(Adams, Bashforth).методы Адамса — схемами(Adams, Moulton).На первый взгляд, метод Адамса (90) кажется привлекательнымтем, что на каждом шаге необходимо вычислять (, ) лишь одинраз, тогда как в методе Рунге-Кутты (85), также имеющем четвёртый порядок точности, (, ) необходимо вычислять четыре раза накаждом шаге.
Это могло бы дать существенный выигрыш в скорости счёта, особенно для тех случаев, когда вычисление (, ) занимает много времени. Однако если оценить коэффициент в остаточномчлене схемы (90) в простейшем случае равномерной сетки, мы получим251 5 IVAdams = 750ℎ () [2, гл. VIII §7], тогда как остаточный член дляметода Рунге — Кутты (85) на равномерной сетке совпадает с остаточ1ным членом для метода Симпсона и равен Simpson = 2880ℎ5 IV . Такимобразом, остаточный член в четырёхчленной схеме Рунге — Кутты (85)в 960 раз меньше, чем в методе√ Адамса, что позволяет использоватьсетку с шагом, увеличенным в 4 960 ≈ 5,57 раз, т. е.
фактически вычислять (, ) меньшее число раз, чем в методе Адамса того же порядка100явнымАдамса — БашфортаАдамса — МултонаНеявныеточности. Относительно большое значение численного коэффициента востаточном члене метода Адамса обусловлено использованием полиномов для экстраполяции функции () ≡ (, ()) за пределы отрезка[− , ], на котором был построен интерполяционный полином (см.п.
5.3).Повышая степень интерполяционного полинома (89), несложно построить численную схему Адамса более высокого порядка точности,хотя даже схема четвёртого порядка (90) на неравномерной сетке достаточно громоздка. К числу «минусов» методов Адамса следует также отнести необходимость задания начальных условий в узлах сетки(для использования метода -го порядка). Это приводит к необходимости выполнения начальных расчётов с использованием других (самостартующих) схем (например, метода Рунге — Кутты четвёртого порядка (85)), что может вдвое увеличить размер программы для ЭВМ.
Данное обстоятельство, наряду с относительно большим коэффициентом востаточном члене, ограничивает практическую применимость методовАдамса.6.7. Жёсткие системы уравненийРассмотрим понятие устойчивости методов численного интегрирования ОДУ на примере уравнения′ = −,>0(91)Решение уравнения (91): () = (0)− . Для любого > 0 справедливо |(+)| ≤ (), т. е. решение уравнения асимптотически устойчиво.Естественно ожидать устойчивости и от численного решения. Для ,построенного по схеме Эйлера (73), имеем:+1 = + ℎ ( , ) = (1 − ℎ) .Для того, чтобы численное решение было устойчивым, необходимо выполнение условия |1 − ℎ| < 1, откудаℎ<2.(92)Следовательно, численное решение уравнения (91), полученное по схеме Эйлера (73), является устойчивым при выполнении условия (92).
Вэтой связи говорят, что численная схема (73) является.условно устой-чивой101Рассмотрим численное решение того же уравнения (91), полученного с помощьюсхемы Эйлера:неявной(93)+1 = + ℎ (+1 , +1 ).яв-Обратим внимание, что если явный метод Эйлера наряду с другимиметодами даёт выражение для вычисления +1 , тосхема (93) представляет собой уравнение, в которое искомое значение +1численного решения в следующем узле сетки входит в аргумент функции (, ). Схема (93), следовательно, является нелинейным конечным уравнением на +1 . Исследуем асимптотическую устойчивостьчисленного решения, построенного по схеме (93) для уравнения (91):ныминеявная+1 = + ℎ (+1 , +1 ) = − ℎ +1 ,откуда.1 + ℎВидно, что для любого шага сетки ℎ численное решение, построенноепо неявной схеме (93), будет устойчивым. В этой связи говорят, чтонеявная численная схема (93) является, в отличие отявной схемы (73). Общий вывод можетбыть сформулирован следующим образом:+1 =абсолютно устойчивойусловно устойчивойявные численные схемы длярешения ОДУ могут обладать лишь условной устойчивостью, тогдакак среди неявных схем есть абсолютно устойчивые.Обратим внимание, что для получения приемлемой точности численного решения уравнения (91) с использованием схем первого порядка (73) и (93) необходимо выполнение условия ℎ ≪ 1/, котороеявляется более сильным, чем полученное выше условие устойчивости(92).
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
