矩陣的特征值與特征向量問題

1、第四章矩陣的特征值和特征向量問題物理、力學和工程技術中的許多問題在數(shù)學上都歸結為求矩陣的特征值和特征向量問 題.計算方陣A的特征值，就是求特征方程A - 入 I =0即-nn_1° n _2_扎+ Pl扎+ P+十Pn = 0的根.求出特征值-后，再求相應的齊次線性方程組(A 丸I厭=0的非零解，即是對應于的特征向量.這對于階數(shù)較小的矩陣是可以的，但對于階數(shù)較大的矩陣來說，求解是十分困難，所以用這種方法求矩陣的特征值是不切實際的我們知道,如果矩陣A和B相似，則A和B有相同的特征值.因此人們就希望在相似變換 下，把A化為最簡單的形式.一般矩陣的最簡單的形式是約當標準形由于在一般情況下，

2、用相似變換把矩陣A化為約當標準形是很困難的，于是人們就設法對矩陣 A依次進行相似變換 使其逐步趨向于一個約當標準形，從而求出A的特征值.本章介紹求部分特征值和特征向量的幕法，反幕法；求實對稱矩陣全部特征值和特征向量的雅可比方法；求特征值的多項式方法；求任意矩陣全部特征值的QR方法第一節(jié)幕法和反幕法一 冪法幕法是一種求任意矩陣 A的按模最大特征值及其對應特征向量的迭代算法該方法最大的優(yōu)點是計算簡單，容易在計算機上實現(xiàn)，對稀疏矩陣較為合適，但有時收斂速度很慢為了討論簡單，我們假設(1) n階方陣A的特征值1，2,，n按模的大小排列為I九1 |>|九2徉糾九n I(1)(2) Vi是對應于特征

3、值i的特征向量i T，2,，n ；(3) Vl，V2，,Vn線性無關.任取一個非零的初始向量Xo ,由矩陣A構造一個向量序列x1 = Ax0 x2 二 Ax1Xk 二 AX"稱為迭代向量由于V1，V2,vn線性無關，構成n維向量空間的一組基，所以，初始向量X。 可唯一表示成x° 二 aMa2V2a.Vn于是Axk 二二 AxkJ =二 A x0-k丄k亠丄_ ka1 1V1 a2 2V2an n£因為比值k * 1a-|V|a2：1 i =23, ,n ,所以k'2V2'1an.Xkkm】=a1v1當k充分大時有n k1 ay從而這說明當k充分大時

4、，兩個相鄰迭代向量Xk 1和Xk近似地相差一個倍數(shù)，這個倍數(shù)便是矩陣A的按模最大的特征值 若用Xki表示向量Xk的第i個分量，則.(X")?L1 屯(xk(8)也就是說兩個相鄰迭代向量對應分量的比值近似地作為矩陣A的按模最大的特征值因為xk 1 " ' 1xk,又Xk 1二AXk,所以有AXk 1 Xk,因此向量Xk可近似地作為對應于1的特征向量這種由已知的非零向量 X。和矩陣a的乘幕構造向量序列'Xk以計算矩陣a的按模最大 特征值及其相應特征向量的方法稱為冪法由(4)式知，幕法的收斂速度取決于比值2'1的大小.比值越小，收斂越快，但當比值接近于1時

5、,收斂十分緩慢人>1,則迭代向量Xk的各個不為零的分量將隨著k無限增用幕法進行計算時，如果九1 <1,則Xk的各分量將趨于零這樣在有限字長的計算機上計大而趨于無窮反之,如果算時就可能溢出停機為了避免這一點，在計算過程中，常采用把每步迭代的向量Xk進行規(guī)范化,即用Xk乘以一個常數(shù)，使得其分量的模最大為 1.這樣，迭代公式變?yōu)閗 = Ax“ mk = max(yk )( k = 1,2,)(9)其中mk是yk模最大的第一個分量相應地取婦止mkXk 或（Yk_21-101A =-12111110-12 J(10)用幕法求其模為最大的特征值及其相應的特征向量（精確到小數(shù)點后三位）。取 x

6、1,1,1 T1的大小.當比值接近于1時，收斂可能kTYkmkXk1101110122-2221-1133-43-4-0.751-0.754-2.53.5-2.53.5-0.7141-0.7145-2.4283.428-2.4283.428-0.7081-0.7086-2.4163.416-2.4163.416-0.7071-0.7077-2.4143.414-2.4143.414-0.7071-0.707解計算結果如表4-1所示。表4-1當k=7時，Xk已經(jīng)穩(wěn)定，于是得到 1 - - m - 3.414及其相應的特征向量 W為M 吒 x? =（-0.707,1,-0.707使用幕法時，應注意

7、以下兩點：（1）使用幕法時，困難在于事先不知道特征值是否滿足（1）式，以及方陣A是否有n個線性無關的特征向量.克服上述困難的方法是：先用幕法進行計算，在計算過程中檢查是否出現(xiàn)了預期的結果.如果出現(xiàn)了預期的結果，就得到特征值及其相應特征向量的近似值；否則，只能用其它方法來求特征值及其相應的特征向量（2）如果初始向量X。選擇不當，將導致公式（3）中W的系數(shù)1等于零.但是，由于舍入誤k差的影響，經(jīng)若干步迭代后，Xk = A X。.按照基向量v，v2，.，vn展開時，v的系數(shù)可能不等于零。把這一向量Xk看作初始向量，用幕法繼續(xù)求向量序列Xk 1，Xk 2，.，仍然會得出預期的結果，不過收斂速度較慢.如

8、果收斂很慢，可改換初始向量.二原點平移法由前面討論知道，幕法的收斂速度取決于比值 很慢.這時，一個補救的方法是采用原點平移法.設矩陣(11)其中P為要選擇的常數(shù)我們知道A和B除了對角線元素外,其它元素都相同,而A的特征值'i和B的特征Ji 之間有關系 "=一 P,并且相應的特征向量相同.這樣，要計算A的按模最大的特征值， 就是適當選擇參數(shù) p,使得一 p仍然是B的按模最大的特征值,且使1對B使用幕法，使得在計算 B的按模最大的特征值1 一 p的過程中得到加速，這種方法稱 為原點平移法.2設4階方陣A有特征值再=15 i,i = 1,2 3 4r比值= 0.9,令P =12作變

9、換B 二 A - pl則B的特征值為卩=12,1,0,-1使用幕法計算B的按模最大的特征值時，確定收斂速度的比值為二 0.50.91 一 P所以對B使用幕法時，可使幕法得到加速。p值，可以使得幕法得到加速，但由于矩陣的特征值的分布情況事先并不知道，所以在計算時，用原點平移法有一定的困難下面考慮當 A的特征值為實數(shù)時，如何選擇參數(shù) p，以使得用幕法計算1時得到加速 的方法設A的特征值滿足1 2 -3 - 'n4 'n則對于任意實數(shù)p, B = A 一 pl的按模最大的特征值1 一 p或n 一 p。如果需要計算及V1時，應選擇p使0 1且確定的收斂速度的比值r = max*二 mi

10、n人1 - p J2'n、 P = 、p ,即2 時，r為最小.這時用幕法計算1及V1時得到加如果需要計算'n及Vn時，應選擇P使且確定收斂速度的比值當.1 _ P = - n 1加速2 時，r為最小這時用幕法計算'1及V1時得到原點平移的加速方法，是一種矩陣變換方法這種變換容易計算，又不破壞A的稀疏性, 但參數(shù)p的選擇依賴于對A的特征值的分布有大致了解 三反幕法反幕法用于求矩陣 A的按模最小的特征值和對應的特征向量，及其求對應于一個給定的近似特征值的特征向量設n階方陣A的特征值按模的大小排列為相應的特征向量為>扎2V1，v2，,vn 則A J的特征值為對應的特征

11、向量仍然為 V1，V2,必.因此,計算矩陣A的按模最小的特征值，就是計算 A'1的按模最大的特征值這種把幕法用到 A 上,就是反幕法的基本思想任取一個非零的初始向量 Xo，由矩陣A，構造向量序列Xk = A"*Xk4， k = 1,2，(12)用(12)式計算向量序列 '兀匚時，首先要計算逆矩陣 A，.由于計算 A時，一方面計算1 1麻煩，另一方面當 A為稀疏陣時，A 不一定是稀疏陣，所以利用 A 進行計算會造成困難 在實際計算時，常采用解線性方程組的方法求X-( 12)式等價于AXk = Xj，k = 1,2，為了防止溢出，計算公式為“ mk = max(yk )N

12、 = yk/mkk =1,2,(14)(13)相應地取Vn(13)式中方程組有相同的系數(shù)矩陣A=LU再解三角形方程組1& mJ” yk或 Vn ： XkA,為了節(jié)省工作量(15)，可先對矩陣A進行三角分解k = 1,2，當A是三對角方陣，或是非零元素較少且分布規(guī)律的方陣時(17)，無論存儲或計算都比較便(16)根據(jù)幕法的討論，我們知道，在一定條件下，可求得AJ的按模最大的特征值和相應的特征向量,從而得到A的按模最小的特征值和對應的特征向量，稱這種方法為反幕法反幕法也是一種迭代算法，每一步都要解一個系數(shù)矩陣相同的線性方程組設p為任一實數(shù),如果矩陣A - pl可逆,則 A pl的特征值為丸

13、1 一 P九2 - P 丸n - P對應的特征向量仍為Vl ,V2，,Vn .如果p是矩陣A的特征值'i的一個近似值，且% - p| v 九j - p ,i 鼻 j1則i P是矩陣 A - pl的按模最大的特征值.因此,當給出特征值的一個近似值p時,可對矩陣A - pl使用反幕法，求出對應于i的特征向量.反幕法迭代公式中的 yk 通過方程組(A pl ”k = Xkjk = 1,2,求得例3 用反幕法求矩陣_2 -1 0 0 1!-1 2 -1 0 !A =|0 -1 2 -10 0 -1 2一的對應于特征值=0.4的特征向量解取X0 = 1,1,11解方程組A - 0.4l 屮二 X

14、0得% = (-40,-65廠65廠40$,g = max(% )= -65,1X1 = (813,1,1,813 f.再解方程組A - 0.4l y2 二 X1得y2 二-445 13,-72013,-720 13,-44513丁，mb 二 max 七 -720 13,x2 =丄 y2 = (89144,1,1,89144 T .X1和x2的對應分量大體上成比例，所以對應于=0.4的特征向量為V 二 89 144,1,1,89 144 T.第二節(jié) 雅可比方法雅可比方法是用來計算實對稱矩陣A的全部特征值及其相應特征向量的一種變換方法.在介紹雅可比方法之前，先介紹方法中需要用到的線性代數(shù)知識和平

15、面上的旋轉變換一 預備知識（1） 如果n階方陣A滿足ata=i即 A = A則稱A為正交陣.（2）設A是n階實對稱矩陣，則 A的特征值都是實數(shù)，并且有互相正交的n個特征向量.（3）相似矩陣具有相同的特征值.（4）設A是n階實對稱矩陣，P為n階正交陣，則B=PtAP也是對稱矩陣.（5）n階正交矩陣的乘積是正交矩陣.（6）設A是n階實對稱矩陣，則必有正交矩陣P，使PtAP 二(1)其中上的對角線元素的是A的n個特征值，向量.由（6）可知，對于任意的n階實對稱矩陣正交陣P的第i列是A的對應于特征值i的特征A，只要能求得一個正交陣P，使PT AP二丄（上的理論基礎. 二旋轉變換 設為對角陣）,則可得到

16、A的全部特征值及其相應的特征向量,這就是雅可比方法為二階實對稱矩陣A二 a a12_a2122 -,即日12 =日21.因為實對稱矩陣和二次型是2f,x2 二 a11x12a12x1x2對應的2a?2 X?,設A對應的二次型為由分析幾何知識知道， 方程f x1,X2 = C表示在x1,X2平面上的一條二次曲線.如果將坐標軸 0冶，0X2旋轉一個角度二，使得旋轉后的坐標軸Oy1 ,°y2和該二次曲線的主軸重合,如圖4-1 所示，n21%這個變換就是xjl變換(4)把坐標軸進行旋轉COST P 二si ncost -siny1sin costy2 一,所以稱為旋轉變換-sinCOST稱為

17、平面旋轉矩陣。顯然有ptapPTP二I ,所以1.其中(5)P是正交矩陣.上面的變換過程即-一2 一.由于a22 sin2 va12sin2寸Ta11 cos2 二PTAP 二 IL0.5 a2 a11 sin2) a12 cos2所以只要選擇"滿足0.5 a22 一 a11 sin2a11 sin2 日 + a22 cos2a12 cos2a12 sin2& 一1a22 一 a sin2)a12cos2)- 022ai2tan2)=a11 a22(6)3T6 =(當a11=a22時，可選取 4 ) PtAP就成對角陣，這時 A的特征值為a11 cosa22sin2ja12

18、sin2-2 = a sin2 丁a22 cos2 丁 - a12 sin2-相應的特征向量為COSV1Sin日- si nV2IL cost雅可比方法雅可比方法的基本思想是通過一系列的由平面旋轉矩陣構成的正交變換將實對稱矩陣逐步化為對角陣，從而得到A的全部特征值及其相應的特征向量n.首先引進R中的平面旋轉變換變換x = yi cost-yi si n片=y sin日+ yi cos日.兀=ykkH i, j記為x = py,其中11COST -:-sin日i-si n日COSTj1+1 一(8)ijx = ( x1,x2,.,Xn $y =( yy,.,yn Jx = P y n x x則稱

19、j y為r中 j平面內的一個平面旋轉變換Px x,ij稱為f ' j平面內的平面旋轉矩P陣.容易證明 j具有如下簡單性質： P為正交矩陣. P的主對角線元素中除第 i個和第j個元素為 COST夕卜，其它元素均為 1 ；非對角線元素 中除第i行第j列元素為一 sin ,第j行第i列元素為sin外,其它元素均為零. PtA只改變A的第i行和第j行元素，AP只改變A的第i列和第j列元素，所以 PtAP只改變A的第i行、第j行、第i列、第j列元素.A=（R）（nA9）a=a 式 0設 'j nn為n階實對稱矩陣，円 ji為一對非對角線元素.令A =PTAP = （af）h則A為實對稱矩

20、陣，且A和A有相同的特征值.通過直接計算知ajj(1)2 2=aii sin 二 ajj cos 二-a0 sin2ija“(1)ij(1)aik(1) 1 /二 aji2( ajj _ ah(1) .二 aki 二 cosii)sin2)cos2ajk sinnk",j(1) a ik I jk(1) aki(1) .二 ajk 二 _ajkSinnajk costk",j二 aki當取二滿足關系式tan 22ajaiajj(10)時，a； =aj' =0，且 a； 2-aj 2ai： 2a, 22aik2aiia2k k = i,ja行 2a2ak1 2二 ak

21、i k,1F,j由于在正交相似變換下，矩陣元素的平方和不變，所以若用平方和，用 SA)表示 A的非對角線元素平方和，則由(11)式得；D(A)=D(A) + 2ai2is(A) = S(A)-2ai2D A表示矩陣(11)A的對角線元素(12)2 2aH=aH cos 二ajj sin 二a0 sin2這說明用ijPj對A作正交相似變換化為 A后，A1的對角線元素平方和比A的對角線元素平方和增加了 2aij , A1的非對角線元素平方和比A的非對角線元素平方和減少了2aij ,且將事先選)=0 A定的非對角線元素消去了 (即ij ).因此，只要我們逐次地用這種變換，就可以使得矩陣 A的非對角線

22、元素平方和趨于零，也即使得矩陣 A逐步化為對角陣.這里需要說明一點：并不是對矩陣A的每一對非對角線非零元素進行一次這樣的變換就能得到對角陣.因為在用變換消去 aij的時候，只有第i行、第j行、第i列、第j列元素在變化， 如果aik或Pj為零，經(jīng)變換后又往往不是零了 .雅可比方法就是逐步對矩陣A進行正交相似變換，消去非對角線上的非零元素，直到將A的非對角線元素化為接近于零為止，從而求得 A的全部特征值，把逐次的正交相似變換矩陣乘起來，便是所要求的特征向量.aij. 一般說來，取絕對值最大的非雅可比方法的計算步驟歸納如下：第一步在矩陣A的非對角線元素中選取一個非零元素 對角線元素；第二步第三步ta

23、n2 氧a 一 ajj求岀日，從而得平面旋轉矩陣 P = Fij ;由公式A = F ARA的元素由公式（9）計算.以A代替A,重復第一、二、三步求岀A2及F，繼續(xù)重復這一過程，直到A（即小于允許誤差）時為止.第四步的非對角線元素全化為充分小R = P；21 2 -1/2#721/V2第五步Am的對角線元素為 A的全部特征值的近似值，只卩2的第j列為對應于_2-10 1A =-12-1110 -12 J的特征值和特征向量.0 =解首先取i =1，j =2,由于a11=a22=2故取4 ,所以個元素）的特征向量.特征值j （ ' j為Am的對角線上第j 例1用雅可比方法求矩陣-1 2-1

24、/V2A = PT AR一-1/丘-1 2再取 i =1,j =3由42)2tan 2門1-2:0.45969 ,cos : 0.88808所以0.88808-0.459690459690.888080.63398-0.32505-0.32505-0.62797-0.627972.36603繼續(xù)做下去A2 二 F2tAF2=IL,直到非對角線元素趨于零 ，進行九次變換后，得中(h) = T 十 P*n+ P2 孫 n_? + Pn(432)0587580.000000.000001A. = 0.000002.000000.000000.000000.000003.41421 一A的對角線元素就

25、是 A的特征值，即0.58758 ,，22.00000 ,，33.41421相應的特征向量為0.50000P.70710 P.50000 1W 二 0.70710 ,2 = 0.00000 飛=-0.70710 .0.50000一0.70710 一0.50000 一相應的特征值的精確值站= 2 2,対=2，丸3 = 2 +72相應的特征向量為12一 12 12_1血,V2 =0,V3 =-1近1/21皿11/2由此可見，雅可比方法變換九次的結果已經(jīng)相當精確了,特別是求得的特征向量正交性很好，所以雅可比方法是求實對稱矩陣的全部特征值及其對應特征向量的一個較好的方法.但由于上面介紹的雅可比方法，每

26、次迭代都選取絕對值最大的非對角線元素作為消去對象，花費很多機器時間.另外當矩陣是稀疏矩陣時，進行正交相似變換后并不能保證其稀疏的性質，所以對階數(shù)較高的矩陣不宜采用這種方法.目前常采用一種過關雅可比方法.這種方法是選取一個單調減小而趨于零的數(shù)列'a即a > a > lim an 06 > a2 >，且n護）作為限值，這些限值稱為 ”關”,對矩陣的非對角線元素規(guī)定一個順序（例如先行后列、自左至右的順序）.首先對限值 日1按規(guī)定的順序逐個檢查矩陣的非對角線元 素,碰到絕對值小于 31的元素就跳過去，否則就作變換將其化為零.重復上述過程，直到所有的非對角元素的絕對值都小

27、于 31為止.再取32 ,a3，類似處理，直到所有的非對角線元素的絕對值都小于3m時,迭代停止.這時的3m應小于給定的誤差限實際運算中常用如下的辦法取限值：對于矩陣A，計算A = A0的非對角線元素平方和S A0，任取N - n,取3k二 S 代 k =1,2,N4. 3多項式方法求特征值問題4.3.1 F-L方法求多項式系數(shù)我們知道，求n階方陣A的特征值就是求代數(shù)方程毋（九）=| a-丸1=0的根。（）稱為A的特征多項式。上式展開為(431)其中pi,p2,pn為多項式（）的系數(shù)。從理論上講，求 A的特征值可分為兩步：第一步 直接展開行列式| A - I求出多項式：（'）；第二步 求

28、代數(shù)方程'（X）=0的根，即特征值。對于低階矩陣，這種方法是可行的。但對于高階矩陣，計算量則很大，這種方法是不適用的。這里我們介紹用 F-L （ Faddeev-Leverrier）方法求特征方程（4.3.2）中多項式，（）的 系數(shù)。由于代數(shù)方程求根問題在第2章中已經(jīng)介紹，所以本節(jié)中解決特征值問題的關鍵是確定矩陣A的特征多項式 （-），所以稱這種方法為多項式方法求特征值問題。記矩陣a= （aj ）n n的對角線元素之和為(4.3.3)trA = aiia 22ann利用遞歸的概念定義以下n個矩陣Bk（k =1,2,，n）:PiB1 -=A,B2-A(B1 - pj ),B3-A(B2

29、- p2 1 ),Bk=A(Bk4 - Pk),Pn= A(BnJ - PnI),P2P3PkPn可以證明,（4.3.4）式中Pk,k - 1,2,，n,即是所求式求矩陣的特征多項式系數(shù)的方法稱為A F-L方法。=trB11 trB221 trB331 trBn n的特征多項式"）相應特征方程為：(4.3.4)的各系數(shù)。用（4.3.4）(4.3.6)(4.3.5)n , nn Jn _2（一1）（丸-PA P2丸一一Pn） = 0 而且可證矩陣A的逆矩陣可表示為41A （Bn4- Pn）Pn例1求矩陣"3 2 41A= 2024 2 3 一的特征值和a1解 用F-L方法求得

30、3 24B!=A=2024 23jPi=trBi= 6_11241B2 = A( Bi Pi I ) = 28242 111P2 = trB 2=152"8 0 0B3 = A( B2 p21) = 0 8 00 0 8 一 1“cp3 = trB 3=8所以A的特征方程為332(-1) ( -6 '-15-8)此方程的根，即特征值為，1 =8, ，2 =_1,，3 =_1214 12 1A(B2-P2l) =P3=01478141 121412從例1中的計算結果可知 B3二P3l. Faddeev曾經(jīng)證明：對n階矩陣A,按(4.3.4) 式計算出的Bn總有Bn = Pn I

31、(4.3.7)4.3.2特征向量求法當矩陣A的特征向量確定以后，將這些特征值逐個代入齊次線性程組(A - I )x=0中,由于系數(shù)矩陣A -訂的秩小于矩陣 A - J的階數(shù)n,因此雖然有n個方程n個未知數(shù)，但實際上 是解有n個未知數(shù)的相互獨立的r個方程(r<n).當矩陣A的所有特征值互不相同時，這樣的問題中要解的齊次方程組中有n-1個獨立方程，其中含有n個特征向量分量，因此特征向量分量中至少有一個需要任意假設其值，才能求出其他特征分量.在計算機中解這樣的齊次線性程組，可用高斯-若當消去法，以便把一組n個方程簡化為等價的一組 n-1個方程的方程組.然而,用高斯-若當消去法簡化一個齊次線性程

32、組時，方程之間不都是獨立的，在消去過程中系數(shù)為零的情況較多 .必需交換方程中未知數(shù)的次序，以避免主元素位置上為零的情況.因此，為了提高精度和避免零元素的可能性，我們總是用主元素措施把絕對值最大的系數(shù)放于主元素位置.例如，假設矩陣A為42-2A =-532一241 一其特征方程為4 -扎 2-2-53-畫 2-2 4 =0展開后為（九一1）（幾一2）（扎一5） =0故特征值分別為F面求特征向量，將 1代入方程組（A-J）x=°中，得3xi +2x2 2x3 = 0« 5xi +2x2 十 2x3 = 0_2xr +4x2 +0x3 =0以-5為主元素,交換上式第一和第二個方程

33、得一5治 +2x2 +2x3 = 0*3為 +2x2 _2x3 = 0_2xr +4x2 _0x3 =0(438)(4.3.9)用高斯-若當消去法消去-5所在列中的xi,并把主元素所在行調到最后，得0x1164小x3 = 05“ 0x1164x3 = 0522X3 = 0% - X2 -55(4.3.10)再以16/5為主元素，消去它所在列中的x2 ,并把主元素所在的行調到最后，得0xi 0x2 0x3 = 0“ x1 +0x2 -1_X320x1 x2 -4%3(4.3.11)因為這就是用高斯若當消去法實現(xiàn)把一組三個方程簡化為等價的一組兩個獨立方程的情形這個等價的方程組包含兩個獨立的方程，而

34、有三個未知數(shù)，所以只要假定其中一個值，則其它兩個值就可以通過兩個獨立方程解出比如，令X3 = T，則得到矩陣a的對應于1的一個特征向量為2_14_ 1對另外兩個特征值的對應特征向量求法和上述對1 =1的推導過程相同計算機中實現(xiàn)求解這樣的齊次線性方程組的消去步驟是，用第3章討論過的高斯-若當消4154(4.3.12)去法的公式，方程組(439)的系數(shù)矩陣經(jīng)過第一次消去后的矩陣B為165(1652'"5以矩陣為方程組(4310)的系數(shù)矩陣，其中省略了有0和1元素的第一列在進行第二次消元之前，要使用完全主元素措施對前兩行進行最大主元素選擇 ，然后再進 行必要的行或列交換每完成一次消

35、元過程，總省略只有0和1元素的第一列，并且計算機僅尋 找矩陣的前n-k行中的最大主元素，其中k是消元過程使用的次數(shù)對(4.3.12)式再進行一次消 元過程，則得到列矩陣B1此矩陣是對應于方程組021，般來說，最后-4-(4.3.13)(4.3.11)的系數(shù)矩陣，不過省略了含0和1元素的前兩列矩陣列的數(shù)目等于矩陣A-'l的階數(shù)和秩的差值.由于方程組(4.3.8)有三個未知數(shù)，兩個獨立方程，所以計算機必須任意給定一個未知數(shù)的 值，以便可以從其他兩個獨立方程中解出另外兩個未知數(shù).為方便，在計算機決定特征向量時要恰當?shù)卦O定任意選取的未知數(shù)的值.例如，令X3 =由方程組(4.3.11)知道，其他

36、兩個分量的值正好能從含x3的非零系數(shù)項得出.為此，從計算機所存儲的最終矩陣中，令B1最上面的0元素為-1,并把它順次調到最下面第三行的位置上，就得到所求的特征向量在工程問題中，從特征方程所求出的特征值，少數(shù)情形也有相同的.一般地，當一個特征方程有k重根時，矩陣A- I的秩可能比其階數(shù)少1,或2,或3,或k，當然對應于的線性無 關的特征向量的個數(shù)也就是1,或2,或3,，或k,下面通過一個特征值對應兩個線性無關特征向量的例子進一步說明計算機求特征向量的方法設矩陣A為_324A=202423j其特征方程為3_k242-人2 0423 -九展開后得('1)2( ' -8)-0所以特征值

37、為為了決定，=-1的特征向量42應1 =，2 = -1,沁=8，將二-1代入方程組(A - J )x=0,得4；X1 12 |X2=04/3使用一次高斯-若當消去法(4314)1 1/2寫成矩陣形式,(4.3.15)式的系數(shù)矩陣為_00X1 10 產2 = 01 X3(4315)0101/2因為方程組(4.3.15)的系數(shù)矩陣的秩為(4.3.16)1,它比矩陣階數(shù)少 2,因此對應于扎=-1有兩個線性無 關的特征向量，必須給兩個未知數(shù)任意規(guī)定值(4.3.15)式可看出，一般總是選擇 X2二7X3 另一個特征向量；這樣有兩個線性無關的特征向量,才能確定這兩個線性無關的特征向量，由=0求一個特征向量

38、；選擇X2 = Q X3 _ - 1求1/2：-1-1 一，在(4.3.16)式的B中，把第一列中第一個0計算機中求兩個線性無關的特征向量的辦法是元素用-1代替,第二列中第二個0元素也用-1代替，然后把第一、第二行順次調到最下面一行 的位置上，第三行自然就成了第一行， 如此調換后矩陣的第一列和第二列就是所求的兩個線 性無關的特征向量。對應于=-1的全部特征向量為刁"11k1 -1 * 0JJ其中K和k2是任意常數(shù)，且不同時為零。為了說明列交換的必要性，避免主元素為零，再舉一個例子，設矩陣8其特征方程為特征值為-21(. '2)1 (九-1) = 0-124，1 = 2 J 2

39、 -0,3 - 1對應于，=2的特征向量可由解下列方程組而求得-1 / 3X3 |1/32/3X2J再進行一次消去過程,2/3得01(4.3.20)0 I%x2 = 01/2_(4.3.21)-4-812xj124|X2=0 100-1X3 一(4.3.17)用一次高斯-若當消去法，得0 01 110 0-1x2=：011 23 一?3 一|(4.3.18)若不進行列交換，則下一個消元過程只能在第一行的第二個元素和第二行的第二個元素中找 最大主元素，而它們都是零，我們不得不對(4317)式進行列交換，即交換未知數(shù)之間的次序，之后再進行消去過程對(4.3.17)式進行列交換，即把絕對值最大系數(shù)放

40、在主元素位置，顯然是第一列和第三 列的交換，交換后成為一12-84 x3 421|x2 =0.T 00 jxi(4 3 19)其中未知數(shù)列矩陣中X1和x3也進行了交換，這樣才能保證(4.3.17)式和(4.3.19)式等價，對(4.3.19)式進行一次高斯-若當消去法，得-2/3在計算機中計算，剩下一個最終的列矩陣_00(4.3.22):1/2J將(4.3.22)式中的列矩陣B中第一個0元素用-1代替，并隨即調到最下面一行，便得到0 1 1/2L -( 4.3.23)這就是對應于方程組(4.3.19)的解，在計算機程序中應把原來進行列交換的列號次序記住，重新把(4.3.23)式中各分量排列一下

41、，即交換第一行和第三行的元素，就得到對應于=2的特征向量-111/20 J對應于的全部的特征向量為'-I!1/2k其中k為不等于零的任意常數(shù)4.4 QR算法QR算法也是一種迭代算法，是目前計算任意實的非奇異矩陣全部特征值問題的最有效的方法之一.該方法的基礎是構造矩陣序列'Ak并對它進行QR分解.由線性代數(shù)知識知道，若A為非奇異方陣，則A可以分解為正交矩陣 Q和上三角形矩陣R 的乘積，即A=QR,而且當R的對角線元素符號取定時，分解式是唯一的.若A為奇異方陣，則零為A的特征值.任取一數(shù)p不是A的特征值，則A-pI為非奇異方陣. 只要求出A-pI的特征值，就很容易求出 A的特征值，

42、所以假設A為非奇異方陣，并不妨礙討論 的一般性.設a為非奇異方陣，令Ai = A，對Ai進行QR分解，即把A1分解為正交矩陣 Qi和上三角形矩陣R1的乘積A1 = Q1 R1做矩陣A2 = R1Q1 =q1 A1Q1繼續(xù)對A2進行QR分解A2二Q2 R2并定義A3 = R2Q2 = Q： A2Q2般地，遞推公式為Ai = A = Q1R1= RQk = Qk AkQk，k = 2,3,.QR算法就是利用矩陣的 QR分解，按上述遞推公式構造矩陣序列%只要A為非奇異方陣，則由QR算法就完全確定'Ak這個矩陣序列 Z具有下列性質.性質1所有人都相似，它們具有相同的特征值 證明 因為Ak 1

43、二 RkQkAkQk=Qk Qk 4 Ak 4Qk jQk =.-Qk Qk 4 - Q1 AQ1Q2 .Qk若令Qk =Q1Q2.Qk，則Qk為正交陣，且有*t*Ak = Qk AkQk因此Ak和A相似，它們具有相同的特征值 性質2 Ak的QR分解式為Ak =QkRk其中 Qk =QQ2Qk，Rk - Rk Rk 4 .R1證明用歸納法.顯然當k=1時，有A = A Q1 R1 Q1R1k假設Ak有分解式A =Qk4Rk4于是Qk Rk - Q1Q2 .Q kd (Qk Rk) Rk.R=Qk j Ak Rk jT因為Ak =Qk丄AQk，所以kQk R = Ak Qk 1 Rk 1 = A

44、因為Q1,Q2,.Qk都是正交陣，所以Qk也是正交陣 解式為，同樣Rk也是上三角形陣，從而Ak的QR分k A 二QkRkT由前面的討論知 Ar 1二Qk AQk .這說明qr算法的收斂性有正交矩陣序列 決定."Qk 的性質定理1如果Qk '收斂于非奇異矩陣 Q :, Rk為上三角形矩陣 三角形矩陣Hu貝H>證明 因為Qk "收斂，故下面極限存在lim Qk = lim QQk =Q；Qk = I k ,k_,R= pRk = jjAk 卅QkT T=lim Qk AQk=QAQqk)：一 一由于Rk(k = 1,2,-)為上三角形矩陣，所以R：:為上三角形矩陣.又因為A.- - lim Ak = lim QkRk_ k 匸k_j：: