Получения квадратурных формул

Итак, квадратурные формулы имеют вид

,

причем, как было установлено, если f (x) – полином степени £ n, то эта формула точна (то есть R= 0).

Используем это обстоятельство для нахождения коэффициентов квадратурной формулы и ее погрешности.

Будем поочередно подставлять полиномы вида[3]

, где 0 £ i,k £ n.

i= 0:	P 0(x) = 1;
i= 1:	P 1(x) =x–xk;
¼	¼	¼	¼
i=n:	Pn (x) = (x–xk) n;

Получили систему n+ 1 линейных алгебраических уравнений для нахождения коэффициентов С_i.

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

Продолжаем процесс подстановки полиномов дальше, но уже в квадратурную формулу с погрешностью R ¹ 0:

i=n+ 1:	Pn+ 1(x) = (x–xk) n+ 1;

Так как коэффициенты C_i уже определены, то из этого уравнения определяем – погрешность квадратурной формулы для полинома вида

P_{n+ 1}(x) = (x–x_k) ^{n+ 1}.

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

P_{n+ 2}(x) = (x–x_k) ^{n+ 2}; P_{n+ 3}(x) = (x–x_k) ^{n+ 3} ;...; P_m (x) = (x–x_k) ^m,

до тех пор, пока не станет отличной от нуля.

Погрешность квадратурной формулы для произвольной подынтегральной функции находится по формуле

.

Пример 6. Получить формулу численного интегрирования на элементарном отрезке по двум точкам:

.

y= 1:	;	C 0+ C 1= h;
y=x–xi:	;	.

Таким образом, С ₀ =C ₁ = h/ 2, и искомая формула имеет вид:

(формула трапеций).

Найдем погрешность.

y= (x–xi)2:

.

Отсюда

. .

Пример 7. Получить формулу численного интегрирования на элементарном отрезке по трем точкам:

.

y= 1:
y= (x–xi– 1/2):
y= (x–xi– 1/2)2:

Таким образом, система уравнений для искомых коэффициентов имеет вид:

C ₀ +C ₁ +C ₂ = h;

–C ₀ +C ₂ = 0;

C ₀ +C ₂ =h/ 3.

Решая ее, получаем: C ₀ =h/ 6; C ₁ = 4 h/ 6; C ₃ =h/ 6.

Найдем погрешность.

y= (x–xi– 1/2)3:

.

Отсюда , и процесс подстановки полиномов продолжаем дальше:

y= (x–xi– 1/2)4:

.

Погрешность , а

.

Таким образом, получили формулу Симпсона:

.

Пример 8. Получить формулу численного интегрирования, не являющуюся квадратурной в смысле данного определения (п. 3.5):

.

В правой части стоят производные в узловых точках – в этом отличие данной формулы численного интегрирования от квадратурной формулы интерполяционного типа (см. п. 3.8). Пусть x ₁ –x ₀ =h. Используем метод неопределенных коэффициентов.

y= 1:	h=C 0 +C 1;	Решение:
y=x–x 0:	h 2/2 =C 1 h+D 0 +D 1;	С 0 =h /2; C 1 =h /2;
y= (x–x 0)2:	h 3/3 =C 1 h2+ 2 D 1 h;	D 0 =h 2/12; D 1 = –h 2/12.
y= (x–x 0)3:	h 4/4 =C 1 h 3 + 3 D 1 h 2.

Погрешность:

y=(x–x ₀ ) ⁴: h ⁵ / 5 =C ₁ h ⁴ + 4 D ₁ h ³ + ; =h ⁵ / 30; .

Следовательно, искомая формула (формула Эйлера) имеет вид:

. (3.11)

Запишем составную формулу Эйлера:

.

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

Более того, можно показать, что при численной реализации этой формулы для аппроксимации производных (не снижая точности формулы Эйлера!) можно использовать формулы второго порядка точности:

3.8 КВАДРАТУРНЫЕ ФОРМУЛЫ НЬЮТОНА–КОТЕСА

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

, (3.12)

где r (х)>0 – заданная интегрируемая функция, называемая весовой функцией.

Если расположение узлов на отрезке интегрирования равномерное, то квадратурные формулы интерполяционного типа называются квадратурными формулами Ньютона-Котеса.

Таким образом, формулы прямоугольников, трапеций, Симпсона являются простейшими частными случаями формул Ньютона-Котеса при r (х)º1.

Основное свойство всех формул Ньютона-Котеса – они точно интегрируют полиномы до степени n включительно. Этот факт мы использовали при работе с методом неопределенных коэффициентов.

Справедливо и обратное утверждение: если квадратурная формула

точна для любого полинома степени n, то она является формулой Ньютона-Котеса (или шире – квадратурной формулой интерполяционного типа).

Получим явные выражения для коэффициентов формул Ньютона-Котеса. Это можно сделать простым интегрированием (вот еще один способ построения квадратурных формул) выражения (см. п. 3.5):

.

Заменим переменную x=a+t×h, где t= 0, 1 ,..., n.

Тогда:

1) узловые точки x ₀ =a; x ₁ =a+h; x ₂ =a+2h;...; x_n=a+nh;

2) бином x–x_i=a+t×h–a–i×h=h(t–i);

3) полином

4) производная

Подставляя все эти выражения в формулу для коэффициентов С_i, имеем:

. (3.13)

Пример 9. Получить коэффициенты квадратурной формулы Ньютона-Котеса при r (х)º1, n =2.

Это коэффициенты квадратурной формулы Симпсона

.

Погрешность формул Ньютона-Котеса оценивается интегралом от соответствующей погрешности интерполяционного полинома. Используя ту же самую замену переменных, получаем:

.

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

Примечание. Формулы Ньютона-Котеса высших порядков.

При n = 1, 2, 3, 4, 5, 6, 7, 9 все коэффициенты C_i (3.13) положительны, а при n =8 и n ³10 среди них имеются как положительные, так и отрицательные. По этой причине формулу Ньютона-Котеса не рекомендуется применять при больших n.

Конец примечания.

3.9 КВАДРАТУРНЫЕ ФОРМУЛЫ ГАУССА

Рассматривая квадратурные формулы интерполяционного типа, можно сделать вывод, что повышение их точности связано с увеличением количества узловых точек. Зададимся вопросом: нельзя ли повысить точность квадратурной формулы не изменяя количества узлов, а лишь перераспределяя их на заданном отрезке? Ранее мы видели, что можно минимизировать погрешность интерполяционного полинома, выбирая в качестве узлов корни полинома Чебышева. Поэтому есть надежда, что и здесь за счет отказа от равномерного расположения узлов можно получить квадратурные формулы, которые будут точны для полиномов степени выше, чем n.

Поставим задачу так: построить квадратурную формулу

,

которая при заданном n была бы точна для алгебраических полиномов возможно больших степеней m>n. [4]

Воспользуемся методом неопределенных коэффициентов (см. п. 3.7), но будем считать неопределенными не только коэффициенты, но и узлы.

Будем поочередно подставлять в это соотношение полиномы вида

f (x) =x^a, где a= 0, 1, 2 ,..., m.

Получим нелинейную систему m +1 уравнений относительно 2 n неизвестных С ₁, С ₂ ,..., С_n, x ₁, x ₂ ,..., x_n:

. (3.14)

Чтобы число уравнений равнялось числу неизвестных, необходимо потребовать m+ 1 = 2 n. Отсюда m= 2 n –1 – искомая наивысшая степень алгебраического полинома. Решая эту систему, находим неизвестные С ₁, С ₂ ,..., С_n, x ₁, x ₂ ,..., x_n.

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

полиномы степени выше, чем m.

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

Пример 10. Рассмотрим частный случай. Пусть r (x) º 1; a=– 1; b= 1; n= 3. Получим квадратурную формулу вида

.

f (x) = 1;	2 =C (1 + 1 + 1)	C=2/3;
f (x) =x;	0 =C ();	;
f (x) =x 2;	2 / 3 =C ();	;
f (x) =x 3;	0=C ();	.

Решая систему нелинейных уравнений[5], найдем

Для нахождения погрешности подставим следующую по порядку функцию f (x).

f (x) =x 4;

;

.

Тогда .

Таким образом, искомая формула имеет вид:

.

Такие формулы называются квадратурными формулами Чебышева. Их общий вид:

. (3.15)

При n =8 и n >9 узловые точки x_i принимают комплексные значения, поэтому квадратурные формулы Чебышева применимы только для n =1, 2, 3, 4, 5, 6, 7, 9.

Квадратурную формулу Чебышева можно применить и к интегралам с произвольными отрезками интегрирования a и b, сделав замену переменных

.

Эта замена переводит отрезок x Î[ a;b ] в отрезок t Î[–1;1]. Тогда

.

В общем же случае решение системы (3.14) довольно затруднительно. На помощь приходит теорема, которая дает рекомендации по построению квадратурных формул Гаусса.

Теорема. Квадратурная формула

(3.16)

точна для любого полинома степени m =2 n– 1 тогда и только тогда, когда выполняются два условия:

1) полином w (x) = (x–x ₁)(x–x ₂) ... (x–x _n), составленный по узловым значениям квадратурной формулы x ₁, x ₂, …, x_n, xÎ [ a, b ], ортогонален с весом r (х) любому полиному q (x)степени меньшей, чем n, то есть

; (3.17)

2) формула (3.17) является квадратурной формулой интерполяционного типа, то есть

(i= 1 ...n). (3.18)

Доказательство.

Необходимость. Пусть формула (3.16) точна для любого полинома степени m= 2 n– 1. Тогда она будет точна и для полинома w (x) q (x), так как его степень не выше 2 n– 1 (полином w (x) имеет степень n, а q (x) имеет степень не больше, чем n– 1). Поэтому

,

так как w (x_i) = 0.

Формула (3.17) доказана. А о справедливости (3.18) мы упоминали ранее (п.3.8).

Достаточность. Пусть f (x) – любой полином степени 2 n– 1. По теореме о делении полиномов его всегда можно представить в виде

f (x) =w (x) q (x) +r (x),

где r (x) – полином степени не выше n– 1.

Тогда

Теорема доказана.

Практический вывод из этой теоремы следующий: чтобы обеспечить наивысший порядок точности квадратурной формулы, необходимо найти систему ортогональных на отрезке [ a;b ] с весом r (х) полиномов Р_n (х) и их корни взять в качестве n узлов х ₁ ... х_n квадратурной формулы. Затем найти коэффициенты квадратурной формулы: либо методом неопределенных коэффициентов, либо непосредственным интегрированием:

.

Для погрешности квадратурных формул Гаусса справедлива формула

.

3.10 НЕКОТОРЫЕ ЧАСТНЫЕ СЛУЧАИ КВАДРАТУРНЫХ

ФОРМУЛ ГАУССА

3.10.1 Формула Гаусса-Чебышева [6]

Известно, что полиномы Чебышева T_n (x)=cos(n ×arccos(x)) обладают свойством ортогональности с весом , (–1£ x £1).

T ₀(x) = 1;

T ₁(x) =x;

T ₂(x) = 2 x ² – 1;

T ₃(x) = 4 x ³ –x;

¼

T_{n+ 1}(x) = 2 xT_n (x)– T_{n– 1}(x).

Поэтому справедлива квадратурная формула

, (3.19)

где x_i – корни полинома Чебышева:

Для коэффициентов С_i,определяя их по формуле

,

получим, что при любом i и n C_i = p/n.

Погрешность формулы можно определить двумя способами: либо методом неопределенных коэффициентов, либо по формуле

.

Таким образом, например, при n= 3 формула Гаусса-Чебышева выглядит следующим образом:

.