la1

4.3.3Свойства ортогональных полиномов

Путь задана система ортогональных с весом ½ полиномов						Pn(x) . Справедлива
Теорема. Все корни		Pn(x) вещественные, простые и принадлежат отрезку						(a; b) .
Доказательство. Пусть Pn(x) имеет k вещественных корней xi							на отрезке (a; b)		нечетной кратности.
Положим		qk(x) = 8 1k;
						k = 0	;
		<
		> jQ			xj);
		>	(x	¡		k > 0
		=1
		:
где корни xi	2 (a; b)	взятые без учета кратности, т.е. входят в произведение только один раз. Тогда
произведение	Pn(x)qk(x) не меняет знак на промежутке (a; b) , и, следовательно,
		Zab Pn(x)qk(x)½(x)dx 6= 0 :
Однако при	k < n интеграл должен равняться		0		в силу ортогональности			Pn	полиномам меньшей
степени. Таким образом k = n.

Теорема. Если алгебраическая степень точности квадратурной формулы c L узлами xk равна 2L¡1 ,

то узлы

xk

суть корни ортогонального полинома PL(x) .

L

Доказательство. Пусть NL(x) = i=1(x ¡ xi) , где

xi

узлы квадратурной формулы и пусть её алгеб-

L 1

. Рассмотрим функцию

f(x) =

L

(x)P

m

(x)

, где

m L 1

,

раическая степень точности равна 2Q

¡

N

·

¡

являющуюся полиномом степени не превосходящей 2L ¡ 1. Для такой функции квадратурная формула

точна по условию, и, следовательно,

L f(xk)¸k = 0 ;

Zb f(x)½(x)dx =

X

a

k=1

b

NL(x)Pm(x)dx = 0 и значит NL ? Pm . Таким образом NL является ортогональным полиномом

то есть

a

и в силуRединственности с точностью до множителя совпадает с

PL .

Пусть теперь корни xi ортогонального полинома PL(x)

являются узлами квадратурной формулы.

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

2L ¡ 1 . Проап-

пpоксимиpуем функцию f(x)

полиномом

g(x) степени

L ¡ 1 по ее значениям в точках

xi :

g(x) = L f(x ) (x) ;

(x) =

L

(x ¡ xj)

:

X

i Li

Li

Y

i=1

j=i

(xi ¡ xj)

6

b

b

Пусть

I = a f(x)½(x)dx ;

J = a g(x)½(x)dx : Тогда

I = J , если f

полином степени до

L ¡ 1

поскольку в этом случае

f = g . Но если f

| полином степени до

2L

1 , то разность

включительно, R

R

¡

полиномов f

и g

также полином степени не превосходящей

2L ¡ 1 , причем

(f ¡ g) jx=xj = 0 , и,

следовательно, справедливо представление f ¡ g = PLqL¡1, где

qL¡1 некоторый полином степени до

L ¡ 1 . Тогда

I ¡ J = Zab ½(x)[f(x) ¡ g(x)]dx = Zab ½(x)PL(x)qL¡1(x)dx = 0 ;

41

то есть алгебраическая степень точности квадратурной формулы равна 2L ¡ 1 , если её узлы корни

ортогонального полинома. Веса при этом равны

¸

=

Za

b

(x)½(x)dx =

b L

(x ¡ xi)

½(x)dx :

Lk

k

Za i=k (xk ¡ xi)

6

Y

Отметим, что корни соседних ортогональных полиномов PL и PL¡1 различны (на самом деле между

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

xi

и

xi+1

полинома

PL

лежит ровно один корень x~i полинома

PL¡1 ). Действительно, пусть fk(x) =

PL(x)PL¡1

(x)

, тогда degfk = 2L ¡ 1 и формула Гаусса-Кристофеля

(x¡xk)

с узлами xi , являющимися корнями полинома

PL , для такой функции точна. Она, как легко увидеть

принимает вид

b

¸kPL0 (xk)PL¡1(xk) = Za

PL(x)PL¡1(x)

½(x)dx :

(x ¡ xk)

Пусть ai старшие коэффициенты ортогональных полиномов

Pi

, тогда

PL(x) = aL Y(x ¡ xi) ; PL¡1(x) = aL¡1 Y(x ¡ x~i) ;

и справедливо представление

PL(x)

=

aL

P

(x) + q

L¡2

(x) ;

x ¡ xk

aL¡1

L¡1

где qL¡2(x) некоторый полином степени не выше

L ¡ 2 . Таким образом с учетом ортогональности

b

aLjjPL¡1jjL2 2;½

1

Za

aL

2

¸k =

PL¡1

(x)½(x)dx =

;

PL0 (xk)PL¡1(xk)

aL¡1

aL¡1PL0 (xk)PL¡1(xk)

но так как ¸k =6 1 (для весов уже получено явное выражение (7), да и кроме того, зная узлы, веса можно однозначно определить через определитель Вандермонда), то PL¡1(xk) =6 0 , и значит ни один из корней полинома PL не может являться корнем полинома PL¡1 . Попутно мы нашли и другое выражение для весов ¸k .

Свойства весов

1) ¸k > 0.

NL(x)

i

2

Доказательство. Пусть fk(x) = h x¡xk

: Это полином степени 2L¡2 , равный 0 во всех узлах, кроме

x = xk , для него формула Гаусса-Кpистофеля точна, поэтому

b

2

2

NL(x)

½(x)dx = ¸

NL(x)

= ¸

[

0

(x )]2

> 0 ;

¸

x=xk ¸

Z ·x

¡

xk

k

·x

¡

xk

k

NL

k

a

¯

следовательно ¸k > 0.

¯

b

2) связь весов ¸k с моментами

cl =

xl½(x)dx :

R

a

XL

xlk¸k = cl ; l = 0 ; 1 ; : : : ; 2L ¡ 1 :

k=1

Свойство становится очевидным, если сосчитать интеграл с весом от степени xl по формуле Гаусса-

Кpистофеля.

3) PL ¸k = Rb ½(x)dx :

k=1 a

Это частный случай свойства 2) пpи l = 0 .

42

4.3.4Примеры ортогональных полиномов

1) Полиномы Лежандра Pn(x) являются ортогональными на промежутке (-1,1) с весом ½(x) = 1. С точностью до нормировки для них справедливо выражение

Pn(x) = (¡n1)n dnn (1 ¡ x2)n :

2 n! dx

Вчастности P0 = 1 ; P1 = x ; P2 = 12 (3x2 ¡ 1) :

2)Полиномы Чебышева первого рода

Tn =

n

[n=2]

(¡1)m(n ¡ m ¡ 1)!

(2x)n¡2m

X

m!(n

¡

2m)!

2 m=0

1

ортогональны на том же промежутке [¡1; 1] ;

с весом

½ =

p

.

1¡x2

с весом ½(x) = e¡x2 . С точностью

3) Полиномы Эрмита

Hn

ортогональны на промежутке (¡1; 1) ;

до нормировки они имеют вид

dn

2

2

Hn(x) = (¡1)nex

e¡x

:

dxn

4) Полиномы Лагеpра

Ln

ортогональны на промежутке

[0; 1) ;

с весом ½(x) = e¡x. Их можно

представить в виде

1

dn

Ln(x) =

ex

(xne¡x) :

n!

dxn

4.3.5Погрешность квадратурных формул

Пусть функция f(x) проинтерполирована по её значениям

f(xi) в

L точках xi

; i = 1; 2; : : : ; L ;

полиномом gL¡1 :

L

L

(x ¡ xj)

f(x) = g

(x) + r(x) ; g

(x) =

f(x )

:

L¡1

L¡1

i=1

j

j=i (xi ¡ xj)

X

6

f(x)

Y

Погрешность интегрирования R

при замене

интерполяционным полиномом

gL¡1 (она же

погрешность соответствующей квадратурной формулы) имеет вид

b

b

b

b

R = Za

f½dx ¡ Za

gL¡1½dx = Za

r(x)½dx = Za

f(L)(»(x))

N(x)½(x)dx ;

L!

L

(x ¡ xi)

N(x) =

=1

iY

и если f полином степени не выше L ¡ 1 , то f(L) ´ 0

и, следовательно, квадратурная формула точна.

Для случая pавноотстоящих узлов

xi ¡ xi¡1 = h

имеем:

jNL(x)j

·

hL max

1

·

hL

L!

k

·CLk ¸

и, значит

Zb

jRj · hLjjf(L)jjC ½(x)dx ;

a

и пpи ½ = 1

jRj · hLjjf(L)jj(b ¡ a) :

Это довольно грубая оценка, однако она показывает порядок точности по h .

43

В случае, если узлы не произвольные, а корни ортогонального полинома PL , то квадратурная формула точна для полиномов степени не превосходящей 2L ¡ 1, хотя полученная оценка этого и не "чувствует". Чтобы улучшить оценку в этом случае поступим следующим образом. Пусть f 2 C2L. Разложим ее в ряд Тейлора в окрестности некоторой точки x¤:

f(x) =

2L¡1

f(k)(x¤)(x ¡ x¤)k

+ f(2L)(x¤)(x ¡ x¤)2L + q(x) ;

X

k!

|

(2L)!

}

k=0

f1(x)

f2{z(x)

тогда

|

{z

}

Zb

Zb

Zb

R = [f ¡ gL¡1(x)]½(x)dx = [f1 ¡ gL¡1(x)]½(x)dx +

		a			a
	½ = 1						=0				(x)
					отбросив от функции f
Пусть вес		, оценим последний интеграл			|		{z			2}
точки разложения x точку			a+b	, тогда
			2
		¤

f2(x)½(x)dx :

a

остаток q(x) и выбрав в качестве

Za

b

f(2L)(x¤)(x ¡ x¤)2L

dx

jjf(2L)jjC

(b

a)2L+1

;

· 22L(2L + 1)!

¡

(2L)!

то есть погрешность R ведет себя как

R

»

jjf(2L)jjC

(b

¡

a)2L+1

:

22L(2L + 1)!

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

В этом пункте мы будем считать, что вес ½ = 1 ;	и, что L число узлов на [a; b] .
4.4.1 Число узлов L = 1	R b
	b
a) Формула левых прямоугольников: x1 = a;	f(x)dx ¼ (b ¡ a)f(a) :
	a
б) Формула правых прямоугольников: x1 = b;	f(x)dx ¼ (b ¡ a)f(b) :
в) Формула средних (прямоугольников)	a
	R

формула наивысшей алгебраической степени точности (она должна быть точной для полиномов не превосходящих степени 2L ¡ 1 = 1). Построим ее в соответствии с изложенными выше соображениями для формул Гаусса-Кристофеля. Для этого сначала с помощью масштабного преобразования и сдвига переве-

дем отрезок

[a; b]

в отрезок

[¡1; 1]

, на котором ортогональными являются полиномы Лежандра Pk .

Тогда

Z

2

Z

2

2

2

2

b

b ¡ a

1

b + a

b ¡ a

b + a

b ¡ a

f(x)dx =

f(

+

y) dy ; x =

+

y :

a

¡1

P

q(y)

y = 0. Вес ¸ (по свойству весов

(y) = y

этого полинома точка

Поскольку

1

, то единственный |корень {z

}

P¸i = Rb ½(x)dx) равен ¸ = R1

dx = 2, таким образом

a¡1

Z	2		Z	¼ 2			¡	2
	b		1						b + a
	f(x)dx =	b ¡ a	q(y)dy		b ¡ a	2q(0) = (b		a)f(		) :

a¡1

44

4.4.2Число узлов L = 2

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

Здесь узлами являются точки x1 = a ; x2 = b. f(x) заменяется интерполяционным полиномом первой степени p1(x), построенным по этим узлам:

f(x)	!	g	(x) =	x ¡ b			f(a) +	x ¡ a			f(b) ;
		1		a	¡	b		b	¡	a

Z		¼	Z	a ¡ b
	b			b
	f(x)dx f(a)			x ¡ b	dx + f(b)

aa

Zb

x ¡ adx = b ¡ a

a

a ¡ b

Za

b ¡ a

Za

¡

(a ¡ b)

2

(b ¡ a) 2

f(a)

b

f(b)

b

f(a)

(a ¡ b)2

f(b)

(b ¡ a)2

=

y dy +

y dy =

+

=

= (b ¡ a)

f(a) + f(b)

:

2

Эта формула разумеется точна для полиномов степени не превосходящей L ¡ 1 = 1 (и не больше).

б) Формула Гаусса-Кpистофеля для L = 2

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

Z	2		Z	2		2
	b		1		b + a
	f(x)dx =	b ¡ a		q(y)dy ; q(y) = f(		+	b ¡ a	y) :

a¡1

1

2

1

1

Полином P2 имеет вид: P2 = 2

(3y

¡ 1) . Его корни y1 = ¡p

; y2

= p

. Веса из симметричности должны

3

3

1

1

1

быть одинаковы: ¸1 = ¸2 ; ¸1 + ¸2 = 2 ) ¸i = 1, следовательно,

g(y)dy

q(

p

) + q(p

). Таким образом

3

3

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

2 ¡ 2 p3

¡R1

¼

¡

Za

¼ 2

·

2

2 p3 ¸

b

b ¡ a

b + a

b ¡ a

1

b + a

b ¡ a

1

f(x)dx

f(

) + f(

+

) :

Алгебраическая степень точности M равна

2L ¡ 1 = 3.

4.4.3Число узлов L = 3

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

Здесь узлами являются точки x1 = a ; x2 = a+2 b ; x3 = b. Для удобства вычислений перейдем к отрезку

[¡1; 1] масштабным преобразованием q(y) = f(b+2a + b¡2a y) :

Z	2		Z	1 ¡ 2	3
	b		1
	f(x)dx =	b ¡ a		q(y)dy ; y = 1 ; y = 0 ; y = 1 :

a¡1

Заменим q(y) интерполяционным полиномом p2(y) :

q(y) ! p2(y) = q(¡1)L1(y) + q(0)L2(y) + q(1)L3(y) ;

где

L1

(y) =

(y ¡

0)(y ¡ 1)

=

y(y ¡ 1)

;

L2

(y) =

(y

¡ (¡1))(y ¡ 1)

=

(y + 1)(y ¡ 1)

;

2

(0

( 1

¡

0)(

¡

1

¡

1)

( 1))(0

¡

1)

¡

1

¡

¡ ¡

45

L3

(y) =

(y ¡ (¡1))(y ¡ 0)

=

(y + 1)y

:

(1

¡

(

¡

1))(1

¡

0)

2

1

Тогда интеграл ¡R1 p2(y)dy равен

Z

¡1

Z

¡

Z

2

2

1

y(y ¡ 1)

1

(y + 1)(y ¡ 1)

1

y(y + 1)

q(

1)

dy + q(0)

dy + q(1)

dy :

¡1

¡1

¡1

Сосчитаем веса

¸3

= ¸1

=

Z

(y2 ¡ 2 )dy = 3 ; ¸2

=

Z (1 ¡ y2)dy =

3 ;

1

2

y

1

1

4

¡1

¡1

1

q(¡1)

4 q(0) + q(1)

таким образом

q(y)dy

+

: Возвращаясь к исходной функции f, получаем формулу

Симпсона

¡R1

¼

3

b

3

3

Za

f(x)dx

¼

b ¡ a

[f(a) + 4f(

a + b

) + f(b)] :

6

2

Заметим, что эта формула точна и для полиномов третьей степени, хотя построение гарантировало точность лишь до значения L ¡ 1 = 2.

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

4.5Составные квадратурные формулы

Разобьем промежуток интегрирования [a; b] на N частей x0 = a ; x1 ; : : : ; xN = b и на каждом промежут-

ке ¢i = [xi; xi¡1] применим ту или иную квадратурную формулу и просуммируем по всем промежуткам.
Пусть hi =Nxi ¡ xi¡1 . Получаем следующие составные квадратурные формулы
M =		hif(xi+xi¡1 );
	i=1				2
	N	f(x		)+f(x			)
T =	P					i¡1		;
	hi		i		2
i=1
S =	N		[f		+ 4f(xi¡1+xi ) + f(x						)].
	P hi
	=1	6	i¡1					2		i
ЛюбопытноiPотметить, что S = 32 M + 31 T .
Удобно составную формулу Симпсона представлять в виде (при четном числе промежутков)
								¹	h
								S(f) =	3	(f0 + 4f1 + 2f2 + : : : + 4fN¡1 + fN ) :

Такая запись называется обобщенной формулой Симпсона.

4.5.1Сходимость квадратурных формул

Устремим в составных квадратурных формулах ранг дробления h = max hi к нулю. Естественным образом возникают вопросы

1)Стpемится ли сумма к интегралу?

2)Если "да", то с какой скоростью?

46

Ответ на первый вопрос положителен. Поскольку и формула средних M и формула трапеций T суть интегральные суммы, а для интегрируемой функции интеграл по определению есть предел интегральных сумм. Поскольку формула Симпсона S является линейной комбинацией (с суммой коэффициентов равной 1) формул средних и трапеций, то при ранге дробления стремящимся к нулю, она также стремится к

интегралу. Нетрудно доказать сходимость и других квадратурных формул.

Теперь обратимся к вопросу о скорости сходимости. Поскольку формулы трапеций T и средних M

точны для полиномов степени не превосходящей 1 , то естественно ожидать, что их погрешность есть

O(h2), а для формулы Симпсона, имеющей алгебраическую степень точности равную трем, погрешность

O(h4).
Рассмотрим ситуацию детально. Пусть x¹i = xi+xi¡1						: Разложим f(x) в ряд Тейлора в окрестности точки
				2
x¹ .								1
		f(x) = f(¹xi) + (x ¡ x¹i)f0(¹xi) +							(x ¡ x¹i)2f00(¹xi)+
								2
+	(x ¡ x¹i)3		f000(¹x ) +	(x ¡ x¹i)4	f(4)(¹x		) +	(x ¡ x¹i)5		f(5)(¹x	) + O(h6) :
	3!		i	24		i			120	i	i

Проинтегрируем это разложение по промежутку [xi¡1; xi]. Заметим, что при этом все члены Тейлоровского разложения с нечетными степенями (x ¡ x¹i) пропадут из-за симметрии расположения точки x¹i. Таким образом

xi		h3			h5			h7
xZ
	f(x)dx = hif(¹xi) +	i		f00(¹xi) +	i		f(4)(¹xi) +	i		f(6)(¹xi) + : : : :	(8)
		3!2	2		5!2	4		7!2	6
i¡1

Из тейлоровского разложения также нетрудно видеть, что

f(x ) + f(x

)

h2

h4

i

i¡1

= f(¹x ) +

i

f00(¹x ) +

i

f(4)(¹x ) + O(h6) ;

2

i

2!22

i

4!24

i

i

откуда

f(x ) + f(x

)

h2

h4

i¡1

f(¹x

) =

i

¡

i

f00

(¹x )

¡

i

f(4)(¹x )

¡

O(h6) :

2

2!22

4!24

i

i

i

i

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

¡ 12 f00(¹xi) +

4!24 f(4)(¹xi)[5 ¡ 1] + O(hi ) :

xZ

f(x)dx = hi

2

¡

xi

f(xi) + f(xi 1)

hi3

hi5

1

6

i¡1

Далее, поскольку S =

2 M +

1 T , то

3

3

f(x)dx = 3i ·2f(¹xi) +

2f(xi¡1

¸+

xZ

i

xi

h

(f(x ) +

))

i¡1

f(4)(¹xi)h5

2

1

1

4

f(4)(¹xi)h5

+

i

·

¢

¡

¢

¸ + O(hi6) = S +

i

[¡2] + O(hi6) :

4!24

3

5

3

5

4!243 ¢ 5

Итого, для pавноотстоящих узлов из (8) погрешность составной формулы средних ±M равна

b

b

N

N

h3

a

a

X

X

±M = Z

f(x)dx ¡ M = Z f(x)fdx ¡ i=1 hif(¹xi) = ¡ i=1

i

f00(¹xi) + O(h5) ;

24

то есть

1

N

N

(b ¡ a)

±

h3 f00

=

h2jjf00jjC

h =

f00

h2 :

X

j

M j ·

24

Xi

24

i

24

jj

jjC

=1

i jj jjC

i=1

(9)

(10)

47

Из (9), аналогично

(b ¡ a)

±

T j ·

f

00

jjC

h2

:

12

j

jj

Из (10)

±

Sj »

jjf(4)jjCh4

(b

¡

a) :

j

6!4

Здесь имеется в виду составная формула Симпсона S . Для обобщенной формулы Симпсона надо h заменить на 2h:

± ¹	jjf(4)jjC	24h4	(b	¡	a) =		M4	h4	(b	¡	a) :
	6! 22
j Sj »						180

4.6Другие формулы

4.6.1Сплайн-квадратура

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

Пусть x 2 ¢i = [xi¡1; xi] ; hi = xi ¡xi¡1 ; ! = 1¡!¹ =

x¡xi¡1

1

для приближенного

hi

. Применим сплайн S3

интегрирования. Заметим, что на промежутке ¢i его можно представить в виде:

1

hi2

3

3

S3 (x) = !fi + !f¹ i¡1 +

[(!

¡ !)Mi

+ (¹!

¡ !¹)Mi¡1] :

6

Здесь Mi = S00(xi). Пусть S(!) = S31(x), тогда

xZxi

S31(x)dx = hi Z01 S(!)d! ; (dx = hid!) :

R

i¡1

R

1

1

При этом !d! = 21 ;

(!3 ¡ !)d! = ¡41 : Таким образом

0

0

xi

fi + fi

Mi + Mi

xZ

S31(x)dx = hi

¡

1

¡ hi3

¡

1

:

2

24

i¡1

Последний член в этой формуле "имитирует"поправку к формуле трапеций. Действительно, вторая производная от сплайна, аппроксимирует вторую производную от функции и

h3	(M	+ M	i¡1	)		h3
i	i				¼	i	f00(¹x ) ;
		24				12
							i

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

	b	N									¢
a
		X¡		fi + fi		1		Mi + Mi		1
Z		f(x)dx ¼ i=1	hi	2	¡		¡ hi3	24	¡		:

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

48

ei!x¹k

4.6.2Формулы Филона

Пусть I = Rb f(x)ei!xdx ; j!j >> 1=jb ¡ aj , а f(x) медленно меняющаяся относительно периода T = 2¼=!

a

колебаний, функция. В этом случае подынтегральная функция f(x)ei!x имеет много осцилляций на про-

межутке (a; b) и использование обычных квадратурных формул весьма затруднено, поскольку приходится делить промежуток интегрирования на большое количество частей. Однако нет необходимости заменять всю подынтегральную функцию интерполяционным полиномом. Достаточно эту процедуру проделать лишь с функцией f(x). Итак, заменим f интерполяционным полиномом p. Тогда

N

N

(x ¡ xk)

f(x)

¼

p(x) =

j=0 Lj

(x)f(x

) ;

Lj

(x) =

;

j

j;k=0

(xj

¡

xk)

X

k6=Y

J = Z

b

N

b

N

p(x)ei!xdx = j=0 f(xj) Z

ei!xLj(x)dx = j=0 Aj(!)f(xj) :

a

X

a

X

b

Интегралы Aj(!) =

a ei!xLj(x)dx берутся в элементарных функциях. Получаемые при этом формулы

приближенного интегрирования называются формулами Филона:

R

I = Zb f(x)ei!xdx ¼

N

Aj(!)f(xj) :

X

a

j=0

x0 = 1 ; x1 = 0 ; x2 = 1.

¡R1

¡R1

1

sin !xf(x)dx ,

1

Задача: Для интегралов

cos !xf(x)dx получить формулу Филона с тремя узлами:

4.6.3Составные формулы Филона

Разобьем промежуток [a; b] на N частей a = x0 < x1 < : : : < xN = b и на каждом промежутке [xk¡1; xk]

заменим f(x) интерполяционным полиномом pk некоторой степени, тогда

b

N

xk

N

xk

I = Z

Z

f(x)ei!xdx = k=1

f(x)ei!xdx » J = k=1 Z

pk(x)ei!xdx :

a

Xxk¡1

Xxk¡1

Пpимеp. Аналог формулы средних.

x Zxk

f(x)ei!xdx »x Zxk

f(¹x)ei!xdx =

k¡1

k¡1

= f(¹xk)

ei!xk ¡ ei!xk¡1

=

2

f(¹xk)ei!x¹k sin

!hk

:

2

i!

!

Оценим погрешность этой формулы. Представим f(x) приближенно: f(x) ¼ f(¹xk) + f0(¹xk)(x ¡ x¹k) . Тогда погрешность R приближенно описывается выражением

xN

N

xk

R = Z

Z

r(x)ei!xdx ¼ k=1 f0

(¹xk)

(x ¡ x¹k)ei!xdx =

x0

(¹xk)

X

¡

xk¡1

= !2 k=1 f0

µsin 2k

2k cos 2k ¶

;

2i

N

!h

!h

!h

X

т.е. если произведение !hk порядка 1, то формула Филона имеет небольшую погрешность, в противном случае погрешность того же порядка, что и значение интеграла.

49

Глава 5

Поиск минимума

5.1Случай одной переменной

5.1.1Метод золотого сечения

Пусть ©(x) : [a; b] ! R и известно, что на промежутке [a; b] функция © имеет хотя бы один локальный минимум. Для применения излагаемого ниже метода золотого сечения, от функции ©(x) не требуется даже непрерывность, достаточно лишь кусочной непрерывности. Будем пока считать, что © имеет на промежутке лишь один локальный минимум (он же и глобальный).

Метод основан на сравнении значений функции в различных точках, с последующим отбрасыванием промежутков, на которых минимум уж точно не может находиться. Ясно, что чтобы осуществлять подобную процедуру, необходимо знать значения функции, вообще говоря, в 4-х точках. Действительно, пусть a = x0 < x1 < x2 < x3 = b , и пусть, скажем, в точке x2 значение функции наименьшее из этих четырех величин. Тогда минимум © заведомо не может находиться на промежутке [x0; x1] и поэтому этот промежуток можно отбросить. Теперь на оставшемся промежутке [x1; x3] нам известны крайние значения функции и значение в одной внутренней точке. Добавляя новую точку x4 мы можем повторить сравнение значений © и вновь сузить допустимый промежуток. Как наиболее разумно размещать добавляемые точки? Представляется естественным, чтобы деление отрезков происходило подобно предыдущему делению.

x0

x1

x2

x4

x3

Это означает, в частности, что внутренние точки должны располагаться симметрично, то есть jx1 ¡ x0j = jx3 ¡ x2j = h . Если длина исходного промежутка равна l , то должно выполняться соотношение

» =

h

=

x1 ¡ x0

=

x2 ¡ x1

=

l ¡ 2h

;

l

x3 ¡ x0

x3 ¡ x1

l ¡ h

откуда

» =

h

=

1 ¡ 2hl

=

1 ¡ 2»

:

l

1 ¡ hl

1 ¡ »

50