《特征值和Jacobi方法.ppt》由会员分享,可在线阅读,更多相关《特征值和Jacobi方法.ppt(49页珍藏版)》请在三一办公上搜索。
1、第9章 矩阵特征值问题的数值方法,9.1 特征值与特征向量9.2 Hermite矩阵特征值问题9.3 Jacobi方法9.4 对分法9.5 乘幂法9.6 反幂法9.7 QR方法,引言,工程实践中有多种振动问题,如桥梁或建筑物的振动,机械机件、飞机机翼的振动,工程实践中有多种振动问题,如桥梁或建筑物的振动,机械机件、飞机机翼的振动,及一些稳定性分析和相关分析可转化为求矩阵特征值与特征向量的问题。,London,England:Millennium(Wobbly)Bridge(1998-2002,Norman Foster and Partners and Arup Associates),I d
2、ecide that I have to write something today,otherwise I would not know how to speak English here.This is a very quick story about a bridge.London launched three major construction projects to celebrate the arrival of the Millennium.After all,Greenwich(pronounced green-ich)is supposed to be(supposed t
3、o be?!)where the prime meridian lies,and the place where the Millennium officially starts in the world.The three projects are the Millennium Dome in North Greenwich,so far the largest single roofed structure in the world,London Eye right across Westminster,which becomes so far the largest observatio
4、n wheel in the world,and the Millennium Bridge that links Southeast London with St.Pauls Cathedral,which is currentlywell.not swinging any more,it is said.,The bridge was designed by Imperial College,a college of my former university.On the very first day that the bridge was open to public,there wer
5、e simply so many people going there to walk from the south bank to St.Pauls that the weight completely exceeded the architects expectation.The slender steel truss bridge began to vibrate with a million people on there.The opening ceremony ended up in an embarrassing vertigo.Millennium left Londoners
6、 a happy adage about swinging bridge,meaning fancy technology that looks good but functions in a funny fashion.Am I using too many Fs here?Or is it simply because my tongue starts to swing in the same direction when I am writing about this wobbly bridge?Next time you visit London,I strongly recommen
7、d this place.After all,with a little swing,this is a shortcut to dash into St.Pauls directly from the southeast!,G:Google Matrix,“the worlds largest matrix computation”.4,300,000,000 x:PageRank(网页级别)vector“The$25,000,000,000 Eigenvector”,搜索引擎,9.1 特征值与特征向量,设A是n阶矩阵,x是非零列向量.如果有数存在,满足,(1)那么,称x是矩阵A关于特征值的
8、特征向量.,如果把(1)式右端写为,那么(1)式又可写为:,记,它是关于参数的n次多项式,称为矩阵A的特征多项式,其中a0=(-1)nA.,(2),显然,当是A的一个特征值时,它必然是 的根.反之,如果是 的根,那么齐次方程组(2)有非零解向量x,使(1)式成立.从而,是A的一个特征值.A的特征值也称为A的特征根.,矩阵特征值和特征向量有如下主要性质:,定理9.1.1 n阶矩阵A是降秩矩阵的充分必要条件是A有零特征值.,定理9.1.2 设矩阵A与矩阵B相似,那么它们 有相同的特征值.,定理9.1.3 n阶矩阵A与AT有相同的特征值.,定理9.1.4 设ij是n阶矩阵A的两个互异特征值,x、y分
9、别是其相应的右特征向量和左特征向量,那么,xTy=0.,9.2 Hermite矩阵特征值问题,设A为n阶矩阵,其共轭转置矩阵记为AH.如果A=AH,那么,A称为Hermite矩阵.,9.2.1 Hermite矩阵的有关性质 设 是Hermite矩阵A的n个特征值.有以下性质:,全是实数.,有相应的n个线性无关的特征向量,它们可以化为一组标准酉交的特征向量组,即,是酉空间中的一组标准酉交基.,记U=(),它是一个酉阵,即UHU=UUH=I,那么即A与以 为对角元的对角阵相似.,A为正定矩阵的充分必要条件是 全为正数.,定理9.2.1 设 是Hermite矩阵A的n个特征值,那么,证:,设x是一个
10、非零向量,A是Hermite矩阵,称 为矩阵A关于向量x的Rayleigh商,记为R(x).,定理9.2.2 如果A的n个特征值为 其相应的标准酉交的特征向量为 那么有,定理9.2.3 设A是Hermite矩阵,那么,9.2.2 极值定理,定理9.2.4(极值定理)设Hermite矩阵的n个特征值为,其相应的标准酉交特征向量为.用Ck表示酉空间Cn中任意的k维子空间,那么,9.2.3 Hermite矩阵特征值问题的性态,矩阵特征值问题与求解线性方程组问题一样,都存在当矩阵A的原始 数据有小变化(小扰动)时,引起特征值问题的变化有大有小的问题,如果引起的变化小,称 该特征值问题是良态的.反之,称
11、为病态的.矩阵特征值问题的性态是很复杂的,通常分别就单个特征值或整体特征值给出条件数进行分 析.对于Hermite矩阵,由于其特征值问题的特殊性质,其特征值都是良态的.下面先证明Hermite矩阵特征值的扰动定理.,定理9.2.5 设矩阵A,E,A+E都是n阶Hermite矩阵,其特征值分别为 那么,证 设矩阵A关于特征值1,2,n 的标准酉交特征向量为u1,u2,un,是由ui,ui+1,un生成的n-i+1维子空间.对 中任意非零向量x,由极值定理,有,由定理,又由定理,对任意x0,有从而有另一方面,A=(A+E)-E.记 为矩阵-E的特征值,那么,重复上面的过程,可得从而有,定理通常又称
12、为Hermite矩阵特征值的扰动定理,定理9.2.6 设矩阵A和A=A+E都是n阶Hermite矩 阵,其特征值分别为 和,那么,这个定理表明,扰动矩阵E使A的特征值的变化不会超过 E2.一般E2小,因此,Hermite矩阵特征值是良态的.,9.3 Jacobi方法,理论上,实对称矩阵A正交相似于以A的特征值为对角元 的 对角阵.问题是如何构造这样的正交矩阵呢?Jacobi方法就是通过构造特殊的正交矩阵 序列,通过相似变换使A的非对角线元素逐次零化来实现对角化的.,9.3.1 平面旋转矩阵与相似约化先看一个简单的例子.,设 是二阶实对称矩阵,即a21=a12,其特征值为1,2.令 使得 记 容
13、易验证BT=B,且,解之得:当 时当 时可选取,为使RTAR为对角阵,要求b12=b21=0,一般的n阶平面旋转矩阵,9.3.2 经典的Jacobi方法,设A是实对称矩阵,记A1=A.Jacobi方法的基本思想是用迭代格式 Ak+1=QTkAkQk,k=1,2,构造一个相似矩阵序列,使Ak收敛于一个对角阵.其中 Qk为平面旋转矩阵,其旋转角k由使Ak的绝对值 最大元a(k)pq=a(k)qp=0 或按列依次使A的非对角元 零化来确定.,定理9.3.1 设A是n阶实对称矩阵,那么由Jacobi方法产生的相似矩阵序列Ak的非对角元收敛于0.也就是说,Ak收敛于以A的特征值为对角元的对角阵.,记 其
14、中Ek是Ak除主对角元外的矩阵.由平面旋转矩阵的性质 中,对于,有因此,又由假设,因此,这样,便有从而,当,9.3.3 实用的Jacobi方法,循环Jacobi方法必须一次又一次扫描,才能使Ak收敛于对角阵,计算量很大.在实际计算中,往往用一些特殊方法来控制扫描次数,减少计算量.下面介 绍一种应用最为广泛的特殊循环Jacobi方法阈Jacobi方法.阈Jacobi方法首先确定一个阈值,在对非对角元零化的一次扫描中,只对其中绝对值 超过阈值的非对角元进行零化.当所有非对角元素的绝对值都不超过阈值后,将阈值减少,再重复下一轮扫描,直至阈值充分小为止.减少阈值的方法通常是先固定一个正整数Mn,扫描一
15、次后,让.而阈值的下界是根据实际问题的精度要求选定的.,9.3.4 用Jacobi方法计算特征向量,假定经过k次迭代得到Ak+1=RTkRT1AR1Rk,(15)这时Ak+1是满足精度要求的一个近似的对角阵.如果记Qk=R1R2Rk=Qk-1Rk,(16)那么,Qk是一个正交矩阵,且(15)式又可表示为Ak+1=QTkAQk.当Ak+1的非对角元素充分小,Qk的第 j列qj可以看成是近似特征值a(k+1)jj相应的特征向量了.,在实际计算中,可以按(16)式在迭代过程中形成Qk,把Qk看成是Qk-1右乘一个平面旋转矩阵得到.不妨记 Q0=I,Qk的元素按下式计算:,9.4 对分法,理论上,一个
16、实对称矩阵正交相似于一个以其特征值为对角元的 对角阵.但是,经典的结果告诉我们,一个大于4次的多项式方程不可能用有限次四则运算 求根.因此,我们不可能期望只用有限次相似变换将一个实对称矩阵约化为一个对角阵.下面先介绍将一个实对称矩阵相似约化为实对称三对角矩阵的方法,再讨论求其特征值的对 分法.,9.4.1 相似约化为实对称三对角矩阵,将一个实对称矩阵正交相似约化为一个实对称三对角矩阵的算法,可归纳如下:记A(1)=A,对k=1,2,n-2 按(4)式、(5)式和(8)式计算;按(9)(12)式,计算A(k+1).,9.4.2 Sturm序列的性质,设实对称三对角矩阵为 其中i0(i=1,2,,
17、n-1)其特征矩阵为T-I.记T-I的第i阶主子式为,这是关于的i次多项式,当i=n时,pn()=T-I是矩阵T的特征多项式.令p0()1,则有p1()=1-,pi()=(i-)pi-1()-2i-1pi-2(),i=2,3,n.(15)多项式序列pi()(i=0,1,,n)称为Sturm序列,定理pi()(i=1,2,,n)的根都是 实根.证 由(14)式,pi()是i阶实对称矩阵的特征多项式,因此,pi()(i=1,2,,n)的根全是实根.,定理定理9.4.2 设是pi()的一个根,那么 pi-1()pi+1()0,即相邻的两个多项式无公共根;pi-1()pi+1()0,即pi-1()与p
18、i+1()反号.,定理9.4.4 pi()的根都是单根,并且将pi+1()的根严格隔离.,9.4.3 同号数和它的应用,定义1 设p0()1,pi()(i=1,2,,n)是一个Sturm序列,称相邻的两个数中符号一致的数目为同号数,记为ai().若某个pi()=0,规定与pi-1()反号.定理9.4.5 设两个实数xy,那么,形如(13)式的实对称 三对角矩阵T的特征多项式在区间(x,y上根的数目为a(x)-a(y).,9.4.4 求Hermite矩阵特征值的对分法,对分法的计算可归纳为以下4个部分 确定(13)式的矩阵T的全部特征值的分布区间.在区间a,b中,用区间对分的方法找出只含T的一个
19、特征值的子区间.在只含一个特征值的子区间上的对分法.同号数的计算.,9.5 乘幂法,设A是n阶矩阵,其n个特征值按模从大到小排序为又假设关于1,2,n的特征向量v1,v2,vn线性无关.,乘幂法是适用于求一般矩阵按模最大特征值及相应特征向量的算法.,xkk1a1v1(k).因此,xk可看成是关于特征值1的近似特征向量.,迭代格式为,按模最大特征值1及其相应的特征向量v1的乘幂法的计算公式:,9.5.2 收缩方法,设矩阵A的n个特征值按模从大到小排序为,其相应的n个线性无关特征向量为v1,v 2,vn.在计算A的最大特征值1及相应特征向量v1后,可以通过收缩方法,继续用乘幂法计算2及其相应的特征
20、向量v2.,定义n阶矩阵把去掉A1的第1行和第1列的n-1阶矩阵记为,那么,B有与A1除1外的相同的n-1个特征值 23n,可以用乘幂法计算2及其相应的 特征向量.在计算1和v后,按(15)式形成n-1阶矩阵B的计算过程称为收缩方法.,9.6 反幂法,反幂法可以求一个非奇异矩阵A的逆矩阵A-1的按模最小的特征值及相应的特征向量,又可以求A的一个近似特征值相应的特征向量.9.6.1 求按模最小特征值及相应特征向量的反幂法,又称为反迭代法.,9.6.2 求近似特征值的特征向量的反幂法,先对矩阵 进行LU分解,记 那么,(7)下面介绍一种选取特殊的初始向量x0的反幂法半迭代法.,假设,选取初始向量x
21、0满足x0=1,这时z0=x0.对照(7)式中的第二个式子.可把z0看成满足Le=z0.(8)这里,e=(1,1,1)T,而z0的各个分量的取值多少是无关重要的.这样,在第一个迭代步的计算中,只需求解(7)式中的上三角方程组Ux1=e.“半迭代法”的命名也由此而得.,9.7 QR方法,定理9.7.1设A是n阶矩阵,其n个 特征值为.那么存在一个酉矩阵U,使UAU是以为 对角元的上三角矩阵.,9.7.1 两个基本定理,定理9.7.2设A是n阶实矩 阵,那么,存在一个正交矩阵Q,使QAQ为一个准上三角矩阵,它的对角元是A的一个特征值,对角元上的二阶块矩阵的两个特征值是A的一对共轭复特征值.,9.7.2 相似约化为上Hessenberg矩阵,对一般n阶矩阵,QR算法的每一个迭代步需要O(n)次乘法运算.如果矩 阵阶数稍大,这个算法几乎没有实际的应用价值.通常采用的方法是先将矩阵相似约化为上Hessenberg形式的矩阵,在此基础上应用QR迭代.这时,一个QR迭代步的乘法运算次数只需O(n)次.,所谓上Hessenberg矩阵是指一个n阶矩阵A,如果当ij+1时,aij=0,称A为上Hessenberg矩阵.例如:一个5阶的上Hessenberg矩阵具有如下的形式:,下面介绍QR方法时,都假设矩阵A是一个上Hessenberg矩阵.,