Вычислительная математика лекции

Отрезок [a,b] разбит на n элементарныx отрезков [xi-1,xi], i=1,…,n точками

n	xi
a=x0 < x1 < … < xn=b . Тогда I Ii , где	Ii	f ( x)dx,i 1, n .
i 1	xi 1

При численном расчёте Ii применяют одну из возможных квадратурных формул. Полученная таким образом конструкция называется составной квадратурной формулой.

8.3.1.Простейшие составные квадратурные

формулы.
Введём обозначения fi =f(xi), fi-1/2	= f(xi-1/2), где x		xi xi 1	.
	i 1/ 2	2

Для простоты принимаем h = xi – xi-1 = const. Однако,

результаты легко обобщаются для переменного шага.

51

8.3.1.1. Составная формула левых прямоугольников.

f

									x
		Xi-1		Xi
Ii hfi 1
				n 1
I h( f0 f1		... fn 1 ) h fi
				i 0
b		n				h	2	n
f (x)dx I Rл.п.i ( f ) I								f '( i )		; ξi	[xi-1,xi].
a		1		2				1

Теорема о среднем. Если f(x) C[a,b], то существует ξ [a,b]

такое, что f ( ) 1 n f (xi ), xi [a,b]. n 1

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

		h	2	n	(b a)	2	1		n	(b a)
Rсост.л.п				f '( i )					f '( i )		f '( )h , ξ [a,b]. Порядок
					2n			n		2
	2			1					1

точности (показатель степени h в формуле погрешности) формулы

равен 1.

52

f1 ... f

8.3.1.2. Составная формула правых прямоугольников.

f

					x
	Xi-1		Xi

Ii hfi

n

I h( f1 f2 ... fn ) h fi ; Rсост.п.п Rсост.л.п . i 1

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

совпадают.

8.3.1.3. Составная формула трапеций.

f

				x
xi-1
		x i

.

Ii h2 ( fi fi 1 )

I h f0

2

Rсост.тр. ( f ) h3 n

12 1

fn

n 1

f0 fn

n 1

h fi

2

2

i 1

' '

(b a)3 1

n

''

(b a)

' '

2

f

( i )

f

( i )

f

( )h

,

12n

2

n

12

1

ξ [a,b], ξi [xi-1,xi].

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

53

2.4. Составная формула центральных прямоугольников.

f

x

Xi-1 Xi-1/2 X i

Ii hfi 1/ 2

n

I h( f1/ 2 f3 / 2 ... fn 1/ 2 ) h fi 1/ 2 i 1

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

аналогичные сделанным в предыдущем параграфе, найдём

R	( f )	(b a)	f ' ' ( )h2 , ξ [a,b].
сост.ц.пр.	24

Если число элементарных отрезков n = 2m чётное, можно избежать использования половинных узлов, для чего частичный интеграл следует вычислять на отрезках [х2i-2, x2i] i=1,…m. Средняя точка этого отрезка – x2i-1, а его длина равна 2h. В этом случае составная формула центральных прямоугольников приобретает вид

Ii 2hf2i 1;i 1, m

m

I 2h( f1 f3

... f2i 1

... fn 1 ) 2h f2i 1

i 1

Так как шаг hнов.=2h,

R

( f )

(b a)

f ''( )h

2

R

( f )

(b a)

f ''( )(2h2 ) R

( f )

(b a)

f ''( )h2

сост.ц.пр.

24

нов

сост.ц.пр.

24

сост.ц .пр.

6

Порядок точности формулы равен 2. ξ [a,b].

54

8.3.1.4. Составная формула Симпсона.

	xi
Ii P2 (x)dx
	xi 1
Здесь P2(x) – интерполяционный полином 2-ой степени,
построенный на узлах xi-1, xi-1/2, xi.
P ( x) f						fi fi 1	( x x	)	fi 2 fi 1/ 2 fi 1	( x x	)2
									h2 / 2
2		i 1/ 2				h	i 1/ 2			i 1/ 2
x
i	P2 ( x)dx		h		( fi 1 4 fi 1/ 2 fi )
			6
xi 1

h=xi-xi-1.

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

n

h

h

n

I Ii

( f0 4 f1/ 2 2 f1 4 f3/ 2 2 f2 ... 2 fn 1 4 fn 1/ 2

fn )

( f0 fn 4 fi 1/ 2

i 1

6

6

i 1

h

5

n

(b a)

. Погрешность формулы Rсост.симпс. ( f )

f (4) ( i )

f (4) ( )h4 ,

2880

2880 i 1

ξ [a,b].

Если число элементарных отрезков n = 2m чётное, можно избежать использования половинных узлов аналогично тому, как это описано в разделе «формула центральных прямоугольников».

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

2 fi )

i 1n 1

			h			m	m 1
	I			( f0 fn 4 f2i 1 2 f2i ) с погрешностью
		3				i 1	i 1
R		( f )			(b a)	f (4) ( )h4 ,	ξ [a,b]. Порядок точности
сост.симпс.				180

квадратурной формулы Симпсона равен 4-ём.

8.4.Квадратурные формулы наивысшей степени точности (формулы Гаусса).

Степень точности квадратурной формулы может быть повышена путем специального способа выбора узлов. В

квадратурных формулах Гаусса весовые коэффициенты и узлы

55

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

b	b
| p(x)xi | dx , i 0,	p(x)dx 0 .(*)
a	a

Тогда параметры Ak и xk квадратурной формулы с n узлами можно выбрать так, чтобы она была точна для полиномов N=2n-1

степени.

Теорема 1. Пусть квадратурная формула имеет стандартный

вид

b	n
p(x) f (x)dx Ak f (xk ) (1).
a	k 1

Тогда, чтобы (1) имела степень точности (2n-1) необходимо и достаточно выполнение следующих условий:

1)Формула (1) должна быть интерполяционной, т. е. её весовые коэффициенты должны вычисляться по формулам

b	(x)
Ak p(x)			dx;
	(x xk )
a		'(xk )

2) Узлы интерполирования xk должны быть такими, чтобы многочлен

b

p(x) (x)Q(x)dx , ω(x)=(x-x1)…(x-xn) был ортогонален с весом p(x) ко

a

всякому многочлену Q(x) степени, меньшей n, т. е. чтобы

b

p(x) (x)Q(x)dx 0 .

a

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

Необходимость. Предположим, что (1) точна для всех полиномов степени (2n-1). Тогда она должна быть точна и для полиномов степени (n-1), а, следовательно, согласно доказанной ранее теореме о степени точности ИКФ эта формула должна быть интерполяционной.

56

Предположим, что Q(x) = Qn-1(x) любой полином степени

(n-1). Тогда ω(x)Q(x) есть полином степени (2n-1). Пусть f(x) = ω(x)Q(x) и, следовательно,

b	b	n
p(x) f (x)dx p(x)(x)Q(x)dx Ak(xk )Q(xk ) 0
a	a	k 1

b

p(x)(x)Q(x)dx 0,

a

что соответствует условию (2) теоремы. Необходимость доказана.

Достаточность. Рассмотрим любой многочлен f(x) степени

(2n-1). Разделим его на ω(x). Тогда f(x) = ω(x)Qi(x)+ri(x), где Qi(x) и ri(x) многочлены степени i (n-1). Поэтому

b	b	b	n
p(x) f (x)dx p(x) (x)Qi (x)dx p(x)ri (x)dx Ak (xk )ri (xk ) .
a	a	a	k 1

Причем, последнее равенство выполнено точно, т. к. первый интеграл равен нулю по условию теоремы, а для второго интеграла ИКФ является точной. Достаточность доказана.

Свойства многочлена ω(x) уточняются дополнительными теоремами.

Теорема 2. Пусть выполнены условия (*). Тогда существует единственный многочлен ω(x) указанного выше вида, ортогональный с весом p(x) ко всякому многочлену степени меньше n.

Теорема 3. Пусть p(x) > 0 на отрезке [a,b], а многочлен ω(x)

ортогонален с весом p(x) любому многочлену Qn-1(x). Тогда корни

ω(x) лежат внутри интервала [a,b] и различны между собой.

Заметим, что требование нахождения корней ω(x) внутри интервала [a,b] становится необходимым, если вне этого интервала f(x) не определена.

57

Возникает вопрос об отыскании многочленов ω(x),

удовлетворяющих условиям, сформулированным в теоремах.

Рассмотрим случай p(x)=1, x [a,b]. Используя замену переменной x=(a+b)/2+t(b-a)/2, задачу интегрирования преобразуют к виду

b		b a 1		a b			b a
	f (x)dx			f				t dt;t [ 1,1]..
		2			2		2
a			1

Известно, что ортогональную на интервале [-1,1] систему многочленов с постоянным весом образуют многочлены Лежандра.

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

(n+1) n+1(t)+(2n+1)t n+n n-1(t) = 0 , 0=1, 1=t.

Другой возможный способ - вычисление по формуле

n (t)	1		d n (t2 1)n
	2n n!		dtn

Полиномы Лежандра

1) ортогональны на интервале [-1,1] любому многочлену Qi(t)

степени, меньшей n:

1
Qi (t)	n (t)dt 0, i<n;
1

2) корни tk действительные, простые и принадлежат интервалу

[-1,1].

Полиномы Лежандра удовлетворяют условиям,

сформулированным в теоремах.

Таким образом, узлы КФ Гаусса могут совпадать с корнями полиномов Лежандра.

Веса Ak вычисляются по формуле

Ak	2	, k=1,...,n.
	(1 tk2 )[ n '(tk )]2

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

58

качестве ω(t) можно использовать нормированный полином Лежандра n/an, где an коэффициент при старшей степени полинома Лежандра.

Погрешность ИКФ Гаусса определяется по формуле

		b a 2n 1								(2n)
Rn	( f ) (n)							f			( ), [a,b]
				2
(n)=		22n 1			(n!)2					2
		(2n 1)(2n)! (2n!)

Напомним, что под n подразумевается число узлов ИКФ.

Величина α(n) быстро уменьшается сростом n, что видно из таблицы

N	2			3			4			5			6
α(n)		1			1			1			1			1
	135			15750			3472875			1.24 *109			649 *109

Значения весов и узлов КФ Гаусса на интервале [-1,1] не зависят от вида подынтегральной функции. В настоящее время построены таблицы узлов и весов квадратур Гаусса по крайней мере до n=4096 c 20 двадцатью десятичными знаками.

Если p(x)= (1 x2 ) 1 , то за узлы xk следует принимать корни

полинома Чебышева.

В этом случае

1

n

2k 1

1

f (x)dx

f (xk ), xk

cos

1 x

2

1

n k 1

2n

R ( f )

f (2n) ( ), [-1,1].

n

22n 1(2n)!

В приведенных формулах x [-1,1]. Произвольный отрезок может быть преобразован к стандартному соответствующей заменой переменных.

59

8.5.Практическая оценка погрешности (правило Рунге).

Пусть Ih приближенное значение интеграла I, вычисленное с

помощью какой либо составной квадратурной формулы.

Предположим, что для погрешности этой формулы справедливо

представление
I - Ih = Chk + o(hk),	(1)

где C ≠ 0 и k > 0 величины, не зависящие от h. Запись o(hk)

означает, что o(hk) пренебрежимо мало по сравнению с hk при h 0.

Тогда величина Chk называется главным членом погрешности

квадратурной формулы.
Заметим, что из (1) следует |I –Ih|	hk с некоторой

постоянной > |C|. Поэтому число k является ни чем иным как

порядком точности соответствующей квадратурной формулы.

Если подынтегральная функция достаточно гладкая, то для

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

Процедуру выделения главного члена погрешности

продемонстрируем на примере составной квадратурной формулы трапеций. Ранее было получено

h

3

n

h

2

n

h

2

b

R

( f )

f ' ' (

)

hf ' ' (

)

f ''(x)dx .

сост.тр.

12

i

12

i

n ,h 0

12

1

1

a

1

b

Таким образом, Rсост.тр. = Ch2 + o(h2), где C

f ''(x)dx .

12

a

При половинном шаге из (1) следует

I - Ih/2 = (1/2k)Chk + o((h/2)k)

(2)

60