\documentclass[12pt,a4paper,oneside,reqno,draft]{amsbook}
\usepackage{amsmath}
\usepackage{amsthm}
\usepackage[T2A]{fontenc}
\usepackage[cp1251]{inputenc}
\usepackage[english,russian]{babel}
\usepackage{amsfonts}
\usepackage{latexsym}

\tolerance 900
\newtheorem{lemma}{Лемма}
\newtheorem{theorem}{Теорема}[section]
\newtheorem{example}{Пример}[section]

\newtheorem{proposition}{Предложение}
\newtheorem{corollary}{Следствие}

\makeatletter
\@addtoreset{equation}{section}
\makeatother
\renewcommand{\theequation}{\arabic{section}.\arabic{equation}}
\renewcommand{\thetheorem}{\arabic{section}.\arabic{theorem}}
\begin{document}

\bigskip

\begin{center}\bf Государственный комитет Российской Федерации\\
по высшему образованию
\vskip20pt
\large

Московский государственный авиационный\\
технологический университет им.\ К.~Э.~Циолковского

\bigskip
\bigskip
 
 
Кафедра ``Высшая математика''

\vskip100pt

\vskip70pt
\Large
КВАДРАТУРНЫЕ ФОРМУЛЫ

\vskip50pt

\rm Методические указания по курсу ``Численные методы''
\vskip70pt

\begin{flushright}
Составитель: Осипенко К.Ю.
\end{flushright}

\vskip 90pt
Москва  1995 г.
\end{center}

\thispagestyle{empty}

\pagestyle{plain}

\section*{\bf Введение}

Если $f(x)$ --- непрерывная на отрезке $[a,b]$ функция и $F(x)$ --- ее
первообразная, то по формуле Ньютона--Лейбница
$$\int_a^bf(x)\,dx=F(a)-F(b).$$
Однако часто первообразная не может быть выражена через элементарные функции
или является слишком сложной. Например, при вычислении интегралов
$$\int_0^1e^{-x^2}\,dx,\qquad\int_0^1\frac{dx}{(1+x^2)^{10}}$$
в первом случае мы сталкиваемся с тем, что первообразная для $e^{-x^2}$ не
выражается через элементарные функции, а во втором --- с тем, что
первообразная для функции $(1+x^2)^{-10}$ является слишком громоздким
выражением.

Кроме того, на практике подынтегральная функция $f(x)$ обычно бывает задана в
дискретном числе точек. В этом случае первообразная $F(x)$ вообще не может
быть найдена точно. Тем самым возникает задача приближенного вычисления
определенного интеграла от функции по информации о значениях этой функции в
некоторой системе точек. Такого рода формулы называются {\it квадратурными
формулами}. В данном пособии рассматриваются основные методы построения
простейших квадратурных формул и оценки их погрешностей.

\section{Формула прямоугольников}

Рассмотрим задачу приближенного вычисления интеграла
$$\int_{-h/2}^{h/2}f(x)\,dx$$
по значению $f(0)=f_0$. Естественно считать, что функция всюду приближенно
равна $f_0$, и заменить вычисление интеграла от исходной функции на
вычисление интеграла от постоянной $f_0$. Таким образом, мы приходим к
простейшей квадратурной формуле
$$\int_{-h/2}^{h/2}f(x)\,dx\approx hf_0,$$
называемой {\it формулой прямоугольников}. Геометрический смысл формулы
прямоугольников заключается в том, что площадь под графиком функции $y=f(x)$
заменяется на площадь прямоугольника с высотой, равной $f(0)$.

Оценим погрешность формулы прямоугольников в предположении, что у функции $f(
x)$ существует непрерывная вторая производная. По формуле Тейлора имеем
$$f(x)=f_0+f_0'x+\frac{f''(\xi)}{2!}x^2,$$
где $f_0'=f'(0)$. Отсюда
$$\int_{-h/2}^{h/2}f(x)\,dx=hf_0+\int_{-h/2}^{h/2}\frac{f''(\xi)}{2!}x^2\,dx.
$$
Положим
$$M_2=\max_{x\in[-h/2,h/2]}|f''(x)|.$$
Тогда
\begin{equation}\label{1.1}
\left|\int_{-h/2}^{h/2}f(x)\,dx-hf_0\right|\le\frac{M_2}2\int_{-h/2}^{h/2}x
^2\,dx=\frac{M_2}{24}h^3.
\end{equation}

\section{Усложненная формула прямоугольников}

Пусть имеется отрезок $[a,b]$ и требуется приближенно вычислить
$$\int_a^bf(x)\,dx.$$
Разобьем отрезок $[a,b]$ на $N$ равных частей точками
$$x_i=a+ih,\quad i=0,1,\ldots,N,\quad h=\frac{b-a}N.$$
На каждом из отрезков $[x_i,x_{i+1}]$, $i=0,\ldots,N-1$, вычислим значение
функции в средней точке $x_{i+1/2}=a+(i+1/2)h$ и применим формулу
прямоугольников
\begin{equation}\label{2.1}
\int_{x_i}^{x_{i+1}}f(x)\,dx\approx hf_{i+1/2},
\end{equation}
где $f_{i+1/2}=f(x_{i+1/2})$.

Поскольку
$$\int_a^bf(x)\,dx=\sum_{i=0}^{N-1}\int_{x_i}^{x_{i+1}}f(x)\,dx,$$
то, сложив приближенные равенства \eqref{2.1}, получим
\begin{equation}\label{2.2}
\int_a^bf(x)\,dx\approx h(f_{1/2}+f_{3/2}+\ldots+f_{N-1/2}).
\end{equation}
Формула \eqref{2.2} называется {\it усложненной формулой прямоугольников\/}
(часто именно эту формулу называют {\it формулой прямоугольников\/}).

Для оценки погрешности формулы \eqref{2.2} предположим, что функция $f(x)$
имеет непрерывную вторую производную на отрезке $[a,b]$ и
$$M_2=\max_{x\in[a,b]}|f''(x)|.$$
Учитывая оценку \eqref{1.1}, имеем
$$\biggl|\int_a^bf(x)\,dx-h\sum_{i=0}^{N-1}f_{i+1/2}\biggr|\le\sum_{i=0}^{N-1
}\frac{M_2}{24}h^3=\frac{M_2}{24}(b-a)h^2.$$
В силу равенства $Nh=b-a$ можно переписать эту оценку следующим образом
$$\biggl|\int_a^bf(x)\,dx-h\sum_{i=0}^{N-1}f_{i+1/2}\biggr|\le\frac{M_2(b-a)
^3}{24N^2}.$$

\section{Использование интерполяционного многочлена Лагранжа для построения
квадратурных формул}

Одним из общих приемов построения квадратурных формул является замена
функции, заданной на отрезке $[a,b]$, некоторой более простой и в то же время
близкой к исходной функцией. Например, если $f(x)$ известна в некоторых
точках отрезка $[a,b]$ \ $x_0,x_1,\ldots,x_n$, то можно заменить ее на
интерполяционный многочлен Лагранжа $L_n(x)$ и положить
$$\int_a^bf(x)\,dx\approx\int_a^bL_n(x)\,dx.$$

Для интерполяционного многочлена Лагранжа имеет место равенство (см.\ [7,
стр.\ 9])
$$L_n(x)=\sum_{i=0}^nl_{ni}(x)f_i,$$
где
$$l_{ni}(x)=\frac{(x-x_0)\ldots(x-x_{i-1})(x-x_{i+1})\ldots(x-x_n)}{(x_i-x_0)
\ldots(x_i-x_{i-1})(x_i-x_{i+1})\ldots(x_i-x_n)},\quad f_i=f(x_i).$$
Таким образом, получаем квадратурную формулу
\begin{equation}\label{3.1}
\int_a^bf(x)\,dx\approx\sum_{i=0}^nA_if_i,
\end{equation}
в которой
\begin{equation}\label{3.2}
A_i=\int_a^bl_{ni}(x)\,dx.
\end{equation}

{\sc Замечание.} Если $f(x)$ --- многочлен степени $k\le n$, то в силу
единственности интерполяционного многочлена Лагранжа $f(x)\equiv L_n(x)$. Тем
самым квадратурная формула \eqref{3.1} точна на многочленах степени $k\le n$. В частности, при всех $k=0,1,\ldots,n$
\begin{equation}\label{3.3}
\int_a^bx^k\,dx=\sum_{i=0}^nA_ix_i^k.
\end{equation}

\begin{example}\rm
Построить квадратурную формулу вида
$$\int_{-1}^1f(x)\,dx\approx A_0f(-1/2)+A_1f(0)+A_2f(1/4).$$
\end{example}

{\sc Решение.} Можно было бы вычислять $A_0$, $A_1$ и $A_2$ по формулам \eqref{3.2}. Например,
$$A_0=\int_{-1}^1\frac{x(x-1/4)}{3/8}\,dx=\frac83\left.\left(\frac{x^3}3-
\frac{x^2}8\right)\right|_{-1}^1=\frac{16}9.$$
Мы же воспользуемся равенствами \eqref{3.3}. Имеем
\begin{align*}
\int_{-1}^1\,dx&=A_0+A_1+A_2,\\
\int_{-1}^1x\,dx&=A_0\left(-\frac12\right)+A_2\frac14,\\
\int_{-1}^1x^2\,dx&=A_0\left(-\frac12\right)^2+A_2\left(\frac14\right)^2.
\end{align*}
Таким образом, получаем следующую систему
$$\begin{cases}3A_0+A_1+A_2&=2,\\
-\dfrac12A_0+\dfrac14A_2&=0,\\
\dfrac14A_0+\dfrac1{16}A_2&=\dfrac23.\end{cases}$$
Решая эту систему, находим
$$A_0=\frac{16}9,\qquad A_1=-\frac{10}3,\qquad A_2=\frac{32}9.$$
Следовательно, искомая квадратурная формула имеет вид
$$\int_{-1}^1f(x)\,dx\approx \frac{16}9f(-1/2)-\frac{10}3f(0)+\frac{32}9f(1/4
).$$

\section{Квадратурные формулы Ньютона--Котеса}

Для приближенного вычисления интеграла
$$\int_a^bf(x)\,dx$$
разобьем отрезок $[a,b]$ на $n$ равных частей точками
$$x_i=x_0+ih,\quad i=0,1,\ldots,n,\quad h=\frac{b-a}n$$
($x_0=a$, $x_n=b$). Заменим функцию $f(x)$ на ее интерполяционный многочлен
Лагранжа, построенный по значениям этой функции в точках $x_0,x_1,\ldots,x_n$.
Для равноотстоящих узлов интерполяции (см.\ [7, стр.\ 12])
$$l_{ni}(x)=\frac{(-1)^{n-i}}{i!(n-i)!}q(q-1)\dots(q-i+1)(q-i-1)\dots(q-n),$$
где
\begin{equation}\label{4.1}
q=\frac{x-x_0}h.
\end{equation}
Сделав замену переменной \eqref{4.1} и учитывая, что $dx=h\,dq$, получим
квадратурную формулу \eqref{3.1}, в которой
$$A_i=h\frac{(-1)^{n-i}}{i!(n-i)!}\int_0^nq(q-1)\ldots(q-i+1)(q-i-1)\dots(q-n)
\,dq.$$

Положим
\begin{equation}\label{4.2}
H_i=\frac1n\frac{(-1)^{n-i}}{i!(n-i)!}\int_0^nq(q-1)\ldots(q-i+1)(q-i-1)
\dots(q-n)\,dq.
\end{equation}
Так как $h=\dfrac{b-a}n$, то $A_i=(b-a)H_i$. Тем самым получаем квадратурную
формулу
\begin{equation}\label{4.3}
\int_a^bf(x)\,dx\approx(b-a)\sum_{i=0}^nH_if_i,
\end{equation}
называемую {\it квадратурной формулой Ньютона--Котеса}. Величины \eqref{4.2}
называются {\it коэффициентами Котеса}. Для них существуют таблицы.

Отметим некоторые простейшие свойства коэффициентов Котеса. Поскольку
квадратурная формула \eqref{4.3} точна для любого многочлена степени $k\le n$, то она, в частности точна для $f(x)\equiv1$. Подставив эту функцию в \eqref{4.3}, находим
\begin{equation}\label{4.4}
\sum_{i=0}^nH_i=1.
\end{equation}
Кроме того, непосредственно из \eqref{4.2} вытекают равенства
\begin{equation}\label{4.5}
H_{n-i}=H_i,\quad i=0,1,\ldots,n-1.
\end{equation}

\section{Формула трапеций}

Рассмотрим квадратурную формулу Ньютона--Котеса при $n=1$. В этом случае $h=b-a$, а сама формула имеет вид
$$\int_a^bf(x)\,dx\approx h(H_0f_0+H_1f_1).$$
Пользуясь свойствами \eqref{4.4} и \eqref{4.5}, находим
\begin{gather*} 
H_0+H_1=1,\\
H_0=H_1.
\end{gather*}
Следовательно, $H_0=H_1=1/2$.

Итак, получена формула
$$\int_a^bf(x)\,dx\approx h\frac{f_0+f_1}2,$$
называемая {\it формулой трапеций}.

Оценим погрешность формулы трапеций. Для погрешности приближения функции
многочленом Лагранжа первой степени имеет место равенство (см.\ [7, стр.\ 10,
12])
$$f(x)-L_1(x)=\frac{f''(\xi)}2h^2q(q-1),$$
где $\xi$ некоторая точка из интервала $(a,b)$, а $q$ определено равенством
\eqref{4.1}. Таким образом, имеем
\begin{multline*}
\left|\int_a^bf(x)\,dx-h\frac{f_0+f_1}2\right|=\left|\int_a^b
(f(x)-L_1(x))\,dx\right|\\
=\left|\int_0^1\frac{f''(\xi)}2h^3q(q-1)\,dq\right|\le\frac{M_2}2h^3\int_0^1q
(1-q)\,dq=\frac{M_2}{12}h^3.
\end{multline*}

Аналогично усложненной формуле прямоугольников можно построить усложненную
формулу трапеций, разбив весь отрезок $[a,b]$ на $n$ равных частей и применив
на каждой из частей формулу трапеций. Тогда $h=\dfrac{b-a}n$, а
соответствующая квадратурная формула имеет вид
$$\int_a^bf(x)\,dx\approx h\left(\frac{f_0}2+f_1+\ldots+f_{n-1}+\frac{f_n}2
\right).$$
При этом для погрешности этой квадратурной формулы (мы обозначаем ее через $R_t$) будут справедливы неравенства
\begin{equation}\label{5.1}
R_t\le n\frac{M_2}{12}h^3=\frac{M_2}{12}(b-a)h^2=\frac{M_2(b-a)^3}{12n^2}.
\end{equation}

\begin{example}\label{p5.1}\rm
В скольких точках надо вычислить функцию $f(x)=e^{x^2}$, чтобы при вычислении интеграла
$$\int_0^1e^{x^2}\,dx$$
по усложненной формуле трапеций погрешность не превосходила $10^{-5}$.
\end{example}

{\sc Решение.} Поскольку $f''(x)=(2+4x^2)e^{x^2}$, то $M_2=6e$. Из \eqref{5.1}
вытекает, что для достижения требуемой точности достаточно выполнения
неравенства
$$\frac{M_2}{12n^2}\le10^{-5}.$$
Отсюда $n\ge100\sqrt{5e}=368,6\ldots\,\,$. Число вычислений функции при
разбиении отрезка на $n$ частей равно $n+1$, поэтому для вычисления
рассматриваемого интеграла с заданной точностью потребуется находить значения
функции в 370 точках.

\section{Формула Симпсона}

Пусть теперь $n=2$, т.е.\ $x_0=a$, $x_1=a+h$, $x_2=b$, а $h=\dfrac{b-a}2$.
Из \eqref{4.2} имеем
$$H_0=\frac14\int_0^2(q-1)(q-2)\,dq=\frac14\left.\left(\frac{q^3}3-\frac{3q^2
}2+2q\right)\right|_0^2=\frac16.$$
Тогда в силу \eqref{4.5} $H_2=H_0=1/6$, а из \eqref{4.4} вытекает, что $H_1=2/3$. Итак, мы получили квадратурную формулу
$$\int_a^bf(x)\,dx\approx\frac h3(f_0+4f_1+f_2),$$
которая называется {\it формулой Симпсона}. Геометрический смысл этой формулы
заключается в том, что вместо интеграла от исходной функции вычисляется
интеграл от параболы, проходящей через значения функции в точках $a$, $\dfrac{a+b}2$ и $b$.

Оценим погрешность формулы Симпсона. Без ограничения общности можно считать,
что начало координат выбрано в точке $\dfrac{a+b}2$. Тем самым будем
оценивать квадратурную формулу
\begin{equation}\label{6.1}
\int_{-h}^hf(x)\,dx\approx\frac h3(f(-h)+4f(0)+f(h)).
\end{equation}

Прежде всего отметим, что формула Симпсона точна не только на многочленах
второй степени, как следует из замечания в п.\ 3, но и на многочленах третьей
степени. Для того, чтобы убедиться в этом, достаточно проверить, что формула
\eqref{6.1} точна для $f(x)=x^3$.

Пусть $c$ --- произвольная точка из интервала $(0,h)$. Рассмотрим
интерполяционный многочлен Лагранжа $L_3(x)$, интерполирующий функцию $f(x)$
в точках $-h$, $0$, $c$ и $h$. Имеем
$$\int_{-h}^hL_3(x)\,dx=\frac h3(L_3(-h)+4L_3(0)+L_3(h))=\frac h3(f(-h)+4f(0)
+f(h)).$$
Поэтому для погрешности формулы Симпсона, которую мы обозначим через $R$,
будем иметь
\begin{multline*}
R=\biggl|\int_{-h}^hf(x)\,dx-\frac h3(f(-h)+4f(0)+f(h))
\biggr|\\
=\biggl|\int_{-h}^hf(x)\,dx-\int_{-h}^hL_3(x)\,dx\biggr|=\biggl|\int_{-h}^h(f
(x)-L_3(x))\,dx\biggr|\\
\le\int_{-h}^h|f(x)-L_3(x)|\,dx.
\end{multline*}
Положим
$$M_4=\max_{x\in[-h,h]}|f^{(4)}(x)|.$$
Тогда для погрешности приближения функции $f(x)$ многочленом Лагранжа
справедлива следующая оценка (см.\ [7, стр.\ 11])
$$|f(x)-L_3(x)|\le\frac{M_4}{4!}|(x+h)x(x-c)(x-h)|.$$
Тем самым для любого $c\in(0,h)$
$$R\le\frac{M_4}{24}\int_{-h}^h|(x+h)x(x-c)(x-h)|\,dx.$$
Переходя к пределу в этом неравенстве при $c\to0$, получаем
\begin{equation}\label{6.2}
R\le\frac{M_4}{24}\int_{-h}^hx^2(h^2-x^2)\,dx=\frac{M_4}{24}\left.\left(h^2
\frac{x^3}3-\frac{x^5}5\right)\right|_{-h}^h=\frac{M_4}{90}h^5.
\end{equation}

Выведем теперь усложненную формулу Симпсона. Разделим отрезок $[a,b]$ на
четное число отрезков $2n$
$$x_i=a+ih,\quad i=0,1,\ldots,2n,\quad h=\frac{b-a}{2n}.$$
На каждом из отрезков $[x_{2i},x_{2i+2}]$ применим формулу Симпсона
$$\int_{x_{2i}}^{x_{2i+2}}f(x)\,dx\approx\frac h3(f_{2i}+4f_{2i+1}+f_{2i+2}).$$
Сложив эти равенства, получим усложненную формулу Симпсона
$$\int_a^bf(x)\,dx\approx\frac h3\biggl(f_0+f_{2n}+4\sum_{i=1}^nf_{2i-1}+2
\sum_{i=1}^{n-1}f_{2i}\biggr).$$

Учитывая оценку \eqref{6.2}, для погрешности усложненной формулы Симпсона,
которую мы обозначим через $R_s$, будем иметь
$$R_s\le n\frac{M_4}{90}h^5=\frac{M_4}{180}(b-a)h^4,$$
где
$$M_4=\max_{x\in[a,b]}|f^{(4)}(x)|.$$
Поскольку $h=\dfrac{b-a}{2n}$, то оценка величины $R_s$ может быть дана через $n$
$$R_s\le\frac{M_4(b-a)^5}{2880n^4}.$$

\begin{example}\rm
Рассмотреть пример~\ref{p5.1} для формулы Симпсона.
\end{example}

{\sc Решение.} Имеем $f^{(4)}(x)=(12+48x^2+16x^4)e^{x^2}$. Тем самым $M_4=76e$. Для
достижения требуемой точности достаточно выполнения неравенства
$$\frac{M_4}{2880n^4}\le10^{-5}.$$
Отсюда
$$n\ge10\sqrt[4]{\frac{19e}{72}}=9,2\ldots\,\,.$$
При использовании формулы Симпсона отрезок разбивается на $2n$ частей.
Поэтому для вычисления рассматриваемого интеграла с заданной точностью
потребуется находить значения функции в 21 точке (напомним, что в аналогичном
примере для формулы трапеций при той же точности требовались вычисления в 370
точках).

\section{Главная часть погрешности квадратурных формул}

Для усложненной формулы прямоугольников было доказано следующее равенство
\begin{equation}\label{7.1}
\int_a^bf(x)\,dx=h\sum_{i=0}^{n-1}f_{i+1/2}+r_1,
\end{equation}
где $h=\dfrac{b-a}n$, а
$$|r_1|\le\frac{M_2}{24}(b-a)h^2.$$
Мы получим более точную оценку в предположении, что функция $f(x)$ имеет
непрерывную четвертую производную.

Рассмотрим сначала случай, когда $[a,b]=[-h/2,h/2]$. По формуле Тейлора
\begin{equation}\label{7.2}
f(x)=f_0+f'_0x+\frac{f_0''}{2!}x^2+\frac{f_0'''}{3!}x^3+\frac{f^{(4)}(\xi)}
{4!}x^4,
\end{equation}
где $f^{(i)}_0=f^{(i)}(0)$, $i=0,1,2,3$, а $\xi$ --- некоторая точка из
интервала $(-h/2,h/2)$. Интегрируя \eqref{7.2}, получаем равенство
$$\int_{-h/2}^{h/2}f(x)\,dx=f_0h+\frac{f_0''}{24}h^3+r_2,$$
в котором
$$r_2=\int_{-h/2}^{h/2}\frac{f^{(4)}(\xi)}{4!}x^4\,dx.$$
Для величины $r_2$ справедлива следующая оценка
$$|r_2|\le\frac{M_4}{4!}\int_{-h/2}^{h/2}x^4\,dx=\frac{M_4}{1920}h^5.$$
Таким образом, для усложненной формулы прямоугольников имеем
\begin{equation}\label{7.3}
\int_a^bf(x)\,dx=h\sum_{i=0}^{n-1}f_{i+1/2}+\frac{h^3}{24}\sum_{i=0}^{n-1}f
''_{i+1/2}+r_3,
\end{equation}
где
$$|r_3|\le n\frac{M_4}{1920}h^5=\frac{M_4}{1920}(b-a)h^4.$$
Применяя \eqref{7.1} для второй производной, находим
$$\int_a^bf''(x)\,dx=h\sum_{i=0}^{n-1}f''_{i+1/2}+r_4,$$
при этом
$$|r_4|\le\frac{M_4}{24}(b-a)h^2.$$
Следовательно,
$$h\sum_{i=0}^{n-1}f''_{i+1/2}=\int_a^bf''(x)\,dx-r_4.$$
Подставляя это выражение в \eqref{7.3}, получаем
$$\int_a^bf(x)\,dx=h\sum_{i=0}^{n-1}f_{i+1/2}+\frac{h^2}{24}\int_a^bf''(x)\,d
x+r_5,$$
где
$$|r_5|=\Bigl|-\frac{h^2}{24}r_4+r_3\Bigr|\le\frac{M_4}{24^2}(b-a)h^4+\frac{M
_4}{1920}(b-a)h^4=\frac{13M_4}{5760}(b-a)h^4.$$

{\sc Определение.} Говорят, что $\varphi(h)=O(h^k)$ (читается ``о большое''),
если существует такая постоянная $C>0$, для которой
$$|\varphi(h)|\le Ch^k.$$

Тем самым мы доказали равенство
\begin{equation}\label{7.4}
\int_a^bf(x)\,dx=h\sum_{i=0}^{n-1}f_{i+1/2}+ch^2+O(h^4),
\end{equation}
в котором
$$c=\frac1{24}\int_a^bf''(x)\,dx.$$
Положим
$$I=\int_a^bf(x)\,dx,\quad I_h^r=h\sum_{i=0}^{n-1}f_{i+1/2}.$$
Тогда равенство \eqref{7.4} запишется в виде
$$I=I_h^r+ch^2+O(h^4).$$
Величина $ch^2$ называется {\it главной частью\/} погрешности формулы
прямоугольников.

Аналогичные равенства можно получить для формул трапеций и Симпсона
\begin{align*}
I&=I_h^t+c_1h^2+O(h^4),\\
I&=I_h^s+c_2h^4+O(h^6)
\end{align*}
(в последнем случае надо требовать, чтобы функция $f(x)$ имела непрерывную
шестую производную).

\section{Правило Рунге практической оценки погрешности}

Пусть $z$ --- неизвестное точное значение некоторой величины, $z_h$ ---
известное ее приближенное значение, зависящее от положительного параметра
$h$, который может принимать сколь угодно малые значения.

Предположим, что установлена связь между точным и приближенным значениями
\begin{equation}\label{8.1}
z=z_h+ch^k+O(h^{k+m}),
\end{equation}
где $c$ --- неизвестная не зависящая от $h$ постоянная. Тогда
\begin{equation}\label{8.2}
z=z_{h/2}+c\left(\frac h2\right)^k+O(h^{k+m}),
\end{equation}
так как для любой постоянной $C$ \ $O((Ch)^n)=O(h^n)$.
Вычитая из \eqref{8.1} равенство \eqref{8.2}, будем иметь
\begin{equation}\label{8.3}
z_{h/2}-z_h=c\left(\frac h2\right)^k(2^k-1)+O(h^{k+m}).
\end{equation}
Отсюда
$$c\left(\frac h2\right)^k=\frac{z_{h/2}-z_h}{2^k-1}+O(h^{k+m}).$$
Следовательно, при $c\ne0$ величина $\dfrac{z_{h/2}-z_h}{2^k-1}$ отличается
от главного члена погрешности $z-z_{h/2}$ на $O(h^{k+m})$, т.е.
\begin{equation}\label{8.4}
z-z_{h/2}=\frac{z_{h/2}-z_h}{2^k-1}+O(h^{k+m}).
\end{equation}

Тем самым при $c\ne0$ оценить погрешность можно так:
$$z-z_{h/2}\approx\frac{z_{h/2}-z_h}{2^k-1}.$$
Такой способ оценки погрешности называется {\it правилом Рунге}.

{\sc Замечание.} На практике считается, что условие $c\ne0$ выполнено, если
\begin{equation}\label{8.5}
\left|2^k\frac{z_{h/2}-z_h}{z_h-z_{2h}}-1\right|<0,1.
\end{equation}
Только в этом случае рекомендуется применение правила Рунге.

Поясним условие \eqref{8.5}. Из \eqref{8.3} следует, что
\begin{align*}
z_{h/2}-z_h&=c\left(\frac h2\right)^k(2^k-1)+O(h^{k+m}),\\
z_h-z_{2h}&=ch^k(2^k-1)+O(h^{k+m}).
\end{align*}
Поэтому
$$\frac{z_{h/2}-z_h}{z_h-z_{2h}}\approx\frac1{2^k}.$$
Таким образом,
$$2^k\frac{z_{h/2}-z_h}{z_h-z_{2h}}\approx1.$$

\section{Уточнение приближенного решения по Ричардсону}

Положим
$$z_h^*=\frac{2^kz_{h/2}-z_h}{2^k-1}.$$
Тогда из \eqref{8.4} получаем
$$z=z_{h/2}+\frac{z_{h/2}-z_h}{2^k-1}+O(h^{k+m})=z_h^*+O(h^{k+m}).$$
При $c\ne0$ \ $z-z_{h/2}$ имеет $k$-ый порядок малости, а $z-z_h^*$ --- $(k+m)$-ый порядок малости, т.е.\ $z_h^*$ --- более точное приближение. Оно носит название {\it уточнение по Ричардсону}.

Таким образом, если вычисляется определенный интеграл по усложненным формулам
прямоугольников или трапеций, то $k=2$, и мы можем оценить приближенно
погрешность по правилу Рунге
$$I-I_{h/2}\approx\frac{I_{h/2}-I_h}3.$$
Кроме того, можно найти уточнение по Ричардсону
$$I_h^*=\frac{4I_{h/2}-I_h}3.$$

Для формулы Симпсона $k=4$, и приближенная оценка погрешности имеет вид
$$I-I^s_{h/2}\approx\frac{I_{h/2}^s-I_h^s}{15}.$$
Для уточнения по Ричардсону имеем
$${I_h^s}^*=\frac{16I_{h/2}^s-I_h^s}{15}.$$

{\sc Замечание.} Формулы трапеций и Симпсона удобны тем, что при переходе от $h$ к
$h/2$ все вычисленные ранее значения функций используются в новой квадратурной формуле.

Правило Рунге и уточнение по Ричардсону можно применять и для других задач приближенного вычисления. Рассмотрим в качестве примера численное дифференцирование.

Пусть известны значения некоторой достаточно гладкой функции $f(x)$ в точках
$x_i=x_0+ih$, $i=0,\pm1$. Положим $f_i=f(x_i)$, $i=0,\pm1$. Хорошо известна
формула численного дифференцирования (см., например, [7, стр.\ 26])
$$f_0''\approx\frac{f_1-2f_0+f_{-1}}{h^2}=f''_{0h}.$$
Найдем главную часть погрешности в этом методе. По формуле Тейлора имеем
$$f_{\pm1}=f_0\pm f_0'h+\frac{f_0''}{2!}h^2\pm\frac{f_0'''}{3!}h^3+\frac{f_0^
{(4)}}{4!}h^4\pm\frac{f_0^{(5)}}{5!}+\frac{f^{(6)}(\xi_{\pm1})}{6!}h^6,$$
где $f_0^{(i)}=f^{(i)}(0)$, $i=0,1,\ldots,5$, $\xi_1\in(0,h)$, а $\xi_{-1}\in(-h,0)$.
Положим
$$r=-\frac{f^{(6)}(\xi_1)+f^{(6)}(\xi_{-1})}{6!}h^4.$$
Тогда
$$f_1+f_{-1}=2f_0+f_0''h^2+\frac{f_0^{(4)}}{12}h^4-rh^2.$$
Тем самым
\begin{equation}\label{9.1}
f_0''=\frac{f_1-2f_0+f_{-1}}{h^2}-\frac{f_0^{(4)}}{12}h^2+r.
\end{equation}
Поскольку
$$|r|\le\frac{2M_6}{6!}h^4,$$
то равенство \eqref{9.1} может быть записано в виде
$$f_0''=f''_{0h}+ch^2+O(h^4),$$
где $c=-f_0^{(4)}/12$. Следовательно, по правилу Рунге
$$f_0''-f''_{0,h/2}\approx\frac{f''_{0,h/2}-f''_{0h}}3,$$
а для уточнения по Ричардсону справедливо равенство
$${f''_{0h}}^*=\frac{4f''_{0,h/2}-f''_{0h}}3.$$

\section*{\bf Литература}

\begin{itemize}
\item[1.] Бахвалов Н.С. {\it Численные методы}. М.: Наука, 1973.
\item[2.] Бахвалов Н.С., Жидков Н.П., Кобельков Г.М. {\it Численные методы}.
М.: Наука, 1987.
\item[3.] Березин И.С., Жидков Н.П. {\it Методы вычислений}. М.: Наука, 1966.
Т.1; Физматгиз, 1962. Т.2.
\item[4.] Волков Е.А. {\it Численные методы}. М.: Наука, 1982.
\item[5.] Демидович Б.П., Марон И.А. {\it Основы вычислительной математики}.
М.: Наука, 1966.
\item[6.] Калиткин Н.Н. {\it Численные методы}. М.: Наука, 1978.
\item[7.] Осипенко К.Ю. {\it  Аппроксимация функций многочленами и численное
дифференцирование}: Методические  указания  по курсу ``Численные методы'';
МГАТУ. М., 1994.
\end{itemize}

\tableofcontents
\end{document}
