第8章矩阵特征值计算-ppt课件

1、上页上页下页下页第第8章章 矩阵特征问题的计算矩阵特征问题的计算 8.1 引言引言 8.2 幂法及反幂法幂法及反幂法 8.3 豪斯霍尔德方法豪斯霍尔德方法 8.4 QR方法方法上页上页下页下页8.1 引引 言言 工程技术中有多种振动问题，如桥梁或建筑物的工程技术中有多种振动问题，如桥梁或建筑物的振动，机械零件、飞机机翼的振动，及一些稳定性分振动，机械零件、飞机机翼的振动，及一些稳定性分析和相关分析在数学上都可转化为求矩阵特征值与特析和相关分析在数学上都可转化为求矩阵特征值与特征向量的问题征向量的问题. 下面先复习一些矩阵的特征值和特征向量的根底下面先复习一些矩阵的特征值和特征向量的根底知识知识

2、.上页上页下页下页 定义定义1 知知n阶矩阵阶矩阵A=(aij)，那么，那么)2()(det)det()(12211212222111211的项的项次数次数 naaaaaaaaaaaaAInnnnnnnnnn 称为称为A的特征多项式的特征多项式. 普通有普通有n个根个根(实的或复的，复根按重数计算实的或复的，复根按重数计算)称为称为A的特征值的特征值. 用用(A)表示表示A的一切特征值的集合的一切特征值的集合. A的特征方程的特征方程)1 . 1(0)det()( AI 上页上页下页下页 设设为为A的特征值，相应的齐次方程组的特征值，相应的齐次方程组 注：当注：当A为实矩阵时，为实矩阵时， (

3、)=0为实系数为实系数n次代数次代数方程，其复根是共轭成对出现方程，其复根是共轭成对出现.)2 . 1(0)( xAI 的非零解的非零解x称为矩阵称为矩阵A的对应于的对应于的特征向量的特征向量. 210131012A 例例1 求求A的特征值及特征向量，其中的特征值及特征向量，其中上页上页下页下页. 0)4)(2)(1(8147210131012)det()(23 AI 解解 矩阵矩阵A的特征方程为的特征方程为求得矩阵求得矩阵A的特征值为：的特征值为：. 4, 2, 1 对应于各特征值矩阵对应于各特征值矩阵A的特征向量分别为：的特征向量分别为：.121,101,111321 xxx上页上页下页下

4、页 定理定理1 设设为为ARnn的特征值的特征值, 且且Ax=x (x0)，那么有，那么有 -p为为A-pI的特征值，即的特征值，即(A-pI)x=(-p)x ； c为的为的cA特征值特征值(c0为常数为常数)； 下面表达有关特征值的一些结论：下面表达有关特征值的一些结论： k为为Ak的特征值，即的特征值，即Akx=kx ； 设设A为非奇特矩阵，那么为非奇特矩阵，那么0 , 且且-1为为A-1的的特征值，即特征值，即A-1x=-1x .上页上页下页下页 定理定理2 设设i(i=1,2,n)为为n阶矩阵阶矩阵A=(aij)的特征的特征值，那么有值，那么有)(Atraniiinii 11 称为称为

5、A的迹；的迹； .nA21 定理定理3 设设ARnn，那么有，那么有. )()(AAT 定理定理4 设设A 为分块上三角矩阵，即为分块上三角矩阵，即, mmmmAAAAAAA22211211其中每个对角块其中每个对角块Aii均为方阵，那么均为方阵，那么. )()(iiniAA1 上页上页下页下页 定理定理5 设设A与与B为类似矩阵即存在非奇特矩阵为类似矩阵即存在非奇特矩阵P使使B=P-1AP，那么，那么 定理定理5阐明，一个矩阵阐明，一个矩阵A经过类似变换，其特征经过类似变换，其特征值不变值不变. 一个亏损矩阵是一个没有足够特征向量的矩阵，一个亏损矩阵是一个没有足够特征向量的矩阵，亏损矩阵在实

6、际上和计算上都存在困难亏损矩阵在实际上和计算上都存在困难. A与与B有一样的特征值；有一样的特征值； 假设假设y是是B的特征向量，那么的特征向量，那么Py是是A的特征向的特征向量量. 定义定义2 假照实矩阵假照实矩阵A有一个重数为有一个重数为k的特征值的特征值, 且对应于且对应于的的A的线性无关的特征向量个数的线性无关的特征向量个数|2|n|，那么对任何非零向量那么对任何非零向量v0(a10)，幂法的算式成立，幂法的算式成立.又设又设A有有n个线性无关的特征向量，个线性无关的特征向量，1对应的对应的r个线个线性无关的特征向量为性无关的特征向量为x1,x2,xr，那么由，那么由(2.2)式有式有

7、 假设假设A的主特征值为实的重根的主特征值为实的重根, 即即1=2=r, 且且 |r|r+1|n|，,)/(11110 nriikiiriiikkkxaxavAv 上页上页下页下页).0(lim111 riiiriiikkkxaxav设设 为为A的特征向量，这阐明当的特征向量，这阐明当A的主特征值是实的重根的主特征值是实的重根时，定理时，定理5的结论还是正确的的结论还是正确的. 运用幂法计算运用幂法计算A的主特征值的主特征值1及其对应的特征向及其对应的特征向量时，假设量时，假设|1|1(或或|1|2|n|，那么对恣意，那么对恣意非零初始向量非零初始向量v0=u0(a10)，有幂法计算公式为，有

8、幂法计算公式为那么有那么有 ,)max(lim11xxukk .lim1 kk上页上页下页下页例例1 用幂法计算矩阵用幂法计算矩阵 1634310232A的主特征值与其对应的特征向量的主特征值与其对应的特征向量.解取解取 v0=u0=(0,0,1)T , 那么那么 TTvuv25. 0 , 1, 5 . 01, 4,1 , 4 , 211111 ), 2 , 1(max1 k/vuvAuvkkkkkkk 上页上页下页下页直到直到k=8 时的计算结果见下表时的计算结果见下表k1 2, 4, 1, 4 0.5, 1, 0.252 4.5, 9, 7.75 90.5, 1, 0.86113 5.72

9、22, 11.4444, 8.36111.44440.5, 1, 0.73604 5.4621, 10.9223, 8.2306 10.92230.5, 1, 0.75365 5.5075, 11.0142, 8.2576 11.01420.5, 1, 0.74946 5.4987, 10.9974, 8.2494 10.99740.5, 1, 0.75017 5.5002, 11.0005, 8.2501 11.00050.5, 1, 0.75008 5.5000, 11.0000, 8.2500 11.00000.5, 1, 0.7500TkvTku Tx7500. 0,0 . 1,5 .

10、 0,0000.1111 从而从而k 见书见书p303-例例3.上页上页下页下页8.2.2 幂法的加速方法幂法的加速方法1、原点平移法、原点平移法 由前面讨论知道，运用幂法计算A的主特征值的收敛速度主要由比值 r=|2/1|来决议，但当r 接近于1时，收敛能够很慢. 这时，一个补救方法是采用加速收敛的方法.其中其中p为参数，设为参数，设A的特征值为的特征值为i，那么对矩阵，那么对矩阵B的的特征值为特征值为i-p ，而且，而且A, B的特征向量一样的特征向量一样. 引进矩阵引进矩阵 B=A-pI .上页上页下页下页 假设要计算假设要计算A的主特征值的主特征值1, 只需选择适宜的数只需选择适宜的数

11、p，使，使1-p为矩阵为矩阵B=A-pI 的主特征值，且的主特征值，且 1212max ppini那么，对矩阵那么，对矩阵B=A-pI运用幂法求其主特征值运用幂法求其主特征值1-p, 收敛速度将会加快收敛速度将会加快. 这种经过求这种经过求B=A-pI的主特征值和的主特征值和特征向量，而得到特征向量，而得到A的主特征值和特征向量的方法叫的主特征值和特征向量的方法叫原点平移法原点平移法. 对于对于A的特征值的某种分布，它是非常的特征值的某种分布，它是非常有效的有效的.上页上页下页下页例例4 设设AR44有特征值有特征值),4 , 3 , 2 , 1(15 jji 比值比值r=|2/1|0.9.

12、做变换做变换B=A-12I (p=12),那么那么B的特征值为的特征值为. 1, 0, 1, 24321 运用幂法计算运用幂法计算B的主特征值的主特征值1的收敛速度的比值为的收敛速度的比值为. 9 . 021121212 pp 虽然经常可以选择有利的虽然经常可以选择有利的p值值, 使幂法得到加速使幂法得到加速, 但设计一个自动选择适当参数但设计一个自动选择适当参数p的过程是困难的的过程是困难的.上页上页下页下页 下面思索当下面思索当A的特征值是实数时，怎样选择的特征值是实数时，怎样选择p使使采用幂法计算采用幂法计算1得到加速得到加速.ppn 1且使收敛速度的比值且使收敛速度的比值.min,ma

13、x112 ppppn 设设A的特征值都是实数，且满足的特征值都是实数，且满足)10. 2(,121nn 那么对实数那么对实数p，使矩阵，使矩阵A-pI的主特征值为的主特征值为1-p或或n-p时，当我们计算时，当我们计算1及及x1时，首先应选取时，首先应选取p使使上页上页下页下页显然，当显然，当2-p=-(n-p )时，即时，即 P=(2+n)/2P* 时时为最小值，这时收敛速度的比值为为最小值，这时收敛速度的比值为.2212112nnnpppp 当希望计算当希望计算n时，应选取时，应选取 p=(1+n-1)/2P* 使得运用幂法计算使得运用幂法计算n得到加速得到加速. 当当A的特征值都是实数，

14、满足的特征值都是实数，满足且且2, n能初步估计出来，我们就能确定能初步估计出来，我们就能确定P*的近似的近似值值.nn 121上页上页下页下页 例例2 用原点平移加速法求例用原点平移加速法求例1中矩阵中矩阵A的主特征值的主特征值与其对应的特征向量与其对应的特征向量.5 . 36345 . 510235 . 4 B,1634310232 A对对B运用幂法，仍运用幂法，仍取取 v0=(0,0,1)T , 那么那么 .875. 0 , 1, 5 . 01, 4,5 . 3 , 4 , 211111TTvuu ), 2 , 1(max1 k/vuvBuvkkkkkkk 解解 取取p=-2.5, 做平

15、移变换做平移变换B=A-pI，那么，那么上页上页下页下页迭代迭代5步的计算结果见下表步的计算结果见下表k1 2, 4, 3.54 0.5, 1, 0.8752 7, 14, 10.5625 140.5, 1, 0.75453 6.76, 13.5179, 10.1406 13.51790.5, 1, 0.75074 6.7503, 13.5007, 10.1256 13.50070.5, 1, 0.75005 6.7500, 13.5000, 10.1250 13.50000.5, 1, 0.7500TkuTkv可得到可得到B的主特征值为的主特征值为 113.5000, 主特征向量为主特征向量

16、为 v1 (0.5 ,1.0, 0.7500)T ,因此，因此，A的主特征值为的主特征值为 1 = 1 +p 11.0000, 主特征向量仍为主特征向量仍为 x1 =(0.5,1,0.7500)T .k 上页上页下页下页 原点位移的加速方法，是一个矩阵变换方法原点位移的加速方法，是一个矩阵变换方法. 这这种变换容易计算，又不破坏矩阵种变换容易计算，又不破坏矩阵A的稀疏性，但的稀疏性，但p的的选择依赖对选择依赖对A的特征值分布的大致了解的特征值分布的大致了解. 见书见书p306-例例5.上页上页下页下页 设设ARnn为对称矩阵，称为对称矩阵，称.),(),()(xxxAxxR 为向量为向量x的瑞

17、利商，其中的瑞利商，其中(x, x)=xTx为内积为内积. 由定理由定理11知道，实对称矩阵知道，实对称矩阵A的特征值的特征值1及及n可用瑞利可用瑞利商的极限值表示商的极限值表示. 下面我们将瑞利商运用到用幂法下面我们将瑞利商运用到用幂法计算实对称矩阵计算实对称矩阵A的主特征值的加速上来的主特征值的加速上来.2、瑞利商、瑞利商(Rayleigh)加速加速上页上页下页下页 定理定理14 设设ARnn为对称矩阵，特征值满足为对称矩阵，特征值满足对应的特征向量对应的特征向量xi满足满足(xi, xj)=ij (单位正交向量单位正交向量) ，运用幂法公式运用幂法公式(2.9)计算计算A的主特征值的主特

18、征值1，那么规范，那么规范化向量化向量uk的瑞利商给出的瑞利商给出1的较好的近似值为的较好的近似值为,121nn kkkkkkuuuAuuR2121, 由此可见，由此可见，R(uk)比比k更快的收敛于更快的收敛于1.上页上页下页下页 证明证明 由由(2.8)式及式及,)max(,)max(00100uAuAAuvuAuAukkkkkkk )11. 2(.),(),(),(),()(2121122112200001 knjkjjnjkjjkkkkkkkkkaauAuAuAuAuuuAuuR 得得上页上页下页下页 幂法的瑞利商加速迭代公式可以写为幂法的瑞利商加速迭代公式可以写为 kkkkkkkkk

19、kvukuuuvAuv /), 2 , 1(),(),(1111其中其中A为为n阶实对称矩阵阶实对称矩阵.,11kkux 对给定的误差限对给定的误差限，当，当| kk-1|时，取近时，取近似值似值上页上页下页下页8.2.3 反幂法反幂法 反幂法是用于求非奇特矩阵A的按模最小的特征值和对应特征向量的方法. 而结合原点平移法的反幂法那么可以求矩阵A的任何一个具有先了解的特征值和对应的特征向量。设矩阵设矩阵A非奇特非奇特,其特征值其特征值i (i=1,2,n) ,满足满足0121 nn 其相应的特征向量其相应的特征向量x1,x2,xn线性无关，那么线性无关，那么 A-1 的的特征值为特征值为1/ i

20、 ，对应的特征向量仍为，对应的特征向量仍为xi (i=1,2,n).iiiiiixxAxAx11 上页上页下页下页此时此时, A-1的特征值满足的特征值满足11111 nn因此因此, 对对A-1运用幂法运用幂法,可求出其主特征值可求出其主特征值1/n k 和特征向量和特征向量 xn uk .从而求得从而求得A的按模最小的按模最小特征值特征值 n 1/k 和对应的特征向量和对应的特征向量 xn uk ，这种求这种求A-1的方法就称的方法就称为反幂法为反幂法.上页上页下页下页为了防止求为了防止求A-1, 可经过解线性方程组可经过解线性方程组Avk=uk-1得到得到vk，采用采用LU分解法，即先对分

21、解法，即先对A进展进展LU分解分解A=LU, 此时反此时反幂法的迭代公式为幂法的迭代公式为 , 2 , 1/max,1 kvuvvzUvzuLzkkkkkkkkkkk 求求出出解解求求出出解解 ), 2 , 1(max11 k/vuvuAvkkkkkkk knknux ,1 反幂法的迭代公式为反幂法的迭代公式为上页上页下页下页对给定的误差对给定的误差，当，当|kk-1|n|0，那么对恣意非零初始向量那么对恣意非零初始向量u0(an0) ，由反幂法计算，由反幂法计算公式构造的向量序列公式构造的向量序列vk，uk满足满足 ,)max(limnnkkxxu .1)max(limnkkv 上页上页下页

22、下页 在反幂法中也可以用原点平移法加速迭代过程，在反幂法中也可以用原点平移法加速迭代过程，或求其它特征值与其对应的特征向量或求其它特征值与其对应的特征向量. 假设矩阵假设矩阵(A-pI)-1存在，显然其特征值为存在，显然其特征值为,1,1,121pppn 对应的特征向量依然是对应的特征向量依然是x1,x2,xn，现对矩阵，现对矩阵(A-pI)-1运用幂法，得到反幂法的迭代公式运用幂法，得到反幂法的迭代公式)12. 2()., 2 , 1()max(/.,)(, 01100 kvuupIAvvukkkkk 初初始始向向量量上页上页下页下页 假设假设p是是A的特征值的特征值j的一个近似值，且设的一

23、个近似值，且设j与与其它特征值是分别的，即其它特征值是分别的，即),(jippij 就是说就是说1/(j-p)是矩阵是矩阵 (A-pI)-1的主特征值，可用反幂的主特征值，可用反幂法法(2.12)计算特征值及特征向量计算特征值及特征向量. 设设ARnn有有 n个线性无关的特征向量个线性无关的特征向量 x1,x2, xn，那么，那么),0(10 jniiiaxau,)max()(0)1(0upIAupIAvkkk 上页上页下页下页,)max()(00upIAupIAukkk 其中其中.)()(10 niikiikxpaupIA 同理可得：同理可得：上页上页下页下页 定理定理16 设设ARnn有有

24、n个线性无关的特征向量，个线性无关的特征向量， 矩阵矩阵A的特征值及对应的特征向量分别记为的特征值及对应的特征向量分别记为i 及及xi (i=1,2,n)，而，而p为为j的近似值，的近似值，(A-pI)-1存在，且存在，且 ,)max(limjjkkxxu jkjkkvppv )max(1,1)max(lim即即那么对恣意非零初始向量那么对恣意非零初始向量u0(aj0) ，由反幂法计算，由反幂法计算公式公式(2.12)构造的向量序列构造的向量序列vk，uk满足满足. )( |jippij . |min/ |pprijij 且收敛速度为且收敛速度为上页上页下页下页 由该定理知，对由该定理知，对A

25、-pI(其中其中pj)运用反幂法，可运用反幂法，可用来计算特征向量用来计算特征向量xj，只需选择，只需选择p是是j的一个较好的的一个较好的近似且特征值分别情况较好，普通近似且特征值分别情况较好，普通r很小，经常只需很小，经常只需迭代一二次就可完成特征向量的计算迭代一二次就可完成特征向量的计算.反幂法迭代公式中的反幂法迭代公式中的vk是经过解方程组是经过解方程组.)(1 kkuvpIA求得的求得的, 为了节省任务量为了节省任务量, 可以先将可以先将A-pI进展三角分解进展三角分解.)(LUpIAP 于是求于是求vk相对于解两个三角形方程组相对于解两个三角形方程组.,1kkkkyUvPuLy 上页

26、上页下页下页实验阐明实验阐明, 按下述方法选择按下述方法选择u0是较好的是较好的: 选选u0使使)13. 2()1 , 1 , 1(011 PuLUv用回代求解用回代求解(2.13)即得即得v1,然后再按公式然后再按公式(2.12)进展迭代进展迭代.反幂法计算公式见书反幂法计算公式见书p311.前面已提到前面已提到.见书见书p311-例例6.上页上页下页下页8.3 豪斯霍尔德方法豪斯霍尔德方法 (1)用初等反射矩阵作正交类似变换约化普通用初等反射矩阵作正交类似变换约化普通实矩阵实矩阵A为上海森伯格矩阵为上海森伯格矩阵.8.3.1 引言引言 本节讨论两个问题本节讨论两个问题 (2)用初等反射矩阵

27、作正交类似变换约化对称用初等反射矩阵作正交类似变换约化对称矩阵矩阵A为对称三对角矩阵为对称三对角矩阵. 于是，求原矩阵特征值问题，就转化为求上于是，求原矩阵特征值问题，就转化为求上海森伯格矩阵或对称三对角矩阵的特征值问题海森伯格矩阵或对称三对角矩阵的特征值问题.上页上页下页下页8.3.2 用正交类似变换用正交类似变换 约化普通实矩阵为上海森伯格矩阵约化普通实矩阵为上海森伯格矩阵 设设ARnn，下面来阐明，可选择初等反射，下面来阐明，可选择初等反射矩阵矩阵U1,U2,Un-2使使A经正交类似变换约化为一个经正交类似变换约化为一个上海森伯格阵上海森伯格阵.(1) 设设,)1(221)1(12112

28、12222111211 AcAaaaaaaaaaaAnnnnnn上页上页下页下页21/212112111 111121sgn()(),(3.1)().niiaaucea 其中其中c1=(a21,an1)TRn-1 ,无妨设无妨设c10，否那么这，否那么这一步不需求约化一步不需求约化. 于是于是, 可选择初等反射阵可选择初等反射阵 使使 ，其中，其中TuuIR11111 1111ecR 令令,111 RU上页上页下页下页那么那么 1)1(221111)1(12111112RARcRRAaUAUA(2)(2)(2)1112131(2)(2)(2)122232(2)(2)1112(2)(2)(2)3

29、2333(2)222(2)(2)(2)23,000nnnnnnnaaaaaaaAAaaacAaaa 其中其中.,),()2()2()2(222)2(2)2(322 nnnTnRARaac上页上页下页下页(2) 第第k步约化：反复上述过程，设对步约化：反复上述过程，设对A已完成第已完成第1步步,第第k-1步正交类似变换，即有步正交类似变换，即有,111 kkkkUAUA或或,11111 kkkUUAUUA且且(1)(2)(1)( )( )( )11121,111,11(2)(1)( )( )( )1222,122,12( )( )( )1,1( )( )( )1,1,11,( )( )( ),1

30、kkkkkkknkkkkkkknkkkkkkkk kknkkkkkkkknkkknkn knnaaaaaaaaaaaAaaaaaaaaa 上页上页下页下页,0)(22)(12)(11knkAcAAknkkkkk 其中其中 为为k阶上海森阶上海森伯格阵，伯格阵，)(11)()(, 1,),(kknTknkkkkkARaac .)()()(22knknkRA 设设ck0, 于是可选择初等反射阵于是可选择初等反射阵Rk使使 其中，其中，Rk计算公式为计算公式为,1ecRkkk 上页上页下页下页)2 . 3(.),(,)()(sgn(1)(, 112/112)()(, 1 Tkkkkkkkkkkkkk

31、nkikikkkkkuuIRaecuaa 令令, kkRIU那么那么)3 . 3(00)1(221)1(12)1(11)(22)(12)1(111 kkkkkkkkkkkkkkkkAcAARARcRRAAUAUA上页上页下页下页221122 nnUUAUUUU其中其中 为为k+1阶上海森伯格阵，第阶上海森伯格阵，第k步约化只需计步约化只需计算算 及及 (当当A为对称矩阵时，只需求计为对称矩阵时，只需求计算算 ).)1(11 kAkkRA)(12kkkRAR)(22kkkRAR)(22(3) 反复上述过程，那么有反复上述过程，那么有.1)(1)1(1, 12)3(332)2(22111 nnnn

32、nnnnnAaaaaa 上页上页下页下页 定理定理17 (豪斯霍尔德约化矩阵为上海森伯格阵豪斯霍尔德约化矩阵为上海森伯格阵) 设设ARnn那么存在初等反射矩阵那么存在初等反射矩阵U1,U2,Un-2 使使.00221122HAUUUUAUUUUTnn 为上海森伯格矩阵为上海森伯格矩阵.总结上述结论，有总结上述结论，有 算法算法1 (豪斯霍尔德约化矩阵为上海森伯格阵豪斯霍尔德约化矩阵为上海森伯格阵) 设设ARnn，本算法计算，本算法计算U0TAU0=H(上海森伯格型上海森伯格型)，其中其中U0=U1U2Un-2为初等反射阵的乘积为初等反射阵的乘积.1. U0I上页上页下页下页2. 对于对于k=1

33、,2,n-2(1) 计算初等反射阵计算初等反射阵Rk使使本算法约需求本算法约需求5n3/3次乘法运算，要明显构成次乘法运算，要明显构成U0还需求附加还需求附加2n3/3次乘法次乘法.,1ecRkkk (2) 约化计算约化计算, kkkkRIUAUUA(3) U0U0Uk上页上页下页下页例例7 用豪斯霍尔德方法将矩阵用豪斯霍尔德方法将矩阵 7242327341AA约化为上海森伯格阵约化为上海森伯格阵. 解解 选取初等反射阵选取初等反射阵R1使使 , 其中其中c1=(2,4)T.1111ecR (1) 计算计算R1:)() 1 , 5 . 0(, 4) 4 , 2max(11规规范范化化Tcc .

34、,472136. 4,809017. 1)5 . 0(,)1,618034. 1(,118034. 125. 11111111111TTuuIRecu 上页上页下页下页.1111ecR 那么有那么有(2) 约化计算约化计算:,111 RU那么得到上海森伯格阵为那么得到上海森伯格阵为.200000. 2399999. 00400000. 0799999. 7472136. 4447214. 0602631. 74112HAUUA 上页上页下页下页8.3.3 用正交类似变换用正交类似变换 约化对称矩阵为对称三对角矩阵约化对称矩阵为对称三对角矩阵 定理定理18 (豪斯霍尔德约化对称矩阵为对称三对豪斯

35、霍尔德约化对称矩阵为对称三对角矩阵角矩阵) 设设ARnn为对称矩阵，那么存在初等反为对称矩阵，那么存在初等反射矩阵射矩阵U1,U2,Un-2使使.111222111221122CcbbcbbcbbcUUAUUUUnnnnnn 为对称三对角矩阵为对称三对角矩阵.上页上页下页下页 证明证明 由定理由定理17, 存在初等反射矩阵存在初等反射矩阵U1,U2,Un-2 使使 为上海森伯为上海森伯格阵，且格阵，且An-1亦是对称阵，因此，亦是对称阵，因此，An-1为对称三对为对称三对角阵角阵.1221122 nnnAHUUAUUUU 由上面讨论可知，当由上面讨论可知，当A为对称阵时，由为对称阵时，由AkA

36、k+1 =Ak Uk Ak一步约化计算中只需计算一步约化计算中只需计算Rk及及Rk A22(k)Rk . 又由于又由于A的对称性，故只需计算的对称性，故只需计算Rk A22(k)Rk的对角线以下元素的对角线以下元素. 留意到留意到).)()(221)(221)(22TkkkkkTkkkkkkuuAAuuIRAR 上页上页下页下页引进记号引进记号)., 1;, 1()(22)(22ikjnkiuttuARARTkkTkkkkkk ,)(221knkkkkRuAr ,)(21knkkTkkkkRururt 那么有那么有对对称阵对对称阵A用初等反射阵正交类似约化为对角三用初等反射阵正交类似约化为对角

37、三对角阵大约需求对角阵大约需求2n3/3次乘法次乘法.用正交矩阵进展类似约化有一些特点，如构造的用正交矩阵进展类似约化有一些特点，如构造的 Uk容易求逆，且容易求逆，且Uk的元素数量级不大，这个算法是的元素数量级不大，这个算法是非常稳定的非常稳定的.算法算法2见书见书p318.上页上页下页下页8.4 QR 方方 法法8.4.1 QR算法算法Rutishauser(1958)利用矩阵的三角分解提出了计利用矩阵的三角分解提出了计算矩阵特征值的算矩阵特征值的LR算法，算法，Francis(1961,1962)利用矩利用矩阵的阵的QR分解建立了计算矩阵特征值的分解建立了计算矩阵特征值的QR算法算法.Q

38、R方法是一种变换方法，是计算普通矩阵方法是一种变换方法，是计算普通矩阵(中小中小型矩阵型矩阵)全部特征值问题的最有效方法之一全部特征值问题的最有效方法之一.上页上页下页下页 (1) 上海森伯格矩阵的全部特征值问题；上海森伯格矩阵的全部特征值问题； (2) 计算对称三对角矩阵的全部特征值问题，计算对称三对角矩阵的全部特征值问题， 目前目前QR方法主要用来计算：方法主要用来计算：且且QR方法具有收敛快，算法稳定等特点方法具有收敛快，算法稳定等特点. 对于普通矩阵对于普通矩阵ARnn (或对称矩阵或对称矩阵)，那么首，那么首先用豪斯霍尔德方法将先用豪斯霍尔德方法将A化为上海森伯格阵化为上海森伯格阵B

39、(或对称或对称三对角阵三对角阵)，然后再用，然后再用QR方法计算方法计算B的全部特征值的全部特征值. 设设ARnn，且，且A对进展对进展QR分解，即分解，即A=QR，上页上页下页下页其中其中R为上三角阵为上三角阵, Q为正交阵为正交阵, 于是可得到一新矩阵于是可得到一新矩阵B=RQ=QTAQ.显然，显然，B是由是由A经过正交类似变换得到，因此经过正交类似变换得到，因此B与与A的特征值一样的特征值一样. 再对再对B进展进展QR分解，又可得一新矩分解，又可得一新矩阵阵,反复这一过程可得到矩阵序列：反复这一过程可得到矩阵序列： 设设A=A1 将将A1进展进展QR分解分解A1=Q1R1 作矩阵作矩阵A

40、2=R1Q1=Q1TA1Q1 QR算法，就是利用矩阵的算法，就是利用矩阵的QR分解，按上述递分解，按上述递推法那么构造矩阵序列推法那么构造矩阵序列Ak的过程的过程. 只需只需A为非奇特为非奇特矩阵，那么由矩阵，那么由QR算法就完全确定算法就完全确定Ak.上页上页下页下页 定理定理19 (根本根本QR方法方法) 设设A=A1Rnn, 构造构造QR算算法法:)1 . 4(), 2 , 1(),(1 kQRARIQQRQAkkkkkTkkkk为为上上三三角角阵阵其其中中则则有有记记,1221RRRRQQQQkkkk 111,;TkkkkkkAAAQ A Q （ ）相相似似于于即即;)()()2(1211211kTkkTkkQAQQQQAQQQA