计算方法第四章矩阵特征值和特征向量的计算课件

《计算方法第四章矩阵特征值和特征向量的计算课件》由会员分享，可在线阅读，更多相关《计算方法第四章矩阵特征值和特征向量的计算课件（32页珍藏版）》请在装配图网上搜索。

1、计算方法第四章矩阵特征值和特征向量的计算1 第四章第四章 矩阵特征值和特征向量的计算矩阵特征值和特征向量的计算 工程实践中有多种振动问题，如桥梁工程实践中有多种振动问题，如桥梁 或建筑物或建筑物的振动，机械机件、飞机机翼的振动，及的振动，机械机件、飞机机翼的振动，及 一些稳定一些稳定性分析和相关分析可转性分析和相关分析可转 化为求矩阵特征值与特征向化为求矩阵特征值与特征向量的问题。量的问题。 但高次多项式求根精度低但高次多项式求根精度低 , 一般不作为求解方一般不作为求解方法法. 目前的方法是针对矩阵不同的特点给出不同的目前的方法是针对矩阵不同的特点给出不同的有效方法有效方法. ( )ijAa

2、AfIAn n nn n矩矩阵阵的的特特征征值值是是 的的特特征征多多项项式式的的 个个零零点点. .计算方法第四章矩阵特征值和特征向量的计算2 常用解法常用解法1 1、 乘幂法和反幂法乘幂法和反幂法2 2、求实对称矩阵特征值的雅可比方法、求实对称矩阵特征值的雅可比方法3 3、求矩阵全部特征值的、求矩阵全部特征值的QR方法方法计算方法第四章矩阵特征值和特征向量的计算3 4.1 4.1 乘幂法和反幂法乘幂法和反幂法一、乘幂法一、乘幂法 乘幂法主要是用来求矩阵的按模最大的特征值与乘幂法主要是用来求矩阵的按模最大的特征值与相应的特征向量。它是通过迭代产生向量序列，由此计相应的特征向量。它是通过迭代产

3、生向量序列，由此计算特征值和特征向量。算特征值和特征向量。设 n 阶实矩阵n nAR的 n 个特征值为(1,2, )iin， 满足120n，所对应的 n 个特征向量 01kkuuAu任取非零的初始向量 ，构造向量序列12, , nxxx线线性性无无关关。 1-1 kkjkjuuu 向向量量逼逼近近A A的的主主特特征征值值（按按模模最最大大的的）对对应应的的特特征征向向量量，计算方法第四章矩阵特征值和特征向量的计算4 11211111111210,(2,3, ) lim()0 lim()0 () inkkiiiiiikkinkkkikiiiinxxkuxxx 由由得得故故只只要要 充充分分大大

4、，就就有有1011,2,iniiiinx存在不全为零的常数 （）,（这里假设0），使得u 101111121()()knnkkkkiiiiiiinkkiiiiuAuA uAxxxx 计算方法第四章矩阵特征值和特征向量的计算5因此，可把 作为与 相应的特征向量的近似，1111 kkkuux 11211111121()()nkiiijjkjinkikjiijjixxuuxx 有有 由111lim()0 =lim,(j1,2,) kjkikkkjunu 11kjkjuu 按上面式子计算矩阵按模最大的特征值与相应的特征向量的方法称为乘幂法。 乘幂法的收敛速度依赖于比值，比值越小，收敛越快。21 A 计

5、算方法第四章矩阵特征值和特征向量的计算61max1max1max1 max(),lim(),lim()( )kkkkkkkkkkkuAVuuuxVuVux，其中是 的绝对值最大的分量 几点说明：）如果 的选取恰恰使得乘幂法计算仍能进行。因为计算过程中舍入误差的影响，迭代若干次后，必然会产生一个向量它在 方向上的分量不为零，可以说实际中出现的可能性几乎为零。0111 10,0kuux 3）重根情形乘幂法也可以计算。）因计算过程中可能会出现溢出或成为的情形。解决方法：每次迭代所求的向量都要规范化。因此，乘幂法实际使用的计算公式是11111 2,(1)0(1)kkux 计算方法第四章矩阵特征值和特征

6、向量的计算7二、乘幂法的加速 因为乘幂法的收敛速度是线性的，而且依赖于因为乘幂法的收敛速度是线性的，而且依赖于比值比值 ,当比值接近于当比值接近于1时，乘幂法收敛很慢。时，乘幂法收敛很慢。乘幂法加速有多种，重点介绍原点平移法乘幂法加速有多种，重点介绍原点平移法。21 11211 12211 ()() ()()iikkkkkkknnnAApIApApIApIuAuuApI upppxxxpp矩阵 与的特征值有以下关系：若是的特征值，则就是的特征值，而且相应的特征向量不变。如果用矩阵按计算，则有 计算方法第四章矩阵特征值和特征向量的计算81211 122111211()() ()() (2,3,)

7、kkkkknnniiuApI upppxxxppppppinp 适当地选取 ，使得且 这样，用乘幂法计算的最大模特征值及相应特征向量的收敛速度比用乘幂法计算要快。这种加速收敛的方法称为原点平移法。1ApIpA 计算方法第四章矩阵特征值和特征向量的计算9123123221222221121121 00)1(), 2 (2,3,2 ) nnninnnnnnppAppppinpp 原原点点平平移移法法使使用用简简便便，但但 的的选选取取困困难难。在在一一些些简简单单情情形形， 可可估估计计。如如当当矩矩阵阵 的的特特征征值值满满足足（或或时时，取取则则有有且且 11 因因此此，用用原原点点平平移移法

8、法求求可可使使收收敛敛速速度度加加快快。计算方法第四章矩阵特征值和特征向量的计算10三、反幂法三、反幂法 反幂法是计算矩阵按模最小的特征值及特征向反幂法是计算矩阵按模最小的特征值及特征向量的方法，也是修正特征值、求相应特征向量的最量的方法，也是修正特征值、求相应特征向量的最有效的方法。有效的方法。设为非奇异矩阵，、为的特征值与相应的特征向量，即此式表明，的特征值是的特征值的倒数，而相应的特征向量不变。因此，若对矩阵用乘幂法，即可计算出的按模最大的特征值，其倒数恰为的按模最小的特征值。 这就是反幂法的基本思想。1111 1 iiiiiiiiAnnxAAxxAxxAAAAA 计算方法第四章矩阵特征

9、值和特征向量的计算111111-11-1kkkkkkkAuA uAAAuuuA uuA对矩阵用乘幂法得， 因为的计算比较麻烦，而且往往不能保持矩阵的一些好性质（如稀疏性），因此，反幂法在实际计算时以求解方程组 ，代替迭代求得，每迭代一次要解一线性方程组。 由于矩阵在迭代过程中不变，故可对先进行三角分解，每次迭代只要解两个三角形方程组。反幂法规范后的计算格式反幂法规范后的计算格式1max()kkkkkkkAuVCuVuC limm ax()nkknxVx 1limkknC 计算方法第四章矩阵特征值和特征向量的计算12 0 () ()iiijjiiijjiiAAIAI 用用带带原原点点平平移移的的

10、反反幂幂法法来来修修正正特特征征值值，并并求求相相应应的的特特征征向向量量是是非非常常有有效效的的。设设已已知知的的一一个个特特征征值值的的近近似似值值为为，因因接接近近，一一般般有有故故是是矩矩阵阵的的按按模模最最小小的的特特征征值值，且且由由上上式式可可知知，比比值值较较小小。因因此此，对对用用反反幂幂法法求求一一般般收收敛敛很很快快，通通常常只只要要经经过过二二、三三次次迭迭代代就就能能达达到到较较高高的的精精度度。四、利用原点平移的反幂法求任一特征值和特征向量四、利用原点平移的反幂法求任一特征值和特征向量计算方法第四章矩阵特征值和特征向量的计算13 4.2 雅可比（雅可比（ Jacob

11、i ）方法）方法 Jacobi方法是用来求实对称矩阵的全部特征方法是用来求实对称矩阵的全部特征值和对应特征向量的一个古典算法。值和对应特征向量的一个古典算法。Jacobi方法方法的基本思想是对矩阵的基本思想是对矩阵A做一系列的正交相似变换，做一系列的正交相似变换，使其非对角元素收敛到零，从而使该矩阵近似为使其非对角元素收敛到零，从而使该矩阵近似为对角矩阵，得到全部特征值和特征向量。对角矩阵，得到全部特征值和特征向量。计算方法第四章矩阵特征值和特征向量的计算14一、古典雅可比方法一、古典雅可比方法11122122cossinsincosaaARaa先考虑二维问题，A为实对称矩阵，R为正交矩阵11

12、12-1T12122ccos ,sincscs=scscsaaARARRARaa做正交相似变换，记222211122212221122221222111112222()()=()()2c acsas acsacs aacsacs aas acsac a2212221122211212()()=010Tcsacs aaaasscac如果，则RAR 成为对角阵，上述方程等价于 计算方法第四章矩阵特征值和特征向量的计算152212221122211212()()=010Tcsacs aaaasscac如果，则RAR 成为对角阵，上述方程等价于 2221112sin10cosaastttca令， =代

13、入则有 2,11tcsc tt 解二次方程可以求出 容易得到1TT111122Tcs0sc0aARARRaR（ ）（ ），的两个列向量是相应的特征向量。计算方法第四章矩阵特征值和特征向量的计算160102021211例 考虑三阶矩阵 A000.70700.7070100.70700.707将A 中（3,1）和（1,3）位置上的元素变成0，取RT10030.7070=0.70720.70700.7071AR AR做正交相似变换后得到计算方法第四章矩阵特征值和特征向量的计算17110.8880.46000.4600.8880001 再将A中（2,1）和（1,2）位置上的元素变成0，取RT2113.

14、36600.325=01.6340.6280.3250.6281AR AR做正交相似变换后得到2210000.9750.22600.2260.975再将A 中（2,3）和（3,2）位置上的元素变成0，取R计算方法第四章矩阵特征值和特征向量的计算18T3223.3660.07350.317=0.07351.78000.31701.145AR AR做正交相似变换后得到kk雅可比方法是一个迭代过程，它生成的是一个矩阵的序列A，当k越大时A 就越接近于对角矩阵，从而得到A特征值更好的近似。计算方法第四章矩阵特征值和特征向量的计算19 11cossin1R1sincos11 ()pqqppq定定义义n

15、n阶阶正正交交矩矩阵阵计算方法第四章矩阵特征值和特征向量的计算20 T(1)1ccos ,sin=ijn nsARARa记，做正交相似变换，得(1)(1)(1)(1)22(1)22(1)22( ,)(,)(,)22()()ijijpjpjqjqjpjqjpppppqqqqqpppqqqpqpqqqppaai jp qacasajp qasacajp qac acsas aac acsas aacsacs aa 22122()()=010pqqqpppqqqpppqcsacs aaaaasscac（ ）令，则=0，上述方程等价于 计算方法第四章矩阵特征值和特征向量的计算212sin10cosqq

16、pppqaastttca令， =代入则有 2,11tcsc tt 解二次方程可以求出 容易得到22122()()=010pqqqpppqqqpppqcsacs aaaaasscac令，则=0，上述方程等价于 计算方法第四章矩阵特征值和特征向量的计算2221()( )2pqD AD Aa2(1)22112nniiiipqiiaaa1TARARJacobi算法的基本思想算法的基本思想: +12112TTTkkkRR RR RRA=A12TnddR Rd=A212D(A)S(A)niiiijijaa记=（对角元的平方和）=（非对角元的平方和）T21kRRR RR记，的每个列向量是对应的特征向量。21

17、()( )2pqS AS Aa计算方法第四章矩阵特征值和特征向量的计算23二、二、 雅可比过关法雅可比过关法（1）循环）循环Jacobi方法：方法： （2）Jacobi过关法：过关法： (,)(1,2)(1,3),(),(2,3),(2,),(),(,)()2p qnnnnp qn n选择依次为1,1,对每个做雅可比变换，完成上述1次变换称做了一次扫描。( )S A 逐次扫描，直到为止。( )11( )( )qkkkpanJacobS AAinS 如果小于某个小参数，就可以让它过关，不做使它零化的变换了，这样可以节约些计算量。每次扫描所取的小参数是一个界限，一般第一个界限取为=，其余取为=。直

18、到这为称止。为过关法。计算方法第四章矩阵特征值和特征向量的计算24多项式运算函数r=roots(p)：求多项式的零点求多项式的零点p=poly(r) : 以以r为零点的多项式为零点的多项式p=poly(A): A的特征多项式的特征多项式PA=polyval(p,S):按数组运算规则，计算多项式的值按数组运算规则，计算多项式的值 其中其中S，PA为矩阵为矩阵PM=polyvalm(p,S):按矩阵运算规则，计算多项式的值按矩阵运算规则，计算多项式的值, 其中其中S，PM为矩阵为矩阵 p=conv(p1,p2)：多项式的乘积多项式的乘积q,r=deconv(p1,p2)：多项式的除法，多项式的除法

19、，p1/p2 p1(x)=p2(x)q(x)+r(x)计算方法第四章矩阵特征值和特征向量的计算25【例】由给定根向量求多项式系数向量。【例】由给定根向量求多项式系数向量。R=-0.5,-0.3+0.4*i,-0.3-0.4*i;P=poly(R)PPR=poly2str(P,x) P = 1.0000 1.1000 0.5500 0.1250PPR = x3 + 1.1 x2 + 0.55 x + 0.125 计算方法第四章矩阵特征值和特征向量的计算26【例】求多项式【例】求多项式 的零点的零点。 r=roots(1 -6 15 -20 15 -6 1)r = 1.0042 + 0.0025i

20、 1.0042 - 0.0025i 1.0000 + 0.0049i 1.0000 - 0.0049i 0.9958 + 0.0024i 0.9958 - 0.0024i665432(1)615201561yxxxxxxx注：尽管利用注：尽管利用MATLAB使得从系数转换到零点或从零使得从系数转换到零点或从零点转换到系数都非常容易，但是使用时一定要注意计点转换到系数都非常容易，但是使用时一定要注意计算的精度。如果存在重根，这种转换可能会降低精度。算的精度。如果存在重根，这种转换可能会降低精度。对于数值计算，计算重根是最困难的问题之一。对于数值计算，计算重根是最困难的问题之一。计算方法第四章矩阵

21、特征值和特征向量的计算27【例】求【例】求3阶方阵阶方阵A的特征多项式。的特征多项式。 A=11 12 13;14 15 16;17 18 19;PA=poly(A) PPA=poly2str(PA,x) PA = 1.0000 -45.0000 -18.0000 -0.0000PPA = x3 - 45 x2 - 18 x - 2.8387e-015 计算方法第四章矩阵特征值和特征向量的计算28【例】求【例】求 的的“商商”及及“余余”多项式。多项式。 p1=conv(1,0,2,conv(1,4,1,1);p2=1 0 1 1;q,r=deconv(p1,p2);cq=商多项式为商多项式为

22、 ; cr=余多项式为余多项式为 ;disp(cq,poly2str(q,x)disp(cr,poly2str(r,x) 商多项式为商多项式为 x + 5余多项式为余多项式为 5 x2 + 4 x + 3 1) 1)(4)(2(32sssss计算方法第四章矩阵特征值和特征向量的计算29n dot(x,y)dot(x,y) 向量的内积向量的内积nnorm norm ： 矩阵或向量范数矩阵或向量范数n det(A) det(A) 方阵的行列式；方阵的行列式；nrank(A) rank(A) 矩阵的秩；矩阵的秩；ntrace(A) trace(A) 矩阵的迹；矩阵的迹；nrref(A) rref(A

23、) 初等变换化矩阵初等变换化矩阵A A为阶梯矩阵为阶梯矩阵ninv(A) inv(A) 矩阵的矩阵的逆逆；即；即 A-1npinv(A) pinv(A) 矩阵的矩阵的广义逆广义逆A A+ +north(A A) 将将A A标准正交化标准正交化n cond(A,flag) cond(A,flag) 矩阵的条件数矩阵的条件数, , flag=2, 1, inf, fro; flag=2, 1, inf, fro;线性代数常用函数计算方法第四章矩阵特征值和特征向量的计算30nd=eig(A) d=eig(A) ： 方阵的特征值；方阵的特征值； V,D=eig(A) : A*V=V*Dn c=cond

24、eig(A) ： 向量向量c c中包含矩阵中包含矩阵A A关于各关于各 特征值的条件数特征值的条件数n V,D,c=condeig(A)：例例: :nA=1 0 0;1 2 0;1 2 3, A=1 0 0;1 2 0;1 2 3, nd=eig(A), d=eig(A), nV,D=eig(A),V,D=eig(A),nC=condeig(A), V,D,C=condeig(A),C=condeig(A), V,D,C=condeig(A),计算方法第四章矩阵特征值和特征向量的计算31例：例：观察观察7 7阶随机矩阵特征值的分布阶随机矩阵特征值的分布a=rands(7,7) %a=rands

25、(7,7) %产生产生7 7阶随机矩阵阶随机矩阵e=eig(a)e=eig(a)title(title(特征值的分布特征值的分布););plot(real(e),imag(e),o)plot(real(e),imag(e),o)xlabel(xlabel(实轴实轴););ylabel(ylabel(虚轴虚轴););注：注：本例验证了如下定理：实方阵的特征值或为实本例验证了如下定理：实方阵的特征值或为实 数或呈共轭对出现。数或呈共轭对出现。计算方法第四章矩阵特征值和特征向量的计算32例：例：观察正交矩阵的特征值分布观察正交矩阵的特征值分布a=rands(7,7);a=rands(7,7);b=o

26、rth(a); %b=orth(a); %构造一个正交矩阵构造一个正交矩阵theta=0:0.01:2theta=0:0.01:2* *pi;pi;e=eig(b);e=eig(b);plot(real(e),imag(e),rplot(real(e),imag(e),r* *,cos(theta),sin(theta);,cos(theta),sin(theta);axis equalaxis equaltitle(title(正交矩阵特征值的分布正交矩阵特征值的分布);); xlabel( xlabel(实轴实轴););ylabel(ylabel(虚轴虚轴););注：注：本例验证了正交矩阵的特征值分布在复平面的单本例验证了正交矩阵的特征值分布在复平面的单 位圆上。位圆上。