книги / Численные методы. Ч. 1

Пусть рассматривается система линейных алгебраических уравнений Ах = f с симметричной положительно определенной матрицей А. Решение будем искать с помощью явного нестационарного метода Ричардсона,

х(п+,)- х (п)		+ Ax(n) = f. п = 0,1,...
	т (п+1)
Попытаемся так определить		набор т(1),т (2),...,x(N\	чтобы ||x(N)-x|| была
минимальной для заданного числа итераций N.
Теорема 2.6. Пусть А	симметричная положительно		определенная матрица,

Xmin > 0, Хтах >0 - наименьшее и наибольшее собственные значения. Пусть задано число

итераций N. Среди всех	наборов т(п),	n = l,N,	наименьшую погрешность ||x(N)-x||
имеет набор, для которого
т	(п ) ’ n = 1,N;
1+Pot'
		Asa!-; t'nl =cosf o		*)*.
^ min	^ m	^		2N
Оценка погрешности в этом случае имеет вид
			Pi = ь К
			i+ V T
Доказательство. Введем, как и ранее, погрешность решения				z(n) = х(п) - х . Схема
Ричардсона позволяет записать систему уравнений относительно погрешностей:
	z(D+,)- z (n) + Az(n) =0,		п = 0,N -1.
	т(пИ)

Отсюда получаем

ztn+l) = 2<n) - х ^ А г '" ' = (E - X("*"A )Z<"1

В частности,

Z<I) = ( E - T",A )Z(0>,

z,2) = (Е - T<2,A )Z(1) = (е - х,2)аХе- t ("A)z,0,1

Z(N)

. п=1

T(N) = (е - t (N)а ) •(е - X(N-"A ) •... • (е - х(,)а ) = П(Е- Т<П>А) •

п-1

Понятно, что T(N)- симметричная матрица. Теперь погрешность на N итерации можно представить выражением

Z(N) =X (N)Z(0\ ||Z(N)||=||T (N)Z(0)||.

Для симметричной положительно определенной матрицы в качестве нормы может быть выбран спектральный радиус v = Xme(T(N)). В самом деле, для собственного вектора р, соответствующего собственному значению Хшах,

T(NV = Xm„ n ,

||т Ы)и||= Цх^ цЦ= Ix^i-pl=v-|HI•

Сучетом свойств нормы получаем

v- H = |r (NV p lT (N1 -M .

			V S ||T (n)|.						(2.19)
С другой стороны,	пусть	р к,	к = 1,ш	- ортонормированная			система векторов,
построенная на основе собственных векторов матрицы T(N),
		(цк.й " )= Ё и !1ц Г = 8 к .-
			i-1
Разложим вектор р	по этому базису:
			k=l
Согласно определению нормы вектора
( m	Го	\	т f т	m	\	т f	m	N	Ь
IN2=(ц,ц)=		У	= 1		J
Ч к=1	n=l		.= lV k=i	n=1		tn = ,V	,=i	V

В компонентной записи вектор T<N)p	с использованием	собственных чисел и
векторов выглядит следующим образом:
{T‘NV}, = £ т Г и , = Е ( Ч М)	= Х л { t o ? }	= i > t Л и *

И

j=l 4

k=l

j

k=l

4 j=j

J

k=l

Здесь использовано определение собственных значений и векторов:

T(NV = X kpk.

Учитывая, что

|T 'N,H||: =(T‘NV .T,NV ) = t ( Z c kx k^	£ с пХпиг1 =			=
1=1 Vk=l	n=l	S k.n=n-	i=l	'

	= t ( v . X kX.ek.) =			k=l	=
	k.n=l	k=l
	можно подсчитать норму оператора	1Г Ч 1.		JU	M|
	= SUP
		н и *	SUP НГн		v *
	йсай*0	||p||	цеН.ц*0	щ

Сравнивая последнее неравенство с выражением (2.19), определяем точное значение нормы:

||T (N)|= V .

С учетом этого можно оценить величину погрешности: ||z(N,|= | r (N,z,0,|S v .|zW|.

Построение доказательства теоремы 2.6 основано на поиске такого набора

т(п), n = 1,N, который минимизирует спектральный радиус v матрицы T(N).
Предположим, что все собственные значения матрицы А упорядочены:
	О Xj = Xmin < Я,2 <	< ^m-1 <	= ^-тшх*
Известно [7],	что если f(A) - матричная функция матричного аргумента А, то
f(X,), f(X2),...,	- полная система собственных значений матрицы f(A).
Поскольку
	T1N)(A) = (е - X(N,A ) ■(е -	X(N‘"A )-... •(е - х(|,а )

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

(l- x ,N,Xl) ( l - x ,N- '\ ) - . . . ( l - x ,llX1)

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

v ^ a x l ( 1- x 'N'X1).(1- x ,- ,X1)....(l-x<"X 1)|.

f(X) = (l- x(N)x) (l - t <N_,,x)•... (l- T(I>X) .

( 2. 20)

Тогда спектральный радиус можно определить следующим образом:

max |f(X)|

и определение набора т(п), n = 1,N сводится к задаче поиска

min	max |f(X)| =	m inj|f|
Очевидно, что функция (2.20) является полиномом степени N, причем f(0) = 1.
Иначе говоря, поиск итерационных параметров		т(п\ n = l,N сводится к задаче об

отыскании полинома степени N, наименее уклоняющегося от нуля на отрезке [linin’А.та*]* которая может быть решена с использованием полинома Чебышёва.

Корни функции (2.20) принимают значения:

1- т,п)Х = 0, Х „ = - ^ , п =Т к

и должны совпадать с корнями полинома Чебышёва:

х „ = Хп,“ + ^ min + Хт" ~ X|nin cos^ ~						, П = Ш .
п		2	2		2N
Теперь очевидно, что итерационные параметры следует выбирать следующим
образом:
	' =Х		—^тах ^~kmm ^ k max		k m|n	О71 _
V	/	n	2		2	2N
_	k mux	k mln [ 1 . k mn ~ k m|n		(n) I _		1 . 1~~^3 t(n) I
			l	J	т < 4	1+ E ) '

n = 1,N.

1+Pot(n) '

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

Для оценки нормы погрешности заметим, что

2РГ v £ ^ al J f(X)H fll= 1+р,

( N ) | b	" z
откуда получаем z’

1+РГ

что и требовалось доказать.

Рассмотрим неявную итерационную схему с положительно определенными матрицами А и В:

	Х(А+‘) _ Y(n)			( 2.21)
	В		I Лх(п) - f
	°	-(n+l) +АХ _ I -
Эта система уравнений для погрешностей принимает вид
	s Z(n+1) - z (n)
	В-		-+Az(n) = 0.
		r (a+D
Указанные свойства матрицы		В	позволяют представить ее в виде В=В,'2В1/2
Тогда предыдущее соотношение можно представить в форме
	B,/2Z(n+l) - B,/2z(n)- + B"I/2Az(n) =0.
	T (n+I)
Обозначим B,/2z(n) = v(n), z(n) = B',/2v(n) Тогда
	v (n+l)		(n)
	~ 5 ^ +Суи,=0*
где С = B”I/2AB”i/2	симметричная		положительно определенная	матрица, причем

минимальное (максимальное) собственное значение матрицы В-1А является одовременно и минимальным (максимальным) собственным значением для матрицы С.

Теорема 2.7. Пусть		А	и		В - симметричные и положительно определенные
матрицы.	Хтах - наименьшее и наибольшее собственные значения матрицы В_,А .
Для заданного	числа	N	итераций неявный			чебышёвский метод (2.21) имеет
минимальную	погрешность		при		наборе т(а),	n = l,N, определенном условиями
предыдущей теоремы, где			^	(	в-'а )
			^	(	в-'а ) '
Удачным выбором матрицы					В можно приблизить значение параметра £ к 1, что
приведет к понижению погрешности |z(N)|.

Погрешность z(n) = х(п) - х решения системы линейных алгебраических уравнений вида1 Ах = f определить невозможно, поскольку точное решение х неизвестно. Однако можно оценить невязку

Г(П) = А2<">= А(х<">_ xj = Ах("> - Ах = Ах(п>- f ,

определяющую, насколько полученное решение не удовлетворяет исходному уравнению.

Рассмотрим явный итерационный метод:

- + Ax<n) = f .

-(п-Н)

Определим из этого соотношения величину х(п+|):

х'"*'1= х<») _ t C»')^Ax(-l _ f) = х<") _ т(-')г(п) .

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

Оценим невязку для следующего шага:

г(п-и) _ д х(п+1) _ f = Ax(n) —x(n+,)Аг(п) - f = г(п) —х(п+1)Аг(п) Определим, как и ранее, квадрат нормы невязки:

J r '- 'f =(г,п) - т '^ ’Аг'"'. г1"1-T (°*"Ar(ol) =

= (г("|,г(",) - ( т (“*,,Аг<",,г<*,) - (г 1”),т(п*1)Аг(“,)+ (т1п*,,Аг(“,)т|"*|)Аг1“’) =

= fr<-f - 2т(“*"(г("\ Аг(0|)+ (т<°*||)!}Аг(”)|!

При выводе последнего выражения учтено, что (Au, V ) = (u, ATv) = (u, A v),

AT = A - в силу симметрии матрицы.

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

1Здесь и далее предполагается, что А - симметричная, положительно определенная матрица.

66

нахождения значения итерационного параметра, при котором невязка г'"*'* минимальна, воспользуемся теоремой Ферма1:

d I r - f

= -2(rl“', Ar‘“’)+ 2х(" ”|Аг1°,|1 = (

						'	K		f
Оценка невязки получаемого решения:
				NxU,'xIsp°Kx<01~4
Здесь	г и I ,	г	£ = —m‘g--		%	Д		наименьшее и наибольшее собственные
			^	ч	mio’ max
				^гшя
значения матрицы		А; п - номер итерации.
Метод минимальных поправок
Неявную итерационную схему
						(п+1) _	(«О
					В	— *— + Axw = f
можно представить			виде
					_ln+l)		(п)
					В*			+Г1*1=0.
					^(0+1) _ ^(п) _
где, как и ранее. г(п^ = Ax*n*—f . Вектор								= В~1г ^ назовем поправкой. Очевидно, что
поправка w‘n) удовлетворяет уравнению
					•д,(П+1) _ «|(®)			+ Aww =0.
					BW	, ,
Предполагая,			что	В	симметричная				положительно определенная матрица,
определим норму в виде

M B =V(Bv' v)-

1Пьер Ферма [17.8.1601-12.1.1665] - французский математик. По профессии был юристом, с 1631 года являлся советником парламента в Тулузе. Основные научные труды изданы лишь после его смерти.

Теорема Ферма [10]: пусть функция у = Г(х), непрерывная в некотором замкнутом интервале [а, Ь], принимает свое наименьшее (наибольшее) значение во внузренней точке £ этого интервала. Если в точке

\ производная функции Г(х) существует, то она равна нулю.

Подсчитаем значение квадрата нормы поправки w(n'fl) = w11” - x^'B^Aw**’:

Jw—”||; =(B{W‘", - T,,,*1,B-'AW,",),WI“' - T("*"B-'AWi"i) =

=(Bw(n), w1"’)- Tln*',(Aw"”, W<”,) - T,"*I)(BWi“), B‘1Aw(",)+(T,"*',)J(Aw,“,,B'lAw(n,) =

=(BW1”’. WII”)-2 T,IW')(AW"", W'",) - ( T("'")!(AW("1, B^AW11”^

=||w""|fB-2 TI""||W'">|[ - ( T,,I*")J||AWI",|B, .

При получении последнего выражения использовались следующие соотношения:

(Bw|n,.B lAwi",) = ((B-'A)TBwln|,wl")) = (АВ 'Bw1"’,w|n|j = (AW1"1,w’"').

AT = A, (B ')T=B-‘,

IW<”,|IA =(AW1"1, W*"’),

||Aw,n)£., =(Aw"’ B-'AW'"1).

Очевидно, что величина ||w(”*l)flBбудет минимальной при условии

—L J s . = - 2||w(n)|A+2T,"‘I)||AW<“’1B.I = ° ,

T,„..) = K i l l =	(Aw(”>, Ww)
|AW(“’||b.,	(AW" 1, B^AW"”) ‘

Для реализации метода минимальных поправок требуется на каждой итерации решать систему уравнений Bw(n) = r(n), откуда находится сама поправка w(n) Кроме этого необходимо определить решение системы уравнений Bv(n) = Aw(n), а именно,

вычислить v(n) = B“,Aw(nJ, требующееся для нахождения итерационного параметра.

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

учетом введенного определения нормы):

lA(x<n)-xIB-^psHxl0,-xL-

Как и ранее, р0 = 7 —\\	=	XminA max	наименьшее и наибольшее
1+ S
собственные значения матрицы	В *А ;	п - номер итерации.

Относительно погрешности zin) = х(п) - х итерационная схема Ричардсона, как это было уже показано ранее, принимает вид:

z (n+l) _ z (n)

-+Аzw =0.

-(n -H )

Отсюда, погрешность

z(n +1) _ z(n) _ т<п + ,)д2(п)

Определим выражение

jjz “■“!£= (Az‘ !. z'**") = (A (Z"" - T‘"-"Az ' z ' n; - т1 Az “') =

=(Az'”'. Z(“')-(A T'"'"AZ'°’. Z'"1)-(A Z"”. X'"'"AZ'°I)+(AX'°'"AZ'“I. T'“‘"Az'”’) =

=(Azln’.z,n,)-2xl"*ll(Az(",,Az("’)+(x,"-")!(A,z'“’.Az'nl) =

=|z'n,|:4 -2x'-"|Az",f +(х'"-"):!|Аг,п,|:л.

При выводе последнего соотношения использована симметричность матрицы А

(A AZ1”1. z1"1) =(Az'"‘. AZ‘“’).

Полученное выражение может рассматриваться как квадратичная функция итерационного параметра т1П+|). Воспользуемся теоремой Ферма для нахождения значения итерационного параметра, доставляющего экстремум этому выражению,

„ |A z'"f	K '. A z " )
|A zf'£	(A ’Z1”1, Az10')
Вторая производная

d1z,n" t =jK t >o

d(x'“-"):

положительна в силу положительной определенности А. то есть выражение

|z‘n-|,|’ - (Az'"*1', Z '"*11) принимает наименьшее значение при найденном х (п+,)

Вспоминая, что Az(n> = г(п) - невязка решения системы уравнений, получаем

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

||х‘"’ -х||А £ р п||х,0,-х||д , р0= Ь 1 , 4

где Xmi ,A.max - наименьшее и наибольшее собственные значения матрицы А; п - номер итерации.

Неявный метод скорейшего спуска

Рассмотрим неявную итерационную схему вида

х(п+1) _ х (п)

В+ Ax(n) = f

х(п+1)

ссимметричной и положительно определенной матрицей В.

Для погрешности z(n) = х(п) - х эта схема принимает следующую форму:

- + Aztn) = 0 .

т (п+1)

Отсюда, z(n+,) = z(n) - T(n+1)B",Az(n).

Как и ранее, с учетом симметрии матрицы А, определим выражение

||zln*"|£ =(Az," " .z ,n*,,) =

= (Az",l.z lnl)-2 T 1"*l,(Az(n,,B 'lAz(n,) + (x,ntll)2(AB''Az(n,, B lAzl"l) =

= l|z,nt - 2 T'"*,>||AZ-”'£ .1+(X‘- ) 2||B"AZ>"'1;.

Благодаря положительной определенности матрицы А, (Az(n+I). z<n+n) >0,

минимум полученного выражения достигается при значении итерационного параметра

Цв-'A z'i; > v t '