数值分析1030求矩阵全部特征值 QR 方法

1、数值分析数值分析数值分析数值分析 第三节第三节 求矩阵全部特征值的求矩阵全部特征值的QR方法方法一、求矩阵全部特征值的一、求矩阵全部特征值的QR方法方法 6060年代出现的年代出现的QRQR算法是目前计算中小型矩阵的算法是目前计算中小型矩阵的全部特征值与特征向量的最有效方法。全部特征值与特征向量的最有效方法。 理论依据：理论依据：任一非奇异实矩阵都可分解成一个正交矩阵任一非奇异实矩阵都可分解成一个正交矩阵Q Q和一个和一个上三角矩阵上三角矩阵R R的乘积，而且当的乘积，而且当R R的对角元符号取定时，的对角元符号取定时，分解是唯一的。分解是唯一的。 11 QRQR (1,2,). kkkkkk

2、AQ RkAR QAAA 方方法法的的基基本本思思想想是是利利用用矩矩阵阵的的分分解解通通过过迭迭代代格格式式将将化化成成相相似似的的上上三三角角阵阵（或或分分块块上上三三角角阵阵），从从而而求求出出矩矩阵阵的的全全部部特特征征值值与与特特征征向向量量。数值分析数值分析数值分析数值分析111111121112, (2 ,3, ) kAAQ RQARAR QQAQAAAAk 由由即即。于于是是即即与与相相似似。同同理理可可得得，。故故它它们们有有相相同同的的特特征征值值。 可证，在一定条件下，基本可证，在一定条件下，基本QRQR方法产生的矩方法产生的矩阵序列阵序列A Ak k “ “基本基本”收

3、敛于一个上三角阵（或收敛于一个上三角阵（或分块上三角阵）。即主对角线（或主对角线子块）分块上三角阵）。即主对角线（或主对角线子块）及其以下元素均收敛，主对角线（或主对角线子及其以下元素均收敛，主对角线（或主对角线子块）以上元素可以不收敛。特别的，如果块）以上元素可以不收敛。特别的，如果A A是实对是实对称阵，则称阵，则A Ak k “ “基本基本”收敛于对角矩阵。收敛于对角矩阵。数值分析数值分析数值分析数值分析QR方法的实际计算步骤方法的实际计算步骤HouseholderAHessenbergB 用用阵阵作作正正交交相相似似变变换换上上第第阵阵一一步步.*:* 1kkkGivenkkkBQ R

4、BBR Q 用用变变换换产产生生迭迭代代序序列列第第二二步步12*n 数值分析数值分析数值分析数值分析HouseholderAB 用阵用阵作正交相似变换作正交相似变换（对称阵）三对角阵（对称阵）三对角阵* 数值分析数值分析数值分析数值分析二、化一般矩阵为上二、化一般矩阵为上Hessenberg阵阵111211121222123233311 (2,3,), Househ old e rnnnnnnnnniihhhhhhhhhhhHhhhinA 称称形形如如 的的矩矩阵阵为为上上海海森森堡堡（H H e es ss se en nb be er rg g) )阵阵。如如果果此此对对角角线线元元全全

5、不不为为零零 则则称称该该矩矩阵阵为为不不可可约约的的上上H H e es ss se en nb be er rg g矩矩阵阵。讨讨论论用用变变换换将将一一般般矩矩阵阵相相似似变变换换成成H H e es ss se en nb be er rg g阵阵数值分析数值分析数值分析数值分析11111111 ,1000 01 HouseholderHHH AHHHHHnHouseholder 首首先先，选选取取矩矩阵阵使使得得经经相相似似变变换换后后的的矩矩阵阵的的第第一一列列中中有有尽尽可可能能多多的的零零元元素素。为为此此，应应取取为为如如下下形形式式其其中中为为阶阶矩矩阵阵。11121111

6、1221121311212131(,) ,(,) ,TTTnnaa HH AHH aH A Haaaaaaaa 于于是是有有 其其中中数值分析数值分析数值分析数值分析222222111111.(, 0) 0,2nnnnTaaAaaHH aH AHn 只只要要取取使使得得就就会会使使得得变变换换后后的的矩矩阵阵的的第第一一列列出出现现个个零零元元。数值分析数值分析数值分析数值分析2211221222211221000*0100*00*0022, .nnnHouseholderHH H AH HHnnHouseholderHHHHH H AH HHHHHessenberg 同同理理，可可构构造造如

7、如下下列列形形式式矩矩阵阵使使得得* *如如此此进进行行次次，可可以以构构造造个个矩矩阵阵使使得得其其中中为为上上矩矩阵阵AH。特特别别地地，当当 为为实实对对称称矩矩阵阵，则则经经过过上上述述正正交交变变换换后后，变变为为三三对对角角阵阵。数值分析数值分析数值分析数值分析1221522 23 2105 2 22 2 021002412,022, (2,2)2(1,0)(22,2) :TTTHouseholderAAHouseholderHHu 用用变变换换将将矩矩阵阵 化化成成上上H H e es ss se e例例解解n nb be er rg g阵阵。求求矩矩阵阵满满：足足，数值分析数值

8、分析数值分析数值分析22 22222210442220122222442210000100220022220022TTuuHIu uH 数值分析数值分析数值分析数值分析2112 100052223201001052 22 222002202102202410022100052510100103222 000223220012220022HH H AH H 于于 是是 有有数值分析数值分析数值分析数值分析用用Household方法对矩阵方法对矩阵A作正交相似变换作正交相似变换, 使使A相似与上相似与上Hessenberg阵，算法如下：阵，算法如下：(1)(1)111221111111(1)112

9、,1,2,.,21(1)()()() ) ,()(0,.,0,.,)kkTkknkkkiki kkkkkkkkkkkknkknHIUUsign aaaUaaa ，计算计算1(2)kHAA计计算算 (1)11(1),1,11()21,nkjlljlkkkijijjijk kntuaiknaat u （ ）（ ）数值分析数值分析数值分析数值分析(1)11(1)1,.,1(1)(2)1,.,nkiilll kkkijijijinta ujknaat u 1(3)kAHA计计算算 数值分析数值分析数值分析数值分析三、上三、上Hessenberg阵的阵的QR分解分解对对上上Hessenberg阵阵只需要

10、将其次对角线上的只需要将其次对角线上的元素约化为零，用元素约化为零，用Given变换比用变换比用Householder变变换更节省计算量。为此先介绍换更节省计算量。为此先介绍Given变换。变换。数值分析数值分析数值分析数值分析,11cossin11sincos11 ()i jiRjjin n阶阶方方阵阵称称为为平平面面旋旋转转阵阵，或或称称为为G Gi iv ve en ns s变变换换阵阵。定定义义 、平面旋转阵、平面旋转阵(Givens(Givens变换阵变换阵) )数值分析数值分析数值分析数值分析1,1,TTi ji ji ji jRRIRRRi i, ,j j平平面面旋旋转转阵阵的的

11、（ ）平平面面旋旋转转阵阵是是非非对对称称交交质质：阵阵性性的的正正。,2Ti ji jRR（ ）也也是是平平面面旋旋转转阵阵。( (3 3) )d de et t( () )= =1 1数值分析数值分析数值分析数值分析1,.,Rxxxxi i, ,j jT T1 12 2n n平平面面旋旋转转阵阵的的作作用用：（ ）将将向向量量 = =的的第第j j个个分分量量约约化化为为零零。,cossinsincos1,., ;,i jiijjijkkyRxyxxyxxyxkn ki j, ,若若令令，有有 111,222121212cossinsincoscossinsincosyxRyxxxxxxx

12、 jy调调整整 ，可可将将 约约化化为为零零。0tanjjixyx令令，得得,.,i jRxx xxxT T12n12n左左乘乘向向量量 = =只只改改变变 的的第第i i个个分分量量和和第第j j个个分分量量。jxix 数值分析数值分析数值分析数值分析0tanjjixyx令令，得得2222cossiniiijjjijxxCrxxxxSrxx 所所以以，取取22,0iijijjyCxSxrxxy于于是是 ,., ,.,0,.,.i jRxxxr xxxxT T1 1i i- -1 1i i+ +1 1j j- -1 1j j+ +1 1n n= =jxix 数值分析数值分析数值分析数值分析2,

13、.,xxxxT T1 12 2n n（ ）将将向向量量 = =的的第第i i+ +1 1个个分分量量到到第第n n个个分分量量约约化化为为零零。22,11,., ,0,.,i iiiRxxxrxxrxxT T1 1i i- -1 1i i+ +2 2n n= =,2,122212,., ,0,0,.,i ii iiiiRRxxxrxxrxxxT T1 1i i- -1 1i i+ +3 3n n= =,2,122,., ,0,.,0,i ni ii iinRRRxxxrrxxT T1 1i i- -1 1= =,(3)i ji jiiji jjRAARAARRARAAT TT T左左乘乘 只只

14、改改变变 的的第第i i，j j行行。右右乘乘 只只改改用用对对矩矩阵阵 作作变变换换得得变变 的的第第i i，j j列列。只只改改变变 的的第第i i，j j行行和和到到的的结结第第论论i i，j j列列。数值分析数值分析数值分析数值分析2,1,4,0,0.xxrT TT T已已知知向向量量 = =，试试用用G Gi iv ve en ns s变变换换将将 约约化化为为例例(1)(1)2,1,4x xxT T：记记 = =，对对计计解解算算C C和和S S。122222121221,55xxCSxxxx1,2(1)(2)1,221055120550015,0,4TRRxx 数值分析数值分析数

15、值分析数值分析(2)4,21xS 5 5对对计计算算C C和和S,C=S,C=2121 (1)1,31,34021010,21,0,04021TRRx 5 52 21 15 52 21 1(2)5,0,4Tx数值分析数值分析数值分析数值分析、用、用 GivensGivens变换对变换对上上Hessenberg阵作阵作QR分解分解(1)(1)(1)11121(1)(1)(1)21222(1)(1)1 1 nnnnnnbbbbbbBbbnGivensBQR 对对上上H H e es ss se en nb be er rg g阵阵, ,通通常常用用个个变变换换阵阵可可将将它它化化成成上上三三角角矩

16、矩阵阵，从从而而得得到到 的的分分解解式式。数值分析数值分析数值分析数值分析(1)211111(2)(2)(2)112131(2)(2)(2)22232(2)(2)(2)1232333(2)(2)1 0(cossin00sincos00(1,2)0011(1,2) nnnnnnnbRrbbbbbbRBBbbbbb 具具体体步步骤骤为为：设设否否则则进进行行下下一一步步），取取旋旋转转矩矩阵阵则则(1)(1)(1)(1)1121111112111 cos, sin, .bbrbbrr 其其中中数值分析数值分析数值分析数值分析2322222( 3 )( 3 )( 3 )( 3 )11213111(

17、 3 )( 3 )( 3 )223212( 3 )( 3 )( 3 )333132( 3 )( 3 )4341 0(10cossinsincos (3, 2)11 (3 , 2)nnnnnnnbRrbbbbrbbbbbbRBbb （）设设否否 则则 进进 行行 下下 一一 步步 ） ， 再再 取取 旋旋 转转 矩矩 阵阵 则则3( 3 )4( 3 )( 3 )1( 2 )( 2 )( 2 )2( 2 )23222222223222 cos, sin, ()() .nnnnnBbbbbbrbbrr 其其 中中数值分析数值分析数值分析数值分析1( )( )( )( )1111111( )( )(

18、)11111( )( )( )1( )( )( )1111( )( )1 (1, ) kkkkkkkknnkkkkkkknknkkkkkknknkkkkkknknkknnnnBR kk Brbbbbrbbbbhhbbbbb 1k 假假设设上上述述过过程程已已进进行行了了步步，有有数值分析数值分析数值分析数值分析()1()()1()2()21 0,11 (1,)cossinsincos1 cos, sin, ()() .kkkkkkkkkkkkkkkkkkkkkkkkbR kkbbrrrbb 设设取取其其中中数值分析数值分析数值分析数值分析(1)(1)(1)11111(1)(1)1(1)(1)(

19、1)111111(1)(1)(1)21212(1)(1)1 (1, )1kkkkknkkkkkknkkkkkkkknknkkkkkknknkknnnnrbbbrbbR kk BBbhbbhhbbn 于于是是因因此此，最最多多做做次次旋旋转转变变换换，即即( )( )( )( )112131( )( )2232( )33 ( ,1) (2,1)(1,2)nnnnnnnnnnnHR n nR nnRBrbbbrbbRrbr 得得数值分析数值分析数值分析数值分析213212132123( ,1),(2,3, ) 4,() TTTnnTTTnnR i iinHR RRRQRQR RRnQRO nHRQ

20、QRQR 因因为为均均为为正正交交矩矩阵阵，故故其其中中仍仍为为正正交交矩矩阵阵。可可算算出出完完成成这这一一过过程程的的运运算算量量约约为为比比一一般般矩矩阵阵的的分分解解的的运运算算量量少少一一个个数数量量级级。可可证证明明仍仍是是上上H H e es ss se en nb be er rg g阵阵，于于是是可可按按上上述述步步骤骤一一直直迭迭代代下下去去，这这样样得得到到的的方方法法的的运运算算量量比比基基本本QR方方法法大大为为减减少少。需需要要说说明明的的是是，通通常常用用方方法法计计算算特特征征值值，然然后后用用反反幂幂法法求求其其相相应应的的特特征征向向量量。数值分析数值分析数

21、值分析数值分析22532644445 (6,4)64 (1,0)(652,4)10 2010.9160250.2773500.8320500.55470020.2773500.08397470.5547000.832050TTTTTQRAAuuuIu u 用用方方法法求求矩矩阵阵的的全全部部特特征征值值。首首先先将将 化化成成上上H H e es ss se en nb be e例例：r rg g：阵阵，取取解解110000.8320500.55470000.5547000.832050H 于于是是 数值分析数值分析数值分析数值分析11221111151.3867503.3282007.211

22、1021.2307688.15384000.1538462.230767, 5( 7.21102)8.774964 cos50.56980. sin0.821781 HH AHHAHQRBHrr 即即为为与与 相相似似的的上上H H e es ss se en nb be er rg g阵阵。将将进进行行分分解解，记记取取0.5698030.8217810(2,1)0.8217810.5698030001R 数值分析数值分析数值分析数值分析122222228.7749641.8015968.597089(2,1)00.4383101.91103000.1538462.230767 (0.438

23、310)( 0.153846)0.464526, cos0.4383100.943564, sin0.1538460.331189 RBrrr 于于 是是再再 取取11100 (3, 2)00.9435640.33118900.3311890.943564 (3, 2)(2,1)8.7749641.8015968.59708900.4645262.541982001.471953RRRBR 于于 是是数值分析数值分析数值分析数值分析12110.5698030.7754030.272165(2,1)(3,2)0.8217810.5376430.18871200.3311890.9435643.5

24、194824.92549110.8401170.3817391.0916272.31065300.4874951.388883,11TTQRRBRQ 第第一一次次迭迭代代得得重重复复上上述述过过程程 迭迭代代次次121232.9920321.0003853 12.013392 0.0074962.0046951.94197100.0003250.9998952.992032,2.004695,0.999895 3,2,1.0.007496BQR 得得精精确确值值下下三三角角非非对对角角元元的的最最大大模模为为。方方法法“基基本本”收收敛敛较较慢慢。数值分析数值分析数值分析数值分析算法：用算法：用Givens方法对上方法对上Hessenberg阵阵A作作 正交分解正交分解,A=QR。,1,22,1,1,1,1,1,1,1,1,1Hessenberg,1,2,.,1(1),(2)1,2.,(3)1,2.,i iiii iiii ii ki kikiki kikTi ik ik ik ik ik ik iA QIinaaraacsrrRAAknacasaasacaQRQknqcqsqqsacqA 输入上阵输入上阵计算计算计算计算输出 为输出 为,Q上三角阵为正交阵。上三角阵为正交阵。数值分析数值分析数值分析数值分析算法：用算法：用Givens方法对上方法对上Hessenberg阵