Апостериорные оценки погрешностей по Рунге и Эйткену

 

Априорные оценки погрешностей (4.7) и (4.10) можно записать в виде

, (4.11)

где А — коэффициент, зависящий от метода интегрирования и вида подынтегральной функции; h — шаг интегрирования; р — порядок ме­тода. Зависимости (4.11) подчиняется главный член погрешности боль­шинства методов численного интегрирования. При численном диффе­ренцировании погрешность также может быть оценена с помощью фор­мулы (4.11), при этом порядок р зависит от количества узловых точек. Пусть вычисляется значение некоторой переменной w с шагом /г, тогда

, (4.12)

где — приближенное значение w; — главный член погрешно­сти; О (h^p+1) — бесконечно малая величина порядка h^p+1. Вычислим ту же самую переменную w с шагом kh

, (4.13)

где коэффициент пропорциональности к может быть как больше, так и меньше единицы. Коэффициент А в выражениях (4.12) и (4.13) будет одинаковым, так как вычисляется одна и та же переменная, одним и тем же методом, а от величины шага h значение А не зависит.

Пренебрегая бесконечно малыми величинами, приравняем правые части соотношений (4.12) и (4.13) с учетом формулы (4.11) и получим

откуда найдем главный член погрешности

(4.14)

Формула (4.14) называется формулой Рунге [2, 1] и позволяет пу­тем двойного просчета величины w с шагами h и kh оценить погреш­ность. Так как оценка осуществляется после вычисления, то она явля­ется апостериорной. Формула (4.14) имеет большое практическое зна­чение, так как позволяет провести оценку погрешности без изменения алгоритма используемого вычислительного процесса. При уменьше­нии шага h главный член погрешности Rо будет стремиться к полной погрешности R.

После определения Rо можно вычислить уточненное значение ис­комой величины

(4.15)

последнее выражение называют второй формулой Рунге. К сожале­нию, погрешность уточненного значения остается неопределенной, хо­тя, как правило, она меньше значения .

Формулы Рунге справедливы для всех вычислительных процессов, для которых выполняется степенной закон (4.11). Для определения по­рядка метода р необходимо проведение априорной оценки погрешно­сти, что не всегда легко осуществить.

Английский математик Эйткен предложил способ оценки погреш­ности для случая, когда порядок р метода неизвестен. Более того, алго­ритм Эйткена позволяет опытным путем определить и порядок метода. Для этого необходимо третий раз вычислить значение величины w с шагом k²h, т.е.

или

(4.16)

Приравнивая правые части выражений (4.13) и (4.16), получим со­отношение

 

 

Подставив в это соотношение значение Rо из первой формулы Рунге (4.14), найдем

(4.17)

Полученное соотношение (4.17) вместе с первой формулой Рунге (4.14) позволяет оценить погрешность при использовании вычис­лительного метода с неизвестным порядком р. Более того, порядок р можно определить, логарифмируя левую и правую части формулы (4.17):

(4.18)

Для выбранного вычислительного процесса алгоритм Эйткена до­статочно применить только один раз при определении порядка метода, а затем использовать формулу Рунге, требующую только двукратного вычисления искомой величины. Формулу (4.14) можно использовать для тестирования программ, реализующих вычислительные методы с известной априорной погрешностью. Априорный и апостериорный по­рядки должны получаться совпадающими для правильных программ. Конечно, это совпадение будет приближенным, так как при получении алгоритмов Рунге и Эйткена учитывались только главные члены по­грешности.

 

Метод трапеций

Подинтегральную функцию заменим на участке полино­мом первой степени . Как и в методах прямоугольников, такая аппроксимация неоднозначна. Одним из возможных способов являет­ся проведение прямой через значения функции на границах интервала интегрирования (рис. 4.4). В этом случае приближенное значение ин-ч теграла определяется как площадь трапеции:

(4.19)

Априорную погрешность R метода трапеций получим путем инте­грирования тейлоровского разложения подынтегральной функции око­ло ТОЧКИ Xi

(4.20)

(4.21)

 

f(xi+h) f(x) Рис. 4.4. Метод трапеций

 

 

 

С помощью разложения (4.20) вычислим подынтегральную функ­цию в точке :

,

откуда

(4.22)

Подставляя произведение (4.22) в выражение (4.21), получим

(4.23)

Следовательно, главный член погрешности метода трапеций на одном интервале будет

(4.24)

Если интегрирование проводится путем разбиения отрезка [ x₀;x_п ] на несколько интервалов, то общую погрешность получим суммирова­нием частичных погрешностей (4.24):

(4.25)

Получили, на первый взгляд, несколько неожиданный результат, оказалось, что метод трапеций имеет погрешность в два раза боль­ше по абсолютной величине, чем метод средних прямоугольников, хо­тя аппроксимация подынтегральной функции проводилась полиномом первой, а не нулевой степени. По-видимому, выбранный вариант ап­проксимации подынтегральной функции прямой, проходящей через ее значения на границах, не является оптимальным. Задача выбора спо­соба аппроксимации полиномом заданной степени с наименьшей воз­можной погрешностью была решена Гауссом, что привело к развитию целого класса методов.

Как видно из выражения (4.25), метод трапеций, как и метод сред­них прямоугольников, имеет второй порядок. Если подынтегральная функция задана аналитически, то предпочтительнее из методов второ­го порядка применять метод средних прямоугольников вследствие его меньшей погрешности.

 

Метод Симпсона

Подынтегральную функцию f(x) заменим интерполяционным по­линомом второй степени Р₂(х) — параболой, проходящей через узлы (рис. 4.5), тогда

где R — погрешность вычисления интеграла.

x₀ x_i x₂ x

Рис. 4.5. Метод Симпсона

Для записи полинома Р₂{х) воспользуемся интерполяционной фор­мулой Ньютона (2.6) для трех узлов:

, (4.26)

где и — разделенные разности, определяемые по формулам -

 

h - расстояние между узлами.

Введем новую переменную , тогда и полином (4.26) принимает вид

. (4.27)

Теперь вычислим интеграл от полинома (4.27):

(4.28)

Последнее соотношение называют квадратурной формулой Симпсона или формулой парабол.

Формулу Симпсона можно получить и с помощью первой и второй формул Рунге, примененных к вычислению интеграла методом трапе­ций. Запишем два приближенных значения интеграла от функции f(x) на интервале [ , ] с шагами h и 2h по формуле трапеций (4.19).:

(4.29)

Интегралы (4.29) подставим в формулы (4.14) и (4.15) и получим уточненное значение интеграла

которое совпадает с формулой Симпсона (4.28).

Для оценки погрешности формулы Симпсона разложим подынте­гральную функцию f(x) в ряд Тейлора около точки х\ и проинтегриру­ем разложение почленно на интервале [ , ]

(4.30)

 

 

Суммируя разложения около точки х_г для функции f(x) в узлах x₀ и x₂, получим, что

тогда интеграл (4.30) принимает вид

(4.31)

 

Первое слагаемое в правой части формулы (4.31) совпадает с форму­лой Симпсона, значит, второе слагаемое является главным членом по­грешности для интеграла на интервале [ , х₂]

(4.32)

Если интеграл вычисляется на интервале [ , ] путем разбиения его на четное число подьштегралов [ ], на каждой паре которых применяется формула Симпсона для узлов , то полная по­грешность будет суммой правых частей соотношения (4.32). При ма­лой величине шага h на основании метода средних прямоугольников получим

,

тогда полная погрешность запишется в виде

(4.33)

следовательно, формула Симпсона имеет четвертый порядок точности с очень малым численным коэффициентом в остаточном члене. Фор­мула Симпсона позволяет получить высокую точность, если четвертая производная подынтегральной функции не слишком велика. В против­ном случае методы второго порядка могут дать большую точность, чем метод` Симпсона. Например,`[17],для`функции формула трапеций при n = 2 для интеграла в пределах [—1,1] дает точный результат, равный 4, тогда как по формуле Симпсона получим результат, не совпадающий даже по знаку (—8/3).

Методы Ньютона-Котеса

В формулах Ньютона-Котеса порядок полинома п, которым за­меняется подынтегральная функция f(x) может принимать различные значения. В этом случае на интервале [а, Ь] интеграл равен

 

 

где — коэффициенты Ньютона-Котеса; — общий знаменатель дробей, выражающих коэффициенты в формуле численного интегри­рования; — относительные координаты (на промежутке [0,1])узло­вых точек, по которым строится интерполирующий полином.

В табл. 4.1 приведены значения и , входящих в формулу (4.34), для п = 1... 5 и выражения для главного члена погрешности . При п = 1 формулы Ньютона-Котеса дают формулу трапеций, п = 2 — формулу Симпсона.

 

Таблица 4.1. Коэффициенты формулы Ньютона-Котеса