Глава 1. Интерполирование 3 страница

где ; .

При грубой оценке погрешности на практике согласно формуле (1.26) можно положить

,

где - среденее арифметическое всех вычисленных конечных разностей порядка (2m+1).

Замечание 1. Если при практических расчетах возникает необходимость в построении интерполяционного многочлена нечетной степени, его можно получить из формулы (1.34) путем отбрасывания последнего слагаемого:

(1.36)

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

(1.37)

Замечание 2. При выводе интерполяционной формулы Гаусса (1.34) предполагалось, что точка интерполирования х лежит ближе к x_k, чем к x_k+1. Если заданная точка х лежит ближе к x_k+1, можно вывести формулу, аналогичную формуле (1.34), используя узлы интерполяции в следующем порядке:

x_k+1, x_k, x_k+2, x_k-1, x_k+3, ….

При отсутствии новой формулы можно воспользоваться формулой (1.34).

Пример. Для функции, заданной таблицей своих значений в примере параграфа 1.7, вычислить приближенно значение в точке х=1,22.

Точка х=1,22 лежит в середине таблицы между узлами x_k=1,2 и x_k+1=1,3.

Воспользуемся интерполяционной формулой Гаусса (1.34) при

m=2; x_k=1,2; h=0,1; .

Значения берем из табл. 1.3:

;

.

1.9. Обратное интерполирование

Пусть функция y=f(x) задана таблицей своих значений . В предыдущих параграфах рассмотрено, как по заданному значению аргумента приближенно найти значение функции. зАдача обратного интерполирования заключается втом, чтобы по заданному значению функции y* найти аргумент x*, при котором . Предположим, что на отрезке [a,b], содержащем узлы интерполяции, функция f(x) монотонна. Тогда существует однозначная обратная функция x=F(y). Она задана той же таблицей значений, что и функция y=f(x), только в роли значений аргумента выступают значения , а - соответствующие им значения функции. Обратное интерполирование в этом случае свелось к обычному интерполировнию для функции x=F(y). Строим интерполяционный многочлен , воспользовавшись, например, интерполяционной формулой Лагранжа. Подставляя в заданное значение y*, получим .

Рассмотрим второй способ нахождения х*, который применим ко всякой функции f(x). Не меняя ролями функцию и аргумент, записываем интерполяционный многочлен по той или иной формуле. Неизвестное значение х* находим приближенно, решая уравнение =у*. Если число узлов велико, то этот способ нахождения х* приводит к решению алгебраического уравнения высокой степени. Различные методы решения таких уравнений будут рассмотрены в гл. 3. Ниже остановимся на методе итерации. Будем рассматривать только случай равноотстоящих значений аргумента, т.е.

.

Пусть для определенности у* находится между у₀ и у₁. Строим интерполяционный многочлен по первой интерполяционной формуле Ньютона. Уравнение =у* принимает вид

. (1.38)

Разрешим это уравнение относительно q, стоящего при разности первого порядка. Получим

(1.39)

Метод итерации (метод последовательных приближений) заключается в следующем. Выбираем начальное приближение .

Подставляя в правую часть равенства (1.39), получаем

Подставляя вычисленное в правую часть равенства (1.39), получаем

Аналогично по находим , а затем и т.д.:

. (1.40)

В значительном числе случаев этот процесс сходится и дает в пределе точное решение уравнения.

На практике интерполяционный процесс заканчивают, когда два соседних приближения совпадают в пределах заданной точности. Последнее вычисленное приближение принимается за решение уравнения (1.38) – q*.

Найдя q*, определяем x* из соотношения

,

т.е. .

Пример. Функция y=f(x) задана таблично:

Х	0,5	0,6	0,7	0,8	0,9
у	1,6487	1,8221	2,0138	2,2255	2,4596

Найти значение аргумента х*, для которого у*=1,7333. Строим табл. 1.5 конечных разностей заданной функции.

Таблица 1.5

0,5	1,6487
0,6	1,8221
0,7	2,0138
0,8	2,2255
0,9	2,4596

Так как , для нахождения q* воспользуемся уравнением (1.39). подставляем в (1.39) у*, у₀ и конечные разности из табл. 1.5:

;

.

Подставляем найденное в правую часть уравнения (1.41):

.

Аналогично находим

.

Четыре знака после запятой у и совпали. Отсюда

;

.

1.10. Численное дифференцирование

Пусть функция f(x) задана таблицей своих значений на отрезке [a,b]. Требуется вычислить производную порядка m в некоторой точке х отрезка [a,b]. Задача эта может быть решена следующим образом. По заданным узлам интерполяции строим интерполяционный многочлен . Тогда

, (1.42)

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

Дифференцируем тождество (1.42) m раз в предположении, что f(x) имеет производную порядка m, и получаем

.

Пренебрегая величиной , получим формулу для приближенного вычисления производной:

. (1.43)

Погрешность этой формулы равна .

Отметим, что численное дифференцирование – операция менее точная, чем интерполирование.

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

, (1.44)

где .

Дифференцируем (1.44) по переменной х. Отметим, что

.

После дифференцирования получаем

. (1.45)

Аналогично находим

(1.46)

Таким же способом в случае надобности можно вычислить и производные функции f(x),более высокого порядка. Если точка х лежит в середине или в конце таблицы, формулы для вычисления производных получим из интерполяционной формулы Гаусса или второй интерполяционной формулы Ньютона соответственно. Формулы численного дифференцировния упрощаются, если требуется вычислить производные в узлах интерполяции. Пусть х=х₀. Тогда q=0. Формулы (1.45) и(1.46) запишутся в виде

;

.

Погрешность в вычислении первой производной в точке х₀ равна

.

Аналогично может быть найдена погрешность для второй производной.

Глава 2. ЧИСЛЕННОЕ ИНТЕГРИРОВАНИЕ

2.1. Общие замечания

Если функция f(x) непрерывна на отрезке [a,b], то определенный интеграл от этой функции в пределах от а до b может быть вычислен по формуле Ньютона – Лейбница

,

где F(x) – первообразная функция f(x).

На практике, как уже указывалось во введении, часто встречаются функции, для которых первообразная не выражается через элементарные функции, а также функции, заданные таблично. вСе это приводит к необходимости применения численных методов. Для получения формул приближенного вычисления интеграла подынтегральную функцию на отрезке [a,b] заменяют интерполирующей функцией (чаще всего интерполяционным многочленом)

.

Тогда

= + .

Отбрасывая второе слагаемое полученного соотношения, приближенно полагаем

. (2.1)

Погрешность равенства (2.1) равна , где - погрешность интерполирования.

Интерполяционный многочлен может быть записан в различных видах в зависимости от применяемой интерполяционной формулы и степени m. После подстановки в (2.1) получаем различные формулы для приближенного вычисления определенного интеграла, которые называются квадратурными формулами.

2.2. Формула трапеций

Заменим подынтегральную функцию интерполяционным многочленом первой степени, построенным по двум узлам интерполяции и . Значения подынтегральной функции в узлах интерполяции равны соответственно и

.

Подставляем в (2.1):

.

Итак,

. (2.2)

Получена формула трапеций. Геометрический смысл формулы (2.2) отображен на рис. 2.1.

B

y=f(x)

у

0 a b x

Рис. 2.1

Левая часть формулы – площадь криволинейной трапеции, правая часть – площадь трапеции аАВb.

Оценим погрешность формулы трапеций. Предположим, что f’’(x) непрерывна на отрезке [a,b]. Имеем

= , (2.3)

где .

К полученному интегралу применим обобщенную теорему о среднем, которая формулируется следующим образом.

Пусть функция непрерывна, а функция q(x) интегрируема на отрезке [a,b]. Пусть, кроме этого, (или ) для . Тогда найдется точка , такая, что

.

Подынтегральная функция интеграла (2.3) представляет собой произведение двух функций: одна функция непрерывна, а на отрезке . После применения обобщенной теоремы о среднем получаем

, (2.4)

где .

Погрешность R пропорциональна (b-a)³, поэтому, если длина отрезка большая, погрешность формулы трапеции также большая. Для увеличения точности формулы трапеций разделим отрезок на n равных частей длины .

Точки деления , , …, .

Применим формулу трапеций (2.2) к каждому отрезку :

. (2.5)

Погрешность этой формулы равна

,

где .

Сумма интегралов (2.5) по всем отрезкам дает общую квадратурную формулу трапеции

. (2.6)

Погрешность этой формулы

.

Величина есть среднее арифметическое, составленное из n значений второй производной; оно лежит где-то между этими значениями. Вторую производную мы предполагали непрерывной на отрезке [a,b], поэтому она принимает все промежуточные значения. Отсюда вытекает, что на отрезке [a,b] найдется такая точка , что

и

(2.7)

Формула трапеций является точной, если подынтегральная функция – многочлен степени, не выше первой.

2.3. Формула Симпсона

Заменим подынтегральную функцию интерполяционным многочленом второй степени, построенным по трем узлам интерполяции:

, , ;

.

Интерполяционный многочлен запишем при помощи первой интерполяционной формулы Ньютона

,

где

; ;

=

= . (2.8)

Полученная формула называется формулой парабол, или формулой Симпсона (Т. Симпсон (1710 – 1761) – английский математик).

Эта формула является точной для любого многочлена второй степени (погрешность интерполирования в этом случае равна нулю). Кроме того, формула Симпсона является точной для (левая и правая части формулы (2.8) обращаются в этом случае в нуль). Отсюда вытекает, что формула Симпсона является точной для любого многочлена третьей степени. Погрешность формулы Симпсона для произвольной подынтегральной функции записывается в следующем виде:

. (2.9)

Предполагается, что непрерывна на [a,b]. Для увеличения точности формулы Симпсона (2.8) разделим отрезок [a,b] на n равных частей длины , причем возьмем n четное, т.е. n=2m. Рассмотрим сдвоенный частичный отрезок и применим к нему формулу Симпсона (2.8):

(2.10)

Погрешность этой формулы, согласно (2.9) равна

= . (2.11)

Всего сдвоенных частичных отрезков .

Сумма интегралов (2.10) по всем сдвоенным отрезкам дает общую формулу Симпсона

. (2.12)

Погрешность этой формулы

= .

В силу непрерывности на отрезке [a,b] на этом отрезке существует такая точка , что

= .

Следовательно, имеем

(2.13)

2.4. Оценка погрешности квадратурных формул

Если подынтегральная функция f(x) задана аналитически и нахождение (для формулы трапеций) или (для формулы Симпсона) не вызывает затруднений, погрешность квадратурной формулы вычисляется по формуле (2.7) или (2.13). но на практике найти и оценить, например, производную четвертого порядка трудно, особенно в тех случаях, когда функция f(x) задана таблично. Поэтому часто прибегают к следующему приему грубой оценки погрешности, предложенному Рунге. Погрешность формулы трапеций может быть записана в виде

Обозначим

; .

Тогда ;

.

Объединяя эти результаты, можно записать, что погрешность квадратурной формулы

,

где k, N – постоянные;

h – длина частичного отрезка при разбиении отрезка [a,b] на n равных частей.

Обозначим: I – точное значение интеграла, I_h – приближенное его значение, вычисленное по квадратурной формуле с шагом h.

Тогда

. (2.14)

Вычислим этот же интеграл по квадратурной формуле с шагом 2h:

. (2.15)

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

Приравниваем правые части равенств (2.14) и (2.15):

.

Отсюда

;

Итак, погрешность квадратурной формулы при вычислении интеграла с шагом h удовлетворяет неравенству

.

Для формулы трапеций к=2, отсюда

. (2.16)

Для формулы Симпсона к=4, отсюда

. (2.17)

Пример. Вычислить приближенно интеграл по формуле трапеций и по формуле Симпсона при n=8. Оценить погрешность. Вычисления оформить в виде табл. 2.1.

Таблица 2.1

			Формула трапеций	Формула Симпсона
	0.125 0.250 0.375 0.500 0.625 0.750 0.875	1.000244 1.003906 1.019775 1.062500 1.152588 1.316406 1.586182
			19.283202	9.765624	28.800780	14.406248
I			Ih=1.2052	I2h=1.22	Ih=1.2000	I2h=1.200

Рассмотрим, как заполняются столбцы, соответствующие вычислению интеграла по формуле трапеций