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

Г1	с, 2	с,э	с.«
0	1	с2з	С2т
и = 0	0	1	С3т
0	0	0	1

“верхняя” треугольная матрица1, у которой равны нулю все элементы, расположенные под главной диагональю. Процедура получения такой матрицы носит название ‘"прямого хода” метода Гаусса. Очевидным условием для успешною

выполнения “прямого хода” является	а**"0 * 0, j = l,m .
“Обратный ход” метода позволяет определить искомые величины:
’*«=У«-
^т-1 —Ут-1	"~^т-1т*т*
, ^т-2 —Ут-1 ” ^т-2т^т —^т-2т-1*т-1»
т
= у . - Е с1Л -
к*2

Таким образом, “прямой ход” метода Гаусса можно	трактовать как
преобразование системы уравнений вида Ах = f в эквивалентную	систему Ux = у,

причем
V - - L		V - f "‘ -	f» ~ a»y.	’
У," а „ ’		У,'а « Г	ай
Последнюю систему соотношений с учетом вышеприведенных выкладок можно
представить в иной форме
fi = а ,|У ,+ 0 у ,			+ 0 у , + ... + 0 у я ,
	= а2|	У, + а<22 *У2 +0-У, + ...+ 0 у т ,
f,	= а „	у, +а','’	у, + а'Д>	у, + ... + 0 ут .
L	=а_	. • У! + ат 2 'У2 + a 'i) -у, + ... + а'";" • у„

1Обозначение матрицы U принято по первой букве английского слова upper - “верхний”.

и записать в виде Ly = f ,	где	L - нижняя треугольная матрица1с отличными от нуля
коэффициентами а*/"0 * 0,	j = 1,ш на главной диагонали.
Все вышесказанное позволяет			трактовать метод Гаусса как последовательное
решение двух систем уравнений:		Ly	= f и Ux = у .

Объединяя эти соотношения, получаем

LUx = f.

Сравнивая последнюю формулу с записью уравнения (2.1), можно сделать заключение, что процедура метода Гаусса эквивалентна разложению исходной

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

L	нижняя треугольная матрица с ненулевыми коэффициентами на главной
диагонали;
U - верхняя треугольная матрица с единицами на главной диагонали.
Обозначим угловые миноры матрицы А символами Aj}				j = l,m :
A,	—a 11»
A:	а,,	а .2
	321	а22
			а п
Аз		а 22	а 2Э
		3 32	а зз
Am= det(A).
Теорема		2.1.	Пусть все угловые миноры матрицы	А отличны от нуля.
А} * 0,	j = 1,ш . Тогда матрицу А можно представить единственным образом в виде
			А = LU,	(2.3)
где L	нижняя треугольная матрица с ненулевыми диагональными элементами, U

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

1Обозначение матрицы L принято по первой букве английского слова lower - “нижний”.

Рассмотрим разложение матрицы для простейшего случая т=2. Пусть в этом случае матрицы имеют вид

А,	’ а п	а 1 2 _		X ,	0 “		1	U . 2 _
			,	Ь 2 =		,	и 2 =
	_ а 2 .	а 2 2 .		^ ■ 2 1	^ 2 2		0	1

Предполагая разложимость для т=2, получаем систему уравнений относительно коэффициентов матриц L2 и U2 :

А, = L, • U2 =	О	м,а]	Г*..		*-1|“.2
	^21 ^22_ о ij		[ к ,	^2lU12+*-22_
	i 12 = а12; ,Х21- а2);			<< +	СО II
			>-2.“l2

Решение этой системы уравнений дает значения коэффициентов:

Хп = а,, = А, ^ 0; Х21= а 21;

.	а1212	_ а“ паП “ 22„ - ав 2 1 а 12 _ А 2 ^ Q
^22 “ а22	а 21 а	“

Это означает, что возможность разложения (2.3), соответствующего условию теоремы, показана для простейшего случая т=2. Теперь предположим, что условия теоремы выполнены для случая (т-1), то есть имеет место разложение Am, = Lm_,Um,. Вводя обозначения

* (a m-l) (а т1 а го2 а т 3

a mn»-l)'

Ulu

U2n

* (^*т-|) (^“ml ^т2 ^inl

представим матрицы Am, Lm, Um в форме

Ащ-I {^m-l}			L n.	{0}	U n = <°)	1 J'
(®ml) ®mm	•	L	= {^ml)	^mm .
Здесь принято, что
	О,
	0:	,	<0) - (о,	0, 0,	0„_,).
м =	0,

	.0-1
Покажем существование разложения для матрицы			Ат . Перемножая		Lm• Umи
сравнивая результат	со структурой матрицы	А т,	получаем	систему	(2m-1)
алгебраических уравнений относительно (2 т -1) неизвестного {um_,},				(Xm_,),	Xmm
^m-l ‘ {^m-1 }
(^m-l	~^mm *
В силу того, что	det(Lm_,) * 0, det(Um_,) * 0	(треугольные матрицы с ненулевыми

значениями на главных диагоналях), существуют единственные решения первых двух

систем уравнений: {um_,} = L^., • {am_,},	(Xm_,) = (a ^ ^ U ;1., . И, кроме того, может
быть вычислено последнее неизвестное	A.mm = атш -(Х т_,)*{ит_,}, значение которого

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

det(A) = d et(A J = det(LmU m) = d et(L J d et(U J = Хшт d e t ^ J					d e tO J ^ ).
В силу того, что det(A)*0,	det(Lm_,) *0,	det(Um4) = l,		имеет	место	равенство
det(A) =Xmmdet(Lm,), из которого	сразу следует,		что Xmm* 0 .		При выполнении
последних преобразований вновь учтено,		что	матрицы	Lm, Um		являются

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

Покажем единственность разложения (2.3). Предположим, что такое разложение возможно не единственным образом, то есть

А = LmU(1) = L(2)U(2).

Проводя простейшие преобразования, получаем

Цп _ L(2)U(2)U(1) «

I -• I	= IT TJ-1
*-(2)^(1)	и (2)и (1) •

Можно заметить, что поскольку исходные и обратные матрицы сохраняют свою “треугольную” форму, в левой части последнего выражения расположена нижняя

треугольная матрица, а в правой части - верхняя треугольная. Но равенство в этом

случае возможно лишь	тогда,		когда		матрицы	Ц^Ь(1)	и U(2)U(",,)	являются
диагональными. Более того, поскольку				U(2) и		содержат единицы на главной
диагонали, полученные	диагональные				матрицы	будут	единичными,	то есть
L;5,L(1) =U (2)U(",,) = Е. А отсюда следует,				что L(l) = Ц 2), U(I)=U (2), но это				и означает
единственность разложения (2.3).
Рассмотрим матрицу коэффициентов системы линейных алгебраических
уравнений следующего вида:
	1	1	О
	А = 1	1	1		det(A) = - l * 0 .
	О	1	1

Первый шаг метода Гаусса приводит матрицу к виду

1 1 О"

А = О 0 1

О1 1

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

угловых миноров = 0, то есть имеется противоречие условиям теоремы 2.1. Для

успешного выполнения процедуры метода Гаусса следует поменять местами строки матрицы А (смена столбцов приведет к необходимости изменить порядок неизвестных), например:

1 1 О'

А* = 0 1 1

1 1 1

Приведение исходной матрицы к новому виду можно интерпретировать как некоторое преобразование вида А* = Р-А. Для рассматриваемого случая

1 0 0

Р = 0 0 1

0 1 0

Теорема 2.2. Если определитель матрицы коэффициентов det(A)*0, то существует матрица перестановок Р такая, что у матрицы А* = Р А все угловые миноры отличны от нуля.

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

Рассмотрим простейший случай - систему двух уравнений:

Если а,, *0, теорема справедлива при Р = Е , то есть в случае тождественного

преобразования.
Пусть	теперь а ,,= 0 . Поскольку	det(A) = а,,а22- а 12а21 = - а |2а21 * 0 , то						а2| *0
Преобразуем исходную матрицу к виду
	а ; = р а 2 =	р=	0	1
			1	о	*
Теперь	А, = а 2 |*0, А2= а21а12- а иа22 * 0 , то		есть		при	ш	2	теорема
справедлива.
Пусть утверждение теоремы также верно и для			Ат_,. Рассмотрим				матрицу Ат
(обозначения введены ранее):

Ат-|

А,

(^т-1) **тт

1.Рассмотрим случай det(Am_,) * 0 . Сформируем матрицу преобразования

Г?.-, («У

"(о)

м построим новую матрицу

	A l = РтА_	^т-1Ащ-i Рт-|{ат-|}
		(а т-|)	а т т

В силу сделанного предположения все угловые миноры матрицы Pm_,Am., отличны от нуля: последний угловой минор Am= det(Am) = det(A) *0 также отличен от нуля. В

рассмотренном частном случае теорема доказана. 2. Пусть теперь det(Aml) = 0.

Представим разложение определителя матрицы А в виде разложения по последнему столбцу:

det(A) = d et(A J = ±a,mdet(A ^,)±aJm det(A(m!’,)± ...± a mm det(A'ml,) .

В силу det(A)*0, существует хотя бы один det(A(J^,) * 0 . Поменяем местами

строки матрица А с номерами ш и k (m * k в силу принятого условия det(Am.,) = 0). В результате получаем преобразованную матрицу А^ = РтАга, у которой угловой минор

Д;., = det(А о т л и ч е н от нуля. Теперь, в соответствии		с доказательством		случая 1,
уже все угловые миноры матрицы А*т отличны от нуля.
Следствие. Если det(A)*0, то существует матрица перестановок Р			такая, что
справедливо разложение
	РА = LU,
где L	нижняя треугольная матрица с ненулевами	коэффициентами	на	главной

диагонали: U - верхняя треугольная матрица с единицами на главной диагонали.

Определение числа операций алгоритма метода Гаусса

Подсчитаем число операций умножения и деления (наиболее длительные операции

при вычислениях на ЭВМ),		необходимое для реализации алгоритма метода Гаусса.
1.	Вычисление	коэффициентов	с , i = l,m, j = i + l,m	(сумма членов
арифметической прогрессии):

2.	Вычисление коэффициентов а‘,*\ k = 2,m, i = k,m, j = k,m (сумма слагаемы

арифметической прогрессии 2-го порядка [6]):

(m - l)-(m - l) + (ro-2)-(m -2) + ... + l-l = I+ 21 + ...+ ( т - 1 )2 = т '(т ~ 1) ( 2 т ~ 0

					6
3. Вычисление значений f,‘k),		k = 2,m,		i = k,m
,	.ч	/	-.ч	,	m -(m -l)
(ш —1)+ (m —2) + ... +1					= — ^

4. Вычисление значений у,. i = l.m производится m раз.

5. Вычисление х(, i = l,m при “обратном ходе” метода Гаусса:

, _ . , гп• (m 1) 1+ 2+ ... + (т -1) = — *----- ' .

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

	_	,	m (m -l) ,	m (m -l) (2m -l)	m / _ 2 .		Л
	ш + 3		------ ------+ -----------	------------= j	(ш	+Зш -1).
Иными словами,			количество	операций умножения		и	деления приблизительно
пропорционально	ПТз						^
		. что и определяет затраты на выполнение метода Гаусса.

Вычисление определителя матрицы

Вычитание строк в методе Гаусса (образование линейных комбинаций уравнений) не изменяет значения определителя матрицы [7]. Знак определителя может измениться лишь при использовании процедуры выбора “главного” элемента, когда производится перестановка строк (столбцов) преобразуемой матрицы. В результате выполнения всех необходимых преобразований метода Гаусса определитель исходной матрицы может быть вычислен достаточно просто:

det(A) = det(LU) = det(L) det(U) = det(L) = П a ,j,■', j=>

Здесь учтено, что L и U треугольные матрицы и, кроме того, матрица U содержит только единицы на главной диагонали.

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

Построение обратной матрицы

Пусть a pq, p,q = l,m коэффициенты обратной матрицы А"1 Согласно определению

1 > и “ к1=8„, 6 = { U = j		символ Кронекера1.
к=1
Теперь q-й столбец	обратной матрицы можно рассматривать как результат
решения системы линейных алгебраических уравнений вида
	а ИЧ)	А »'
	а 2(Ч)
	. = •
	ttp(q)

Ат>.

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

1Кронекер Леопольд [7.12.1823 - 29.12.1891] - немецкий математик. С 1861 года - член Берлинской академии наук; с 1872 года стал члеаом-корреспондентом Петербургской академии наук. В 1883 году занял должность профессора Берлинского университета.

Метод квадратного корня

Метод квадратного корня предназначен для решения систем линейных

алгебраических уравнений

вида

Ах =

f

с симметричной

матрицей

коэффициентов

„ = а„, i,j= 1,ш

Метод основан на разложении матрицы коэффициентов А в произведение

А = ST • D • S .

(2.4)

где S

верхняя треугольная

матрица с

положительными

значениями на

главной

диагонали;

D - диагональная матрица со значениями +1 или -1.

Согласно теореме 2.1

при неравенстве нулю всех угловых миноров матрицу Л

можно разложить в произведение А = LU.

Представим нижнюю треугольную матрицу L

с ненулевыми коэффициентами на

главной диагонали в виде произведения

N K ,

где

N - нижняя треугольная матрица с

единицами на главной диагонали, К - диагональная матрица, причем kn = Хм,

i = l,m :

'*п

0

0

0

"1

0

0

0'

0

0

0

* 2 .

>■2»

0

0

П2,

1

0

0

0

0

0

L = К

*22

0

=

1

0

0

0

* 3 3

0

Хт1

K z

К г

Кш

Пт ,

Пт2

П т З

1

0

0

0

Jt„BJ

После перемножения

матриц

N

и

К

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

относительно величин n4j,

i = 2,m,

j = l,m -l:

= п 21Х п ;

Х 3 | = П 3 1 * 1 Р

* 3 2 - П 3 2 * 2 2 ’

*•4,

= П 4 1 * 1 Р

* 4 2

= П 4 2 * 2 2 ’

* 4 3 — П 4 3 * 3 3 *

Очевидно, что п2) = ^ - ,

п31 = —

,

п „

= —

,

*п

*п

*22

С учетом этого коэффициенты n(j, kSj матриц N и К можно представить в общем виде