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

《第8章-矩阵特征值计算..ppt》由会员分享，可在线阅读，更多相关《第8章-矩阵特征值计算..ppt（107页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。

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

2、矩阵的特征值和特征向量的基础知识及性质知识及性质.上页上页下页下页8.1.1 特征值问题及性质特征值问题及性质求求A的特征值问题的特征值问题(1.1)等价于求等价于求A的特征方程的特征方程设矩阵设矩阵A Rnn，特征值问题是求，特征值问题是求 C和非零向和非零向量量x Rn，使，使的根的根.在在5.1.3节已经给出特征值的一些性质，下面节已经给出特征值的一些性质，下面再补充一些基本性质再补充一些基本性质.Ax=x (1.1)其中其中 是矩阵是矩阵A的的特征值特征值，x是矩阵是矩阵A属于特征值属于特征值 的的特征向量特征向量.本章讨论计算矩阵特征值的数值方法本章讨论计算矩阵特征值的数值方法,在在

3、科学和工程技术中很多问题在数学上都归结为求矩科学和工程技术中很多问题在数学上都归结为求矩阵的特征值问题阵的特征值问题.上页上页下页下页 定理定理1 设设 为为ARnn的特征值的特征值,且且Ax=x,x 0，则有则有 (2)-为为A-I的特征值，即的特征值，即(A-I)x=(-)x；(1)c 为的为的cA特征值特征值(c 0为常数为常数)；(3)k为为Ak的特征值，即的特征值，即Akx=kx；(4)设设A为为非奇异矩阵，那么非奇异矩阵，那么 0,且且-1为为A-1的特的特征值，即征值，即A-1x=-1x.上页上页下页下页 定理定理2 (1)设矩阵设矩阵ARnn可对角化，即存在非可对角化，即存在非

4、奇异矩阵奇异矩阵P使使的充分必要条件是的充分必要条件是A具有具有n个线性无关的特征向量个线性无关的特征向量.(2)如果矩阵如果矩阵ARnn有有m个个(mn)不同的特征不同的特征值值 1,2,m，则对应的特征向量则对应的特征向量 x1,x2,xm 线性无关线性无关.上页上页下页下页 定理定理3(对称矩阵的正交约化对称矩阵的正交约化)设设ARnn为对称为对称矩阵，则矩阵，则 (3)存在一个正交矩阵存在一个正交矩阵P使的使的且且 1,2,n为为A的特征值，而的特征值，而P(u1,u2,un)的的列向量列向量uj为为A的对应于的对应于 j 的单位特征向量的单位特征向量.(1)A的特征值均为实数；的特征

5、值均为实数；(2)A有有n个线性无关的特征向量；个线性无关的特征向量；上页上页下页下页 定义定义 设设ARnn为为对称矩阵对称矩阵，对于任一非零向对于任一非零向量量x 0，称，称为对应于向量为对应于向量x的的瑞利瑞利(Rayleigh)商商.定理定理4 设设ARnn为为对称矩阵对称矩阵(其特征值依次记其特征值依次记为为 1 2 n)，则则 (1).(对任何非零对任何非零xRn)；(2).；(1.3)上页上页下页下页 证明证明 只证只证(1)，关于，关于(2)自己作练习自己作练习.由于由于A为实对称矩阵，可将为实对称矩阵，可将 1,2,n 对应的对应的特征向量特征向量 x1,x2,xn 正交规范

6、化，则有正交规范化，则有(xi,xj)=ij，设，设x 0为为Rn中中任一向量，则有任一向量，则有于是于是从而从而(1)成立成立.结论结论(1)说明说明瑞利商瑞利商必位于必位于 n和和 1之之间间.上页上页下页下页8.1.2 特征值估计与扰动特征值估计与扰动 定义定义1 设设n阶矩阵阶矩阵A=(aij)，令，令 (1)；(2)集合集合 称为复平称为复平面上以面上以aii为为圆心，以圆心，以ri为半径的为半径的n阶矩阵阶矩阵A的的n个个格什格什戈林戈林(Gerschgorin)圆盘圆盘.上页上页下页下页 定理定理5(Gerschgorin圆盘定理圆盘定理)特别地，如果特别地，如果A的一个圆盘的一

7、个圆盘Di是与其它圆盘分离是与其它圆盘分离(即即孤立圆盘孤立圆盘)，则，则Di中中精确地包含精确地包含A的一个特征值的一个特征值.(1)设设n阶矩阵阶矩阵A(aij)，则，则A的每一个特征值必的每一个特征值必属于下面某个圆盘之中属于下面某个圆盘之中 (2)如果如果A有有m个个圆盘圆盘组成一个连通的并集组成一个连通的并集S，且且S与余下与余下n-m个个圆盘圆盘是分离的，则是分离的，则S内恰包含内恰包含A的的m个特征值个特征值.或者说或者说 A的特征值都在的特征值都在n个圆盘的个圆盘的并集并集中中.上页上页下页下页 证明证明 只就只就(1)给出证明给出证明.设设 为为A的特征值，即的特征值，即 A

8、x=x，其中，其中 x=(x1,x2,xn)T 0.或或 记记 ，考虑，考虑Ax=x的第的第k个方程，个方程，即即于是于是即即上页上页下页下页这说明，这说明，A的每一个特征值必位于的每一个特征值必位于A的一个圆盘的一个圆盘中，并且相应的特征值中，并且相应的特征值 一定位于第一定位于第k个圆盘中个圆盘中(其中其中k是对应特征向量是对应特征向量x绝对值最大的分量的下标绝对值最大的分量的下标).利用相似矩阵性质，有时可以获得利用相似矩阵性质，有时可以获得A的特征值进的特征值进一步的估计，即适当选取非奇异对角阵一步的估计，即适当选取非奇异对角阵并做相似变换并做相似变换 .适当选取适当选取 可使可使某些

9、圆盘半径及连通性发生变化某些圆盘半径及连通性发生变化.上页上页下页下页 例例1 估计矩阵估计矩阵A的特征值范围，其中的特征值范围，其中 解解 矩阵矩阵A的的3个圆盘为个圆盘为 由定理由定理8，可知，可知A的的3个特征值位于个特征值位于3个圆盘的并个圆盘的并集中，由于集中，由于D1是是孤立圆盘，所以孤立圆盘，所以D1内恰好包含内恰好包含A的一的一个特征值个特征值 1(为实特征值为实特征值)，即，即A的其它两个特征值的其它两个特征值 2,3包含在包含在D2,D3的并集中的并集中.上页上页下页下页现在取对角阵现在取对角阵做相似变换做相似变换 矩阵矩阵A1的的3个圆盘为个圆盘为上页上页下页下页显然，显

10、然，3个圆盘都是孤立圆盘，所以，每一个圆个圆盘都是孤立圆盘，所以，每一个圆盘都包含盘都包含A的一个特征值的一个特征值(为实特征值为实特征值)，且有估计，且有估计上页上页下页下页 定理定理6(Bauer-Fike定理定理)设设 是是A+IRnn的一个的一个特征值，且特征值，且P-1AP=D=diag(1,2,n)，则有，则有其中其中|p为矩阵的为矩阵的p范数，范数，p=1,2,.证明证明 由于由于 (A)时时显然成立显然成立,故只考虑故只考虑 (A).这时这时D-I非奇异，设非奇异，设x是是A+I对应于对应于 的特征向量，由的特征向量，由(A+I-I)x=0左乘左乘P-1可得可得而对角矩阵而对角

11、矩阵(D-I)-1的范数为的范数为所以有所以有这就得到这就得到(1.5)式式.这时总有这时总有(A A)中的一个中的一个 取到取到m值值.上页上页下页下页 由定理由定理6可知可知|P-1|P|=cond(P)是特征值扰动的是特征值扰动的放大系数，但将放大系数，但将A对角化的相似变化矩阵不是唯一的，对角化的相似变化矩阵不是唯一的，所以取所以取cond(P)的下确界的下确界称为特征值问题的条件数称为特征值问题的条件数.只要只要v(A)不很大，矩阵微不很大，矩阵微小扰动只带来特征值的微小扰动小扰动只带来特征值的微小扰动.但是但是v(A)难以计算难以计算,有时只对一个有时只对一个P，用，用cond(P

12、)代替代替v(A).特征值问题的条件数和解线性方程组时的矩阵是特征值问题的条件数和解线性方程组时的矩阵是两个不同的概念，对于一个矩阵两个不同的概念，对于一个矩阵A，两者可能一大一，两者可能一大一小，例如二阶矩阵小，例如二阶矩阵A=diag(1,10-10)，有，有v(A)=1，但解，但解线性方程组的矩阵条件数线性方程组的矩阵条件数cond(A)=1010.上页上页下页下页 关于计算矩阵关于计算矩阵A的特征值问题，当的特征值问题，当n=2,3时，我们时，我们还可按行列式展开的办法求特征方程还可按行列式展开的办法求特征方程p()=0的根的根.但但当当n较大时，如果按展开行列式的办法较大时，如果按展

13、开行列式的办法,首先求出首先求出p()的系数，再求的系数，再求p()的根，工作量就非常大，用这的根，工作量就非常大，用这种办法求矩阵特征值是不切实际的，由此需要研究种办法求矩阵特征值是不切实际的，由此需要研究A的特征值及特征向量的数值方法的特征值及特征向量的数值方法.本章将介绍一些计算机上常用的两类方法，一类本章将介绍一些计算机上常用的两类方法，一类是幂法及反幂法是幂法及反幂法(迭代法迭代法)，另一类是正交相似变化的，另一类是正交相似变化的方法方法(变化法变化法).上页上页下页下页8.2 幂法及反幂法幂法及反幂法 幂法与反幂法都是求幂法与反幂法都是求实矩阵实矩阵的特征值和特征向的特征值和特征向

14、量的量的向量迭代法向量迭代法，所不同的是幂法是计算矩阵的，所不同的是幂法是计算矩阵的主主特征值特征值(矩阵矩阵按模最大的特征值按模最大的特征值称为称为主特征值主特征值，其模，其模就是该矩阵的就是该矩阵的谱半径谱半径)和和相应特征向量相应特征向量的一种向量迭的一种向量迭代法，特别适用于大型稀疏矩阵代法，特别适用于大型稀疏矩阵.而反幂法则是计而反幂法则是计算非奇异算非奇异(可逆可逆)矩阵矩阵按模最小的特征值按模最小的特征值和和相应特征相应特征向量向量的一种向量迭代法，特别是计算海森伯格阵或的一种向量迭代法，特别是计算海森伯格阵或三对角阵的对应一个给定近似特征值的特征向量的三对角阵的对应一个给定近似

15、特征值的特征向量的有效方法之一有效方法之一.下面分别介绍幂法与反幂法下面分别介绍幂法与反幂法.上页上页下页下页8.2.1 幂法幂法(又称又称乘幂法乘幂法)现讨论求现讨论求 1及及x1的方法的方法.设实矩阵设实矩阵A=(aij)有有一个一个完全的特征向量组完全的特征向量组，即，即A有有n个线性无关的特征向量个线性无关的特征向量，设矩阵，设矩阵A的特征值为的特征值为 1,2,n,相应的特征向量为相应的特征向量为x1,x2,xn.已知已知A的主特征值的主特征值 1是实根，且满足条件是实根，且满足条件上页上页下页下页 幂法的幂法的基本思想基本思想是是:任取任取非零非零的初始向量的初始向量v0,由由矩阵

16、矩阵A构造一向量序列构造一向量序列vk称为迭代向量，由假设，称为迭代向量，由假设，v0可唯一表示为可唯一表示为上页上页下页下页于是于是其中其中由假设由假设 故故 从而从而为为 1的特征向量的特征向量.上页上页下页下页所以当所以当k充分大时，有充分大时，有即为矩阵即为矩阵A的的对应特征值对应特征值 1 的一个近似特征向量的一个近似特征向量.用用(vk)i 表示表示vk的第的第i个分量，则当个分量，则当k充分大时，有充分大时，有即为即为A的的主特征值主特征值 1的近似值的近似值.由于由于 这种由已知非零向量这种由已知非零向量v0及矩阵及矩阵A的乘幂的乘幂Ak构造向构造向量序列量序列vk以计算以计算

17、A的的主特征值主特征值 1(利用利用(2.7)式式)及相及相应特征向量应特征向量(利用利用(2.5)式式)的方法就称为的方法就称为(乘乘)幂法幂法.上页上页下页下页 迭代公式实质上是由矩阵迭代公式实质上是由矩阵A的乘幂的乘幂 Ak与非零向量与非零向量v0相乘来构造向量序列相乘来构造向量序列vk=Akv0，从而计算主特征从而计算主特征值值 1及其对应的特征向量，这就是及其对应的特征向量，这就是幂法幂法的思想的思想.的收敛速度由比值的收敛速度由比值来确定，来确定，r 越小收敛越快，但当越小收敛越快，但当r 1 1时收敛可能很慢时收敛可能很慢.总结上述讨论，有下面的定理总结上述讨论，有下面的定理.上

18、页上页下页下页 定理定理7 设设ARnn有有n个线性无关的特征向量，个线性无关的特征向量，主特征值主特征值 1满足满足|1|2|n|，则对则对任何非零向量任何非零向量v0(a1 0)，(2.4)式和式和(2.7)式成立式成立.又设又设A有有n个线性无关的特征向量，个线性无关的特征向量，1对应的对应的r个线性个线性无关的特征向量为无关的特征向量为x1,x2,xr，则由则由(2.2)式有式有 如果如果A的主特征值为实的重根的主特征值为实的重根,即即 1=2=r,且且|r|r+1|n|，上页上页下页下页为为 1的特征向量，这说明当的特征向量，这说明当A的主特征值是实的重根的主特征值是实的重根时，定理

19、时，定理7的结论还是正确的的结论还是正确的.应用应用幂法幂法计算计算A的主特征值的主特征值 1及其对应的特征向及其对应的特征向量时量时，如果如果|1|1(或或|1|2|n|，则对，则对任意非零任意非零初始向量初始向量v0=u0(a1 0)，有，有幂法计算公式为幂法计算公式为(uk,vk)则有则有 (1)(2)上页上页下页下页例例2 用幂法计算矩阵用幂法计算矩阵的主特征值和相应的特征向量的主特征值和相应的特征向量.解解取取 v0=u0=(1,1,1)T,则则上页上页下页下页计算结果如下表计算结果如下表k0(1,1,1)1(0.9091,0.8182,1)2.75000005(0.7651,0.6

20、674,1)2.558791810(0.7494,0.6508,1)2.538002915(0.7483,0.6497,1)2.536625616(0.7483,0.6497,1)2.536584017(0.7482,0.6497,1)2.536559818(0.7482,0.6497,1)2.536545619(0.7482,0.6497,1)2.536537420(0.7482,0.6497,1)2.5365323上页上页下页下页这个结果是用这个结果是用8位浮点数字进行运算得到的位浮点数字进行运算得到的,uk的分量值是舍入值的分量值是舍入值.于是得到于是得到及相应的特征向量及相应的特征向量

21、(0.7482,0.6497,1)T.1和相应的和相应的特征向量的真值特征向量的真值(8位数字位数字)为为上页上页下页下页例例 用幂法计算矩阵用幂法计算矩阵的主特征值与其对应的特征向量的主特征值与其对应的特征向量.解解取取 v0=u0=(0,0,1)T ,则则上页上页下页下页直到直到k=8 时的计算结果见下表时的计算结果见下表k1 2,4,1,4 0.5,1,0.252 4.5,9,7.75 90.5,1,0.86113 5.7222,11.4444,8.36111.44440.5,1,0.73604 5.4621,10.9223,8.2306 10.92230.5,1,0.75365 5.5

22、075,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.7500从而从而上页上页下页下页8.2.2 加速方法加速方法1、原点平移法原点平移法 由前面讨论知道，应用幂法计算由前面讨论知道，应用幂法计算A的主特征值的的主特征值的收敛速度主要由比值收敛速度主要由比值 r=|2/1|来决定，但当来决定，但当r 接近于接近于1时，收敛可能

23、很慢时，收敛可能很慢.这时，一个补救办法是采用加这时，一个补救办法是采用加速收敛的方法速收敛的方法.其中其中p为参数，设为参数，设A的特征值为的特征值为 i，则，则对矩阵对矩阵B的特征的特征值为值为 i-p，而且而且A,B的特征向量相同的特征向量相同.引进矩阵引进矩阵 B=A-pI.上页上页下页下页 如果要计算如果要计算A的主特征值的主特征值 1,只要只要选择合适的数选择合适的数p，使使 1-p为矩阵为矩阵B=A-pI 的主特征值，且的主特征值，且 那么，对矩阵那么，对矩阵B=A-pI应用幂法求其主特征值应用幂法求其主特征值 1-p,收收敛速度将会加快敛速度将会加快.这种通过求这种通过求B=A

24、-pI的主特征值和特的主特征值和特征向量，而得到征向量，而得到A的主特征值和特征向量的方法叫的主特征值和特征向量的方法叫原原点平移法点平移法.对于对于A的特征值的某种分布，它是十分有的特征值的某种分布，它是十分有效的效的.上页上页下页下页例例3 设设AR44有特征值有特征值比值比值r=|2/1|0.9.做变换做变换B=A-12I (p=12),则则B的的特征值为特征值为应用幂法计算应用幂法计算B的主特征值的主特征值1的收敛速度的比值为的收敛速度的比值为 虽然常常能够选择有利的虽然常常能够选择有利的p值值,使幂法得到加速使幂法得到加速,但设计一个自动选择适当参数但设计一个自动选择适当参数p的过程

25、是困难的的过程是困难的.上页上页下页下页 下面考虑当下面考虑当A的特征值是实数时，怎样选择的特征值是实数时，怎样选择p使采使采用幂法计算用幂法计算 1得到加速得到加速.且使且使收敛速度的比值收敛速度的比值 设设A的特征值都是实数，且满足的特征值都是实数，且满足则对实数则对实数p，使矩阵使矩阵A-pI的主特征值为的主特征值为 1-p或或 n-p时，时，当当我们计算我们计算 1及及x1时，时，首先应选取首先应选取p使使上页上页下页下页显然，当显然，当 2-p=-(n-p)时，即时，即P=(2+n)/2P*时时 为最小值，这时为最小值，这时收敛速度的比值收敛速度的比值为为 当希望计算当希望计算 n时

26、，应选取时，应选取 p=(1+n-1)/2P*使得应用幂法计算使得应用幂法计算 n得到加速得到加速.当当A的特征值都是实数，满足的特征值都是实数，满足(2.10)式式且且 2,n能初步估计出来，我们就能确定能初步估计出来，我们就能确定P*的近似值的近似值.上页上页下页下页 例例4 用原点平移加速法求用原点平移加速法求例例2中矩阵中矩阵A的主特征值的主特征值与其对应的特征向量与其对应的特征向量.对对B应用幂法，仍应用幂法，仍取取 v0=(1,1,1)T ,则则 解解 取取p=0.75,做平移变换做平移变换B=A-pI，则则上页上页下页下页k0(1,1,1)1(0.8750,0.7500,1)2.

27、00000005(0.7516,0.6522,1)1.79140116(0.7491,0.6511,1)1.78884437(0.7488,0.6501,1)1.78733008(0.7484,0.6499,1)1.78691529(0.7483,0.6497,1)1.786658710(0.7482,0.6497,1)1.7865914计算结果如下表计算结果如下表由此得由此得B的主特征值为的主特征值为 1 1.7865914，A的主特征值的主特征值 1为为 1=1+0.75 2.5365914.上页上页下页下页 原点位移的加速方法，是一个矩阵变换方法原点位移的加速方法，是一个矩阵变换方法.这

28、这种变换容易计算，又不破坏矩阵种变换容易计算，又不破坏矩阵A的稀疏性，但的稀疏性，但p的的选择选择依赖依赖对对A的的特征值分布的大致了特征值分布的大致了解解.及相应的特征向量为及相应的特征向量为 x1(0.7482,0.6497,1)T.与例与例2结果比较，上述结果比例结果比较，上述结果比例2迭代迭代15次还好次还好.若若迭代迭代15次次 1 1.7865258(相应的相应的 1 2.5365258)就是真就是真值值.上页上页下页下页 设设ARnn为为对称矩阵对称矩阵，称，称为向量为向量x的的瑞利商瑞利商，其中，其中(x,x)=xTx为内积为内积.由定理由定理4知道，实对称矩阵知道，实对称矩阵

29、A的特征值的特征值 1及及 n可用可用瑞利商瑞利商的的极限值表示极限值表示.下面我们将下面我们将瑞利商瑞利商应用到用幂法计算应用到用幂法计算实对称矩阵实对称矩阵A的主特征值的加速上来的主特征值的加速上来.2、瑞利商、瑞利商(Rayleigh)加速加速上页上页下页下页 定理定理9 设设ARnn为为对称矩阵对称矩阵，特征值满足，特征值满足对应的特征向量对应的特征向量xi满足满足(xi,xj)=ij(单位正交向量单位正交向量),应用幂法公式应用幂法公式(2.9)计算计算A的主特征值的主特征值 1，则规范化，则规范化向量向量uk的的瑞利商瑞利商给出给出 1的较好的近似值为的较好的近似值为 由此可见，由

30、此可见，R(uk)比比 k更快的收敛于更快的收敛于 1.上页上页下页下页 证明证明 由由(2.8)式及式及得得上页上页下页下页 幂法的幂法的瑞利商加速迭代公式瑞利商加速迭代公式可以写为可以写为其中其中A为为n阶实对称矩阵阶实对称矩阵.对给定的误差限对给定的误差限，当，当|k k-1|时，取近似时，取近似值值上页上页下页下页8.2.3 反幂法反幂法 反幂法是用于求非奇异矩阵反幂法是用于求非奇异矩阵A的的按模最小的特征值按模最小的特征值和对应特征向量和对应特征向量的方法的方法.而结合原点平移法的反幂法而结合原点平移法的反幂法则可以求矩阵则可以求矩阵A的任何一个的任何一个给定近似特征值对应的特征给定

31、近似特征值对应的特征向量向量.设矩阵设矩阵A非奇异非奇异,其特征值其特征值 i (i=1,2,n),满足满足其相应的特征向量其相应的特征向量x1,x2,xn线性无关，则线性无关，则 A-1 的特的特征值为征值为1/i，对应的特征向量仍为对应的特征向量仍为xi(i=1,2,n).上页上页下页下页此时此时,A-1的特征值满足的特征值满足因此因此,对对A-1应用幂法应用幂法,可求出其可求出其主特征值主特征值 1/n k 和和特征向量特征向量 xn uk.从而求得从而求得A的的按模最小特征值按模最小特征值 n 1/k 和对应的和对应的特征向量特征向量 xn uk，这种求这种求A-1的方法就称为的方法就

32、称为反幂法反幂法.上页上页下页下页为了避免求为了避免求A-1,可通过解线性方程组可通过解线性方程组Avk=uk-1得到得到vk，采用采用LU分解法，即先对分解法，即先对A进行进行LU分解分解A=LU,此时此时反幂反幂法的迭代公式法的迭代公式为为 反幂法的迭代公式反幂法的迭代公式为为:任取任取v0=u0 0 0，构造，构造上页上页下页下页对给定的误差对给定的误差，当，当|kk-1|n|0，则对则对任意非零初始向量任意非零初始向量u0(an 0)，由反幂法计算公由反幂法计算公式构造的向量序列式构造的向量序列vk，uk满足满足 (1)(2)收敛速度收敛速度的比值的比值r=|n/n-1|.上页上页下页

33、下页 在反幂法中也可以用在反幂法中也可以用原点平移法原点平移法加速迭代过程加速迭代过程，或或求其它特征值与其对应的特征向量求其它特征值与其对应的特征向量.如果矩阵如果矩阵(A-pI)-1存在，显然其特征值为存在，显然其特征值为对应的特征向量仍然是对应的特征向量仍然是x1,x2,xn，现对矩阵现对矩阵(A-pI)-1应用幂法，得到反幂法的迭代公式应用幂法，得到反幂法的迭代公式上页上页下页下页 如果如果p是是A的特征值的特征值 j的一个近似值，且设的一个近似值，且设 j与其它与其它特征值是分离的，即特征值是分离的，即就是说就是说1/(j-p)是矩阵是矩阵(A-pI)-1的主特征值，可用反幂的主特征

34、值，可用反幂法法(2.12)计算特征值计算特征值 j及特征向量及特征向量xj.设设ARnn有有n个线性无关的特征向量个线性无关的特征向量x1,x2,xn,则则上页上页下页下页其中其中同理可得下面的定理同理可得下面的定理.上页上页下页下页 定理定理11 设设ARnn有有n个线性无关的特征向量，个线性无关的特征向量，矩阵矩阵A的特征值及对应的特征向量分别记为的特征值及对应的特征向量分别记为 i 及及xi(i=1,2,n)，而，而p为为 j的近似值，的近似值，(A-pI)-1存在，且存在，且 (1)(2)则对则对任意非零初始向量任意非零初始向量u0(aj 0)，由反幂法计算公由反幂法计算公式式(2.

35、12)构造的向量序列构造的向量序列vk，uk满足满足且收敛速度为且收敛速度为上页上页下页下页 由该定理知，对由该定理知，对A-pI(其中其中p j)应用反幂法，可应用反幂法，可用来计算特征向量用来计算特征向量xj，只要选择只要选择p是是 j的一个较好的近的一个较好的近似且特征值分离情况较好，一般似且特征值分离情况较好，一般r很小，常常只要迭代很小，常常只要迭代一二次就可完成特征向量的计算一二次就可完成特征向量的计算.反幂法迭代公式中的反幂法迭代公式中的vk是是通过解方程组通过解方程组求得的求得的,为了节省工作量为了节省工作量,可以先将可以先将A-pI进行三角分解进行三角分解于是求于是求vk相对

36、于解两个三角形方程组相对于解两个三角形方程组上页上页下页下页实验表明实验表明,按下述方法选择按下述方法选择u0是是较好的较好的:选选u0使使用回代求解三角形方程组用回代求解三角形方程组(2.13)即得即得v1，然后再按公，然后再按公式式(2.12)进行迭代进行迭代.上页上页下页下页反幂法计算公式：反幂法计算公式：1.分解计算分解计算P(A-pI)=LU,且保留且保留L,U及及P信息信息2.反幂迭代法反幂迭代法(1)解解Uv1=(1,1,1)T求求v1(2)解解k=2,3,解解Lyk=Puk-1求求yk 解解Uvk=yk求求vk k=max(vk)计算计算uk=vk/k上页上页下页下页例例5 用

37、反幂法求用反幂法求的对应于计算特征值的对应于计算特征值=1.2679(精确特征值为精确特征值为 )的的特征向量特征向量(用用5位浮点数进行计算位浮点数进行计算).解解 用部分选主元的三角分解将用部分选主元的三角分解将A-pI(其中其中p=1.2679)分解为分解为 P(A-pI)=LU,其中其中上页上页下页下页由由Uv1=(1,1,1)T,得得由由LUv2=Pu1,得得 3对应的特征向量是对应的特征向量是由此看出由此看出u2是是x3的相当好的近似的相当好的近似.特征值特征值 3的真值为的真值为上页上页下页下页8.3 正交变化与矩阵分解正交变化与矩阵分解正交变化是计算矩阵特征值的有力工具，本节介

38、正交变化是计算矩阵特征值的有力工具，本节介绍豪斯绍豪斯霍尔德霍尔德(Householder)变换变换和和吉文斯吉文斯(Givens)变换变换，并利用它们讨论，并利用它们讨论矩阵分解矩阵分解.主要讨论主要讨论实矩阵实矩阵和和实向量实向量.上页上页下页下页8.3.1 豪斯霍尔德变换豪斯霍尔德变换定义定义2 设向量设向量w Rn，且，且wTw=1，称矩阵，称矩阵为为初等反射矩阵初等反射矩阵，这个矩阵也称为，这个矩阵也称为豪斯霍尔德变换豪斯霍尔德变换.如果记如果记w=(w1,w2,wn)，则有，则有上页上页下页下页定理定理12 设有初等反射矩阵，其中设有初等反射矩阵，其中H=I-2wwT，其，其中中w

39、Tw=1，则，则(1)H是对称矩阵，即是对称矩阵，即HT=H.(2)H是正交矩阵，即是正交矩阵，即H-1=H.(3)设设A为对称矩阵，那么为对称矩阵，那么 A1=H-1AH=HAH 即亦即亦是对称矩阵是对称矩阵.证明证明 只证只证H的正交性，其中显然的正交性，其中显然.HTH=H2=(I-2wwT)(I-2wwT)=I-4wwT+4w(wTw)wT=I.设向量设向量u 0，则显然，则显然是一个初等反射矩阵是一个初等反射矩阵.上页上页下页下页vwSxv 图图 8-18-1Oy下面考察下面考察初等反射矩阵初等反射矩阵的的几何意义几何意义.参见参见图图8-1，考，考虑以为虑以为w法向量且过原点法向量

40、且过原点O的的超平面超平面S:wTx=0.设任意向量设任意向量v Rn,则则v=x+y,其中其中x S,y ST.于是于是 Hx=(I-2wwT)x=x-2wwTx=x.对于对于y ST，易知，易知Hy=-y,从而对从而对任意向量任意向量v Rn，总有，总有 Hv=x-y=v.其中其中v 为为v关于平面关于平面S的镜面反射的镜面反射(见图见图8-1).初等反射矩阵在计算上的意义是它能用来约化矩初等反射矩阵在计算上的意义是它能用来约化矩阵，例如设向量阵，例如设向量Hx 0，可选择一初等反射矩阵，可选择一初等反射矩阵H，使使Hx=y.上页上页下页下页定理定理13 设设x,y为两个不相等的为两个不相

41、等的n维向量维向量,|x|2=|y|2,则存在一个初等反射矩阵则存在一个初等反射矩阵H，使，使Hx=y.证明证明 令令 ,则得到一个初等反射矩阵则得到一个初等反射矩阵而且而且因为因为所以所以容易说明，容易说明，w是使是使Hx=y成立的唯一长度等于成立的唯一长度等于1的的向量向量(不计符号不计符号).上页上页下页下页定理定理14(约化定理约化定理)设设x=(x1,x2,xn)T 0,则存则存在初等反射矩阵在初等反射矩阵H，使，使Hx=-e1，其中，其中证明证明 记记y=-e1，设，设x y，取，取=|x|2，则有，则有|x|2=|y|2，于是由定理，于是由定理13存在存在H变换变换 H=I-2w

42、wT，其中其中 ，使，使Hx=y=-e1.上页上页下页下页记记u=x+e1=(u1,u2,un)T，于是，于是其中记其中记u=(x1+,x2,xn)T，.显然显然如果如果和和x1异号，那么计算异号，那么计算x1+时有效数字可能时有效数字可能损失，我们取损失，我们取和和x1有相同的符号，即取有相同的符号，即取在计算在计算时，可能上溢或下溢，为了避免溢出，时，可能上溢或下溢，为了避免溢出，将将x规范化规范化上页上页下页下页则有则有H 使使 H x =e1，其中，其中例例6 设设x=(3,5,1,1)T，则，则|x|2=6.取取=6,可直接验证可直接验证Hx=-e1=(-6,0,0,0)T.上页上页

43、下页下页8.3.2 吉文斯变换吉文斯变换设向量设向量x,y R2，则变换，则变换或或是平面上向量的一个旋转变换，其中是平面上向量的一个旋转变换，其中为正交变换为正交变换.Rn中变换中变换其中其中x=(x1,x2,xn)T,y=(y1,y2,yn)T，而，而上页上页下页下页称为称为Rn中平面中平面xi,xj的的旋转变换旋转变换，也称为，也称为吉文斯变换吉文斯变换.P=P(i,j,)=P(i,j)称为平面旋转矩阵称为平面旋转矩阵.上页上页下页下页显然，显然，P=P(i,j,)具有性质：具有性质：(1)P与单位阵与单位阵I只是在只是在(i,i),(i,j),(j,i),(j,j)位置元素不一样，其它

44、相同位置元素不一样，其它相同.(2)P为正交矩阵为正交矩阵(P-1=PT).(3)P(i,j)A(左乘左乘)只需计算第只需计算第i行与第行与第j行元素，即行元素，即对对A=(aij)mn有有其中其中c=cos,s=sin.(4)AP(i,j)(右乘右乘)只需计算第只需计算第i列与第列与第j列元素列元素利用平面旋转变换，可得向量利用平面旋转变换，可得向量x中的指定元素变为零中的指定元素变为零.上页上页下页下页定理定理15(约化定理约化定理)设设x=(x1,xi,xj,xn)T,其中不全为零，则可选择平面旋转阵其中不全为零，则可选择平面旋转阵P=P(i,j,)，使，使其中其中于是，由于是，由c,s

45、的取法得的取法得证明证明 取取c=cos=xi/x i,s=sin =xj/x i，由，由P(i,j,)x=x=(x 1,x i,x j,x n)T,利利用矩阵乘法用矩阵乘法,显然有显然有上页上页下页下页8.3.3 矩阵的矩阵的QR分解与舒尔分解分解与舒尔分解 定理定理16 设设ARnn非奇异，则存在正交矩阵非奇异，则存在正交矩阵P，使使PA=R，其中，其中R为上三角矩阵为上三角矩阵.证明证明 我们先用吉文斯变换给出构造我们先用吉文斯变换给出构造P的方法的方法 (1)第第1步约化，由设有步约化，由设有j(j=1,2,n)使使aj1 0，则，则可选择吉文斯变换可选择吉文斯变换P(1,j)，将，将

46、aj1处的元素化为零处的元素化为零.即即若若aj1 0(j=2,3,n)，则存在，则存在P(1,j)使得使得可简记为可简记为P1A=A(2)，其中，其中P1=P(1,n)P(1,2).上页上页下页下页 (2)第第k步约化，设上述过程已完成第步约化，设上述过程已完成第1步至第步至第k-1步，于是有步，于是有 由设有由设有j(n j k)使使ajk(k)0(j=k+1,n)，则可，则可选择吉文斯变换选择吉文斯变换P(k,j)(j=k+1,n)使使 PkA(k)=P(k,n)P(k,k+1)A(k)=PkPk-1P1A=A(k+1),其中其中 Pk=P(k,n)P(k,k+1).上页上页下页下页 (

47、3)继续上述约化过程，最后则有继续上述约化过程，最后则有Pn-1P2P1A=A(n)=R (上三角形矩阵上三角形矩阵).令令P=Pn-1P2P1,它是一个正交矩阵它是一个正交矩阵,有有PA=R.也可以用豪斯霍尔德变换构造出正交矩阵也可以用豪斯霍尔德变换构造出正交矩阵P，记，记A(0)=A，它的第一列记为，它的第一列记为u1，不妨设，不妨设a1(0)0，可按公，可按公式式(3.2)找到矩阵找到矩阵 使使于是于是其中其中上页上页下页下页 一般地，设一般地，设其中其中D(j-1)为为j-1阶方阵，其对角线以下元素均为阶方阵，其对角线以下元素均为0，矩，矩阵阵 A(j-1)为为n-j+1阶方阵，设其第

48、一列为阶方阵，设其第一列为a1(j-1)，可选，可选择择n-j+1的豪斯霍尔德矩阵变换的豪斯霍尔德矩阵变换使使根据根据 Hj构造构造nn阶的变换矩阵阶的变换矩阵Hj为为上页上页下页下页于是有于是有它和它和A(j-1)有类似的形式，只是有类似的形式，只是D(j)为为j阶方阵，其对阶方阵，其对角线以下元素是角线以下元素是0，这样经过，这样经过n-1步运算得到步运算得到Hn-1H2H1A=A(n-1)=R (上三角形矩阵上三角形矩阵).其中其中P=Hn-1H2H1是一个正交矩阵是一个正交矩阵,从而有从而有PA=R.上页上页下页下页定理定理17(QR分解定理分解定理)设设ARnn为非奇异矩阵为非奇异矩

49、阵,则存在正交矩阵则存在正交矩阵Q与上三角矩阵与上三角矩阵R，使，使A=QR，且当且当R的对角元素为正时，分解是唯一的的对角元素为正时，分解是唯一的.证明证明 从定理从定理16可知，只要令可知，只要令Q=PT就有就有A=QR.下下面证明分解的唯一性，设有两种分解面证明分解的唯一性，设有两种分解A=Q1R1=Q2R2，其中其中Q1,Q2为正交矩阵，为正交矩阵，R1,R2为对角元素均为正的为对角元素均为正的上三角矩阵，则上三角矩阵，则由假设及对称正定矩阵由假设及对称正定矩阵ATA的楚列斯基分解的唯一性，的楚列斯基分解的唯一性，则得则得R1=R2，从而可得，从而可得Q1=Q2.证毕证毕.上页上页下页

50、下页定理定理16保证了保证了A可分解为可分解为A=QR.若若A非奇异，则非奇异，则R也非奇异也非奇异.如果不规定如果不规定R的对角元素为正，则分解不的对角元素为正，则分解不是唯一的是唯一的.一般按吉文斯变换或豪斯霍尔德变换方法一般按吉文斯变换或豪斯霍尔德变换方法作出的分解作出的分解A=QR，R的对角元不一定是正的，设上三的对角元不一定是正的，设上三角矩阵角矩阵R=(rij)，只要令，只要令则则 为正交矩阵，为正交矩阵，为对角元是为对角元是|rii|的上三的上三角矩阵，这样角矩阵，这样 便是符合定理便是符合定理17的唯一的唯一QR分解分解.上页上页下页下页例例7 用豪斯霍尔德变换作矩阵用豪斯霍尔