Спецкурс численные методы.

101

ЧМ-8. АЛГЕБРАИЧЕСКАЯ ПРОБЛЕМА СОБСТИЕННЫХ ЗНАЧЕНИЙ И СОБСТВЕННЫХ ВЕКТОРОВ КВАДРАТНОЙ МАТРИЦЫ

8.1. Входная информация для самопроверки.

Для изучения данной темы Вам необходимо восстановить в памяти:

-из курса прикладной математики – решение систем линейных однородных и неоднородных алгебраических уравнений, понятия линейной и векторной алгебры, собственные значения и собственные вектора квадратной матрицы, системы линейно независимых векторов;

-из настоящего спецкурса - понятия: линейного пространства, подпространства и нормированного пространства.

8.2. Содержание темы.

Напомним некоторые сведения из курса линейной алгебры.

Если А — квадратная матрица п-го порядка и Ax x при x 0, то число называется собственным, значением матрицы, а ненулевой вектор х — соответствующим ему собственным вектором. Перепишем задачу в виде

A E x 0, x 0 (8.1.)

Для существования нетривиального решения задачи (8.1) должно выполняться условие

det A E 0	(8.2)
Этот определитель является многочленом	n-й степени относительно ; его называют

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

Если найдено некоторое собственное значение, то, подставляя его в однородную систему (8.1), можно определить соответствующий собственный вектор. Будем нормировать собственные векторы. Тогда каждому простому (не кратному) собственному значению соответствует один (с точностью до направления) собственный вектор, а совокупность всех собственных векторов, соответствующих совокупности простых собственных значений, линейно – независима. Таким образом, если все собственные значения матрицы простые, то она имеет п линейно-независимых собственных векторов, которые образуют базис пространства .

Пример.

Заданы осевые и центробежные моменты инерции некоторого пространственного тела в центральных осях: Jxx =6 кн м2 , Jyy =5 кн м2 , Jzz =7 кн м2 , Jxy =-2 кн м2 ,

Jxz =2 кн м2 , Jyz =0. Найти значения главных центральных моментов инерции и их

направления. Решение.

Составим матрицу центральных моментов инерции

102

J			J		J		6		2	2
		xx		xy		xz			5
A= Jyx			Jyy		Jyz = 2					0 .
	Jzx		Jzy					2	0	7
					Jzz

Запишем характеристическое уравнение (8.2) и найдем собственные значения этой матрицы:

	6	2	2
	2	5	0	=0.
	2	0	7
Раскрывая определитель , получим:
(6- )(5- )(7- ) – 4(5- )-4(7- )=0,
( 3)( 6)( 9) 0.
Отсюда находим собственные значения матрицы: 1 =3,					2 =6,	3 =9.
Найдем собственный вектор x1,		соответствующий			первому собственному

значению 1 =3 из уравнения (8.1)

(A-3E) X1=0.

Запишем это уравнение в координатной форме

3x1 2x2			2x3		0
		2x2			0,
3x1
2x			4x	3	0
	1

Преобразуем данную систему линейных алгебраических уравнений по методу Гаусса:

	3	2	2		3	2	2	1	0	2	1 0		2
		2				1			1			1
	2		0 1				0 0			2 0			2
		0	4			0			2			0
	2				1		2	0		4	0		0
Таким образом,		система имеет решение:						x1	2x3 ,		x2 2x3 ,		а собственный
вектор x1 2,	2,	1 .				x2 ,
Найдем собственный				вектор			соответствующий				второму		собственному
значению 2 =6 из уравнения (8.1)

(A-6E) X2 =0.

103

Запишем это уравнение в координатной форме

		2x2 2x3		0
				0,
2x1 x2
2x		x	3	0
	1

Преобразуем данную систему линейных алгебраических уравнений по методу Гаусса:

0		2 2			4		2	0		0		0	0
		1					1					1
2			0 2					0 2					0 .
	2	0	1			2	0	1			2	0	1

Таким образом,		система имеет решение: x2 2x1 ,	x3 2x1 ,	а собственный
вектор x2 1,	2,	2 .
Найдем	собственный вектор x3 , соответствующий		второму	собственному

значению 3 =9 из уравнения (8.1)

(A-9E) X3 =0.

Запишем это уравнение в координатной форме

2x1 4x2					0,
2x		2x2	2x	3	0
	1				0
3x1			2x3

Преобразуем данную систему линейных алгебраических уравнений по методу Гаусса:

3	2	2	3	2	2	0	2	1	0		0	0
	4			2							2
2		0 1			0 0 2			1 . 0				1 .
				0						1 0
2	0 2		1		1	1 0		1				1
Таким образом, система имеет решение:						x1	2x2	,	x3 2x2 ,		а собственный
вектор x3 2, 1,	2 .
Однако, если	степень		многочлена достаточно					высокая,		то		аналитически

определить корни достаточно сложно. И в этом случае их определяют численными методами.

Если собственное значение известно, то собственный вектор .удовлетворяет системе (8.1). Но любой численный метод дает вместо точного собственного значения

і приближенное значение ~і , так что det( A ~і E ) 0, хотя отличается от нуля

104

очень мало. В таком случае задача ( A ~і E )x 0 при использовании приближенного

собственного значения имеет только тривиальное решение x=0. Поэтому при численных расчетах находить собственные векторы непосредственно из системы (8.1) нельзя.

1. Метод обратных интеграций.

Для нахождения собственных векторов удобен метод обратной итерации, заключающийся в следующем. Выберем наудачу вектор b и рассмотрим линейную неоднородную систему

	E x b	(8.3)
A i

Определитель этой системы отличен от нуля, так что она имеет единственное решение. Покажем, что найденный из нее вектор x окажется почти равным собственному вектору xі соответствующему данному собственному значению і .

Для простоты ограничимся случаем, когда матрица n-го порядка имеет n линейнонезависимых собственных векторов xі — например, матрица нормальная (случай

произвольных матриц рассмотрен ниже). Тогда собственные векторы образуют базис, по которому можно разложить векторы x и b:

n	n
x jxj,	b jxj	(8.4)
j 1	j 1

Подставляя это разложение в систему (8.3), перенося все члены влево и учитывая, что Axj j xj , получим

n	j
						x		0	(8.5)
			j				j
				i		j

j 1

Поскольку собственные векторы линейно – независимы, то их линейная комбинация обращается в нуль только в том случае, когда все ее коэффициенты равны нулю. Поэтому из (8.5) следует

j

j (8.6)

j i

Видно, что если j ~i - мало, то коэффициент j будет очень большим; в противном случае он невелик. Рассмотрим следствия из этого в трех основных случаях.

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

Это означает, что найденный вектор x почти совпадает с собственным вектором xi (с

точностью до нормировочного множителя), что и требовалось доказать. Заметим, что

105

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

Из (8.6) видно, что при обратной итерации (т. е. при переходе от b к x) компонента i усиливается по сравнению с другими компонентами примерно во столько раз, во сколько погрешность данного собственного значения меньше разности соседних собственных значений. Поэтому чем точнее найдено ~i (очевидно, хорошая точность особенно важна при наличии близких собственных значений), тем ближе x будет к xi . Если собственные значения найдены слишком грубо, или случайно вектор b выбран неудачно, так что i - очень мало, то разница между x и xi может оказаться

заметной. Тогда подставляют найденный вектор x в правую часть уравнения (8.3) вместо b и организуют итерационный процесс

	E x	S	x	S 1	, x	0		b	(8.7)
A i

Обычно он сходится настолько быстро, что двух итераций вполне достаточно. Напомним, что на каждой итерации обязательно надо нормировать найденные x( S ) чтобы не получать в расчетах слишком больших чисел, вызывающих переполнение на ЭВМ.

Замечание 8.1. Очень эффективен один простой способ выбора b. В качестве его компонент в декартовых координатах возьмем последовательные многоразрядные псевдослучайные числа k . Тогда вероятность того, что i окажется очень малым,

будет ничтожна.

Второй случаи	— собственное		значение i		кратно; например, 1	2 ... p ,
1 p n. Напомним,	что в этом случае собственные векторы x1 ,x2 ,...,xp					определены
неоднозначно; любая их линейная комбинация удовлетворяет уравнению
		p		p	p
	A	jxj		j Axj	1 jxj
	j 1			j 1	j 1

и являются собственным вектором. Т.е. они порождают p – мерное подпространство, любой базис которого можно взять в качестве системы искомых собственных векторов.

Теперь из (8.6) следует, что большими оказываются коэффициенты 1, 2, , p ,

причем степень их усиления одинакова; остальные коэффициенты остаются малыми. Значит, найденный из (8.3) вектор х будет приближенно линейной комбинацией x1,x2, ,xp , а тем самым — искомым собственным вектором. Если точность полученного

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

столько линейно-независимых векторов b k , какова кратность корня. Обратными

итерациями получим столько же векторов x k , которые и будут искомыми; они будут линейно-независимыми, поскольку преобразование (8.3) невырожденное. Остается только ортогонализовать найденные векторы, если это требуется по условиям задачи.

det A iE 0, то при

106

Напомним, что в качестве декартовых координат векторов b k целесообразно брать псевдослучайные числа; здесь это имеет то дополнительное преимущество, что векторы автоматически получаются линейно-независимыми.

Третий случай — когда матрица имеет кратные корни, но число ее собственных векторов меньше п — выходит за рамки нашего доказательства. Однако метод обратных итераций здесь также применим в той форме, которая описана для кратных корней. Разница лишь в том, что если p-кратному собственному значению соответствуют всего

q собственных векторов q p , то из полученных обратной итерацией векторов

x k только q будут линейно-независимыми. Это выясняется при их ортогонализации:

первые q векторов ортогонализуются «без приключений», а при ортогонализации

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

Каков объем расчетов в методе обратной итерации? Нахождение собственного вектора требует (при одной итерации) не более 23n3 действий, так что для нахождения всех

их надо около n4 арифметических действий. Таким образом, при больших порядках матрицы метод неэкономичен, но при n 10 вполне удовлетворителен. Особенно употребителен этот метод из-за своей простоты, универсальности и хорошей устойчивости алгоритма.

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

одного собственного вектора в этом случае требуется 10n, а для всех—10n2 арифметических действий.

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

то нахождение каждого собственного вектора требует 23n3 действий (тот же прием для трехдиагональной матрицы позволяет сократить число действий до 7n). А для нахождения всех собственных векторов требуется соответственно 23n3 арифметических действий.

Отметим одну существенную деталь. Поскольку

нахождении собственных векторов в формулах прямого хода метода исключений (прогонки) на главной диагонали появится хотя бы один очень малый элемент. Чтобы формально можно было вести расчет, диагональные элементы не должны обращаться в нуль; для этого надо, чтобы погрешность собственного значения была не слишком мала, т. е. составляла бы 10—15 последних двоичных разрядов числа на ЭВМ. Если корни характеристического многочлена находят методом парабол (или секущих), то такая погрешность получается естественно, ибо из-за ошибок округления эти методы перестают сходиться в очень малой окрестности корня. Но если корни определялись методом Ньютона, то при этом могли быть найдены верно все знаки собственного значения; тогда,

107

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

2. Степенной метод .

Степенной метод (счет на установление) применяется для получения наибольшего по

модулю собственного значения. Пусть		1		2		3	. Построим такой
итерационный процесс:
x S 1	Ax S				(8.8)

Он не сходится в обычном смысле. Разложим нулевое приближение по собственным

векторам матрицы:

x 0

ixi

. Тогда легко убедиться,

что x S

iS ixi

и при

i

i

достаточно большом

числе итераций x S S x , т.

е.

вектор

x S

сходится к

i

i i

собственному вектору по направлению. Очевидно, при этом x S 1 x S .

q

1

Процесс сходится

линейно

со знаменателем

2

1

. Считается,

что

процесс

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

	1				x S 1				x S 1 ,x S 1	(8.9)
									x S ,x S
						x S

Отметим, что при расчетах на ЭВМ на каждой итерации после вычисления 1 вектор x S 1

надо нормировать, чтобы не получать переполнений или исчезновений чисел.

Формально при 1 0 итерации сходятся к следующему собственному значению.

Однако из-за ошибок округления 1 не может быть точно нулем, а при малом 1

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

Если наибольшее собственное значение кратное, но соответствующий элементарный делитель матрицы линеен, то итерации сходятся обычным образом. Но если 1 2 , а их

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

Если 1 2 , то сходимость очень медленная; этот случай нередко встречается в

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

Одна итерация для матрицы общего вида требует 2n2 арифметических действий, а для ленточной матрицы — 2тп действий. Из-за медленной сходимости степенной метод применяют только к матрицам, содержащим очень много нулевых элементов (и даже к ним— довольно редко).

108

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

квадратичную сходимость: x S A x 0	, где	A A		A	и A A. Однако если
S		S	S 1	S 1	0

матрица А имеет много нулевых элементов, то ее степени уже такими не будут. Поэтому этот вариант обычно не экономичен.

3. Обратные итерации со сдвигом.

Напишем итерационный процесс, обратный по отношению к степенному процессу:

x S 1 A 1x S (8.10)

Очевидно, он сходится в указанном в предыдущем пункте смысле к наибольшему по модулю собственному значению матрицы A 1, т. е. к наименьшему по модулю собственному

значению матрицы А (ибо собственные значения матриц А и A 1 обратные друг другу). Все, что говорилось в предыдущем пункте о характере сходимости, разумеется, справедливо и в этом случае; сходимость будет довольно медленной.

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

наименьшее, собственное значение i Тогда так называемая сдвинутая матрица

A iE будет иметь собственные значения i . У этой матрицы интересующее нас

собственное значение i будет намного меньше по модулю, чем остальные. Поэтому

обратные итерации со сдвинутой матрицей (которые мы запишем в несколько иной форме)

	E x	S 1	x	S	,	(8.11)
A i

будут быстро сходиться и определят требуемое нам собственное значение i . Напомним,

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

S

x

S

,

A iE y

S

xkS

x

S 1

y S

(8.12)

i

i

ykS

,

y S

Здесь индекс k относится к компонентам векторов, а скобки означают некоторое

усреднение по всем компонентам: например, среднеарифметическое.

Если исходное приближение было хорошим, то иногда процесс сходится за несколько итераций; тогда выгодно непосредственно решать линейную систему (8.10). Если же

требуемое число итераций велико, то лучше обратить матрицу A iE . Выгодней всего

при решении линейной системы (8.11) методом исключения Гаусса использовать полученные

109

на первой же итерации вспомогательные коэффициенты cmk на каждой последующей

итерации; но это не предусмотрено в обычных стандартных программах.

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

S

E y

S

x

S

,

A i

S

y

S

(8.13)

S 1

S

xk

,

x S 1

ykS

y S

i

i

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

допускать слишком точное совпадение i с собственным значением нельзя, ибо матрица системы (8.13) становится плохо обусловленной; в связи с нахождением собственных

векторов. Поэтому, когда в ходе итераций у величины S	устанавливаются (т. е.
i

перестают меняться) 5—7 знаков, то итерации следует прекращать.

Замечание 8.2.. Переменный сдвиг собственного значения (8.13) нельзя включать с первой итерации; сначала надо получить грубую сходимость итераций с постоянным сдвигом.

Замечание 8.3. Обратные итерации особенно удобны, если матрица заранее приведена преобразованием подобия к почти треугольной форме. Тогда одна обратная

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

Выводы. Обратные итерации с постоянным и особенно с переменным сдвигом — очень эффективный метод расчета. Для нахождения собственных векторов этот метод

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

Все перечисленные методы реализованы на ЭВМ и их можно найти в пакетах прикладных программ, например: MathCAD, MathLAB.

8.3. Критерий усвоения

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

-математическую постановку задачи о нахождении собственных значений и собственных векторов квадратной матрицы;

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

В результате изучения данной темы Вы должны знать:

-существуют численные методы для нахождения собственных значений и собственных векторов квадратной матрицы;

-их достоинства и недостатки;

110

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

-исследуя физическую постановку задачи, составить ее математическую модель;

-на основе анализа математической модели, выбрать наиболее оптимальный численный метод ее решения;

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

8.4. Выход темы в другие темы и дисциплины

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

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

8.5.Тест - контроль для самопроверки

8.1.Если выполняется соотношение Ax x, где A - квадратная матрица, то число называется:

А. собственным вектором матрицы;

. Б собственным значением матрицы; В. элементом матрицы.

8.2.Если выполняется соотношение Ax x, где A - квадратная матрица, то x называется:

А. собственным вектором матрицы; Б. собственным значением матрицы;

В. элементом матрицы 8.3. С помощью метода обратных итераций можно определить:

А. собственным значением матрицы; Б. собственный вектор матрицы;

В. собственное значение и собственный вектор матрицы. 8.4. С помощью степенного метода можно определить:

А. наибольшее по модулю собственное значение матрицы; Б. собственный вектор матрицы; В. все собственные значения матрицы.

8.5. С помощью метода обратных итераций со сдвигом можно определить: А. наибольшее по модулю собственное значение матрицы; Б. собственный вектор матрицы; В. все собственные значения матрицы

Ответы на тест – контроль

8.1.«Б». - собственным значением матрицы;

8.2.«А» собственным вектором матрицы.

8.3.«Б» собственный вектор матрицы.

8.4.«А» наибольшее по модулю собственное значение матрицы.

8.5.«В» все собственные значения матрицы.