第3章矩陣特征值與特征向量的計算

1、 引言引言 在科學(xué)技術(shù)的應(yīng)用領(lǐng)域中，許多問題都?xì)w為求解一個在科學(xué)技術(shù)的應(yīng)用領(lǐng)域中，許多問題都?xì)w為求解一個特特征系統(tǒng)。如動力學(xué)系統(tǒng)和結(jié)構(gòu)系統(tǒng)中的振動問題，求系統(tǒng)的征系統(tǒng)。如動力學(xué)系統(tǒng)和結(jié)構(gòu)系統(tǒng)中的振動問題，求系統(tǒng)的頻頻率與振型；物理學(xué)中的某些臨界值的確定等等。率與振型；物理學(xué)中的某些臨界值的確定等等。3.1 3.1 乘冪法及其變體乘冪法及其變體3.1.1 乘冪法乘冪法定理定理 設(shè)設(shè)A Rn n有完全特征向量系有完全特征向量系，若若 1, 2, n為為A的的n個特征值且滿足個特征值且滿足n21 對任取初始向量對任取初始向量x(0) Rn，對乘冪公式對乘冪公式)1()(kkAxx確定的迭代序列確定的

2、迭代序列xk，有下述結(jié)論有下述結(jié)論： （1）當(dāng)當(dāng) 時時，對對i = 1, 2, , n211)()1(lim kikikxx收斂速度取決于收斂速度取決于 的程度的程度，r 越小收斂越快越小收斂越快，r 1收斂慢收斂慢，112r且且x(k)（當(dāng)當(dāng)k充分大時充分大時）為相應(yīng)于為相應(yīng)于 1的特征向量的近似值的特征向量的近似值。321（2）當(dāng)當(dāng) 時時a）若若 1 = 2，則主特征值則主特征值 1及相應(yīng)特征向量的求法同（及相應(yīng)特征向量的求法同（1）；）；(2)21( )limkikkixx收斂速度取決于收斂速度取決于 的程度的程度。向量向量 、113r)(1)1(kkxxc）若）若 ，則連續(xù)迭代兩次則連

3、續(xù)迭代兩次，計算出計算出x(k+1)，x(k+2)，210)()1()2( kjkjkjqxpxx)(1)1(kkxx分別為主特征值分別為主特征值 1、 2相應(yīng)的特征向量的近似值相應(yīng)的特征向量的近似值。然后對然后對j = 1, 2, , n 解方程解方程b）若若 1 = - 2，對，對i = 1, 2, , n求出求出 、 后后，由公式由公式pq2122 pqip 2222 pqip 解出主特征值解出主特征值 1、 2。此時收斂速度取決于此時收斂速度取決于 的程度的程度。113r)(2)1(kkxx)(1)1(kkxx向量向量 、 分別為相應(yīng)于分別為相應(yīng)于 1， 2的特征向量的近似值的特征向量

4、的近似值。令令max(x)表示向量表示向量x分量中絕對值最大者。即如果有某分量中絕對值最大者。即如果有某i0，使，使iniixx 1max0則則 max (x) = xi對任取初始向量對任取初始向量x(0)，記，記)max()0()0()0(xxy 則則)0()1(Ayx 一般地，若已知一般地，若已知x(k)，稱公式，稱公式 ), 1, 0()max()()1()()()(kAyxxxykkkkk定理定理 設(shè)設(shè)A Rn n具有完全特征向量系具有完全特征向量系， 1, 2, , n為為An21則對任初始向量則對任初始向量x(0)，由規(guī)范化的乘冪法公式確定的由規(guī)范化的乘冪法公式確定的1)()max

5、(lim kkx(1)（2）y(k)為相應(yīng)于主特征值為相應(yīng)于主特征值 1的特征向量近似值的特征向量近似值的的n個特征值個特征值，且滿足且滿足向量序列向量序列y(k)，x(k)滿足滿足3.1.2 反冪法反冪法若若 A 有有| 1 | | 2 | | n |，則，則 A 1 有有11111 nnA 1 的主特征根的主特征根 A的絕對值最小的特征根的絕對值最小的特征根)()1(kkxAx )(1)1(kkxAx 如何計算如何計算解線性方程組解線性方程組對應(yīng)同樣一組特征向量。對應(yīng)同樣一組特征向量。設(shè)設(shè)A Rn n可逆，則無零特征值，由可逆，則無零特征值，由)0( xxAx 有有 xxA 11 規(guī)范化反

6、冪法公式為規(guī)范化反冪法公式為 ), 1, 0()max()()1()()()(kyAxxxykkkkk為節(jié)省工作量，可先對為節(jié)省工作量，可先對A進行進行LU分解，再解三角形方程分解，再解三角形方程( )(1)(0,1,)kkLxyUxxk3.1.3 乘冪法的加速（乘冪法的加速（原點位移法）原點位移法） njjkjjkkkvAxx111)1()( ，希望，希望 | 2 / 1 | 越小越好。越小越好。取取 0（常數(shù)），用矩陣（常數(shù)），用矩陣B = A - 0I 來代替來代替A進行乘冪迭代。進行乘冪迭代。0 ii (i = 1, 2, , n)iiiiiivvAvvIABv)()(000 設(shè)設(shè) i

7、 (i = 1, 2, , n)為矩陣為矩陣B B 的特征值，則的特征值，則B與與A特征值之間特征值之間應(yīng)有關(guān)系式：應(yīng)有關(guān)系式：關(guān)于矩陣關(guān)于矩陣B的乘冪公式為的乘冪公式為 )0(0)0()()(xIAxBxkkk 為加快收斂速度，適當(dāng)選擇參數(shù)為加快收斂速度，適當(dāng)選擇參數(shù) 0，使使00210()maxjj n 達(dá)到最小值。達(dá)到最小值。 0101 1210()knjkjjjvv jkjnjjkvv12111 當(dāng)當(dāng) i (i = 1, 2, , n)為實數(shù)，且為實數(shù)，且 1 2 n時，取時，取)(212*0n 則為則為 ( 0) 的極小值點。這時的極小值點。這時122122122*01*022212

8、12121 nnnn若已知矩陣若已知矩陣A的特征值的特征值 i的一個近似值的一個近似值 0 ，要求，要求 i和對和對記記B = A - 0I，對任取初始向量對任取初始向量x(0) Rn， ), 1, 0()max()()1()()()(kyBxxxykkkkk應(yīng)的特征向量，可考慮利用原點移位加速的反冪法。應(yīng)的特征向量，可考慮利用原點移位加速的反冪法。3.2 子空間迭代法子空間迭代法斯密特斯密特（Schmidt）正交化過程：正交化過程： 設(shè)設(shè) 1， 2， 3 為為R3上的三個線性無關(guān)的向量，上的三個線性無關(guān)的向量，令令 ，則，則 1為單位長度的向量，再令為單位長度的向量，再令2111222211

9、222,),( 可以驗證可以驗證（ 1, 2）= 0，即，即 1與與 2正交。若令正交。若令22311333),(),( 則則0),(),(2313 2333 即與即與 1, 2正交，將其單位化為正交，將其單位化為于是向量組于是向量組 1, 2, 3構(gòu)成構(gòu)成R3上一組標(biāo)準(zhǔn)正交基，且上一組標(biāo)準(zhǔn)正交基，且232322131221321321),(),(),(,QR其中其中Q = 1, 2, 3為正交矩陣，為正交矩陣，R是上三角陣。是上三角陣。對對n維向量空間，設(shè)維向量空間，設(shè) 1, , n為為Rn上上n個線性無關(guān)的向量，個線性無關(guān)的向量，類似有類似有11 2111 11222),( 2222 22

10、311333),(),( 2333 jjnnjnn ),(11 2nnn 232322322113122111),(),(),(),(),(),(,nnnnnn 即即QR Q為正交陣為正交陣，R 為上三角陣為上三角陣將將n個線性無關(guān)向量變換為個線性無關(guān)向量變換為n個兩兩正交向量的方法稱為個兩兩正交向量的方法稱為 斯密特正交化方法。斯密特正交化方法。斯密特正交化過程將可逆陣斯密特正交化過程將可逆陣A分解為正交陣與上三角陣的乘積。分解為正交陣與上三角陣的乘積。3.3 雅克比雅克比 (Jacobi) 旋轉(zhuǎn)法旋轉(zhuǎn)法 1預(yù)備知識預(yù)備知識1）若若B是上（或下）三角陣或?qū)顷嚕巧希ɑ蛳拢┤顷嚮驅(qū)顷?，則

11、則B的主對角元素即是的主對角元素即是B的特征值的特征值。2）若矩陣若矩陣P滿足滿足PTP = I，則稱則稱P為正交矩陣為正交矩陣。顯然顯然PT = P-1，且且P1, P2, 是正交陣是正交陣時時，其乘積其乘積P = P1P2Pk仍為正交矩陣仍為正交矩陣。3）稱矩陣稱矩陣jiPij11cossin11sincos11為旋轉(zhuǎn)矩陣為旋轉(zhuǎn)矩陣 2雅克比方法雅克比方法設(shè)矩陣設(shè)矩陣A Rn n是對稱矩陣是對稱矩陣，記記A0 = A，對對A作一系列旋轉(zhuǎn)相似變換作一系列旋轉(zhuǎn)相似變換TCP AP其中其中A與與C的的 元素有下面的關(guān)系元素有下面的關(guān)系11c o ss i n11s i nc o s11Ti jP

12、Pij列列ij行行P是一個正交陣是一個正交陣，我們稱它是我們稱它是（i, j）平面上的旋轉(zhuǎn)矩陣平面上的旋轉(zhuǎn)矩陣 PTAP只改變只改變A的第的第i行行、j行行、i列列、j列的元素列的元素；A和和C的元素僅在第的元素僅在第i行（列）和第行（列）和第j行（列）不同行（列）不同，它們之間有如下的關(guān)系它們之間有如下的關(guān)系：cossin,cossinikkiikjkjkkjjkikccaaki jccaa22(1)22(1)cossinsin2sincossin21sin2(cos22kiiiijjijkjjiijjijijjiiijjijcaaacaaaccaaa( , )lklkcal ki j其它元

13、素不變22,1,1( )( )nnijiji ji jijCaAa引入記號22( )( )( )( )22(3.3.8)ijijCACAac易知我們選取我們選取Pij，使得使得 ，因此需使因此需使 滿足滿足0ijjicc2tg2(3.3.9)ijiijjaaa將將 限制在下列范圍內(nèi)限制在下列范圍內(nèi)44 如果如果 0iijjaa0ija 4 0ija 4 直接從三角函數(shù)關(guān)系式計算直接從三角函數(shù)關(guān)系式計算sin 和和cos ，記，記)1()1()1()1()1(2sinkijkjjkiikjikiiaaaxaay則則yx2tg當(dāng)當(dāng) 時，有下面三角恒等式：時，有下面三角恒等式：422222tg112

14、cos1cos2yxy于是于是 2221cos2yxy采用下面公式計算采用下面公式計算 sin 222cos2tgcossin22sinyxx3.3.2 Jacobi旋轉(zhuǎn)法 由于一次正交相似變換AC=PTAP可將A的兩個非對角元素化為零。因此可選一系列正交變換矩陣Pk，對A進行正交相似變換，直至將化為近似對角矩陣首先，在首先，在A=A1=(als(1)的非對角元素中選取絕對值最大元素，記為的非對角元素中選取絕對值最大元素，記為11111 1111,1211 1(2)(2)(2)2,=()0ijTlsijj iijPPAAP APAaaa對確定的，用公式(3.3.9)確定出 ，從而產(chǎn)生平面旋轉(zhuǎn)矩

15、陣約化 為且的非對角元素1111(1)(1),(0)maxijijl saa然后，再在然后，再在A2的非對角元素中選取絕對值最大元素的非對角元素中選取絕對值最大元素22222,222223222211 12,12112( ,),( ,),( , ),1,2,kkijTTTkijTTTkkkijPPAijj iAP A PP P APPPPkAPP P A PPP確定出其位置，按公式(3.3.9)確定 角，從而產(chǎn)生平面旋轉(zhuǎn)矩陣并將 中位置的元素化為零 得到如此做下去 可得到一系列平面旋轉(zhuǎn)矩陣使得(1)2( )2(1)2111111(3.3.8)()()2()2()()2()22()()()1()

16、(1)(1)21()0 ()(1)kkkkkijijkijkkkkkAAaaAaAAAAn nn nAkn n 由關(guān)系式知并有P78 例 Jacobi 過關(guān)法 Jacobi方法中 ,每變換一次都要在非對角元素中掃描選取絕對值最大元素,這難免增加很大的工作量.為減少工作量和提高運行效率,下面介紹一種改進的方法,即Jacobi過關(guān)法 具體做法:0112131232(1)111212( ),( ),.nnnnijijrrAAAnaaaaaaaanAn計算矩陣 的非對角元素平方和預(yù)取閥值對 的非對角元素進行掃描對的元素進行旋轉(zhuǎn)點還 將化為零 其余元素視為過關(guān) 不作變換 當(dāng)所有的非對角元

17、素絕對值都小于 后,縮小閥值,如取重復(fù)前面的步驟.按此方法,閥值不斷縮小 直到滿足之后停止計算這種方法稱為Jacobi過關(guān)法.3.4 Household方法 Householder方法是計算實對稱矩陣A的部分或全部特征值及其特征向量,計算過程是:先利用正交相似變換將A約化為對稱的三對角矩陣C,其次應(yīng)用對分法計算C的特征值,最后計算特征向量.3.4.1 實對稱矩陣的三對角化22(, ).cossin0 (, )1/,sinsjkkjjkkjjkikikkjjkikccki jccaaki jaaacoa 選取 使得某個變?yōu)榱阒灰? , , ,2,3,12,4,12, ,13,4,23,5,23

18、, ,21, ,2, ,1,0(, ),;,;,2,3,1,1,2,.Ti j ki j ki j kjkkjnnnn ni j iPPCPAPccki jPPPPPPPinPjiinAAC 記上述旋轉(zhuǎn)矩陣為則相似變換使元素取或?qū)τ靡淮螌?進行正交相似是變換,便可將 化為三對角矩陣Household變換 設(shè) 為單位向量,即 nRu1uuT2(3.4.2)THIuu稱此為Household矩陣，簡稱H矩陣。容易驗證 H 矩陣是對稱的且正交的矩陣，即IHHHHHTT2,由于 uuuuIHuT)2(且對如何與u 直交的向量v, 都有 0Tu v vvuuIHvT2因此，對任意 ，可設(shè) ，則其H變換為

19、nRxvuxvuHvHuHx若將 中所用與向量u正交的方向視為一個鏡面，有上述公式看到，H 變換不改變向量在鏡面上的投影，并將向量沿法向量的投影改變?yōu)榉捶较虻乳L度的向量，因此H變換也稱為鏡面反射變換nR由H變換的性質(zhì)不難知道，對任意非零向量 ，如果 則必存在則必存在H矩陣，使得矩陣，使得 nRyx,22yxyHx 事實上，當(dāng)取 時，即可驗證由（3.4.2）式所定義的矩陣滿足 的要求 2yxyxuyHx 121/221211,()(,).(1),(,() ,0,0) ,TnnTrrii rHAaa aanraca aasign assarn 利用變換將 化為三對角矩陣的過程 可逐列 行 的方式進行.設(shè)為矩陣的某一列向量如想將此向量的后個分量化為零 即將 變成其中1221122112,(0,0,() ,) ,21()rrnTrrrnrra aaacacHHacwacasign as aas ssign aaw假設(shè)不全為零由于且故存在矩陣 滿足若令11200() ,) ,Tn rTrrrnHIHIwwwasign as aa 則 矩陣可取為其中P83 例3.5 3.4.2 求對稱三對角矩陣特征值的對分法1122210,1,2,1nnniabbabAbbabin 考慮對稱的三對角矩陣這里約定00112112( ),( )1.(