矩阵特征值与特征向量的计算

1、第第5 5章章 矩阵特征值与特征向量的计算矩阵特征值与特征向量的计算 方程(5.1)可改写为： Ax-Ax-x=0, x=0, 或(A A-I I)x x=0 0 因此，n阶方阵A A的特征值是特征方程 det(A A-I I)=0 (5.2) 的n个根：1,2,n .A A的属于i的特征向量是如下齐次线性方程组的解 (A A-iI I)x x=0 0 (5.3) 设A A为n阶矩阵，数和非零向量x x满足： AxAx=x x (5.1)则称为A A的特征值，x x为A A的属于的特征向量。圆盘定理圆盘定理 设矩阵A=(aij)nn ,记复平面上以aii为圆心,以 ri=nijjija1|则

2、(1)A A的任一特征值至少位于其中一个圆盘内; (2)在m个圆盘相互连通(而与其余n-m个圆盘互不连通)的区域内,恰有A A的m个特征值(重特征值按重数记).显然对较大的n，通过解方程(5.2)和(5.3)求A的特征值是不实际的，因此需要采用数值方法。注意：若x x是A A的特征向量，则y y=cx x,c0,仍为特征向量。这是因为： Ay=AAy=A(cx x)= =cAxAx= =c x=x= (cx x)= = y y 为半径的n个圆盘为 Ri=aiiri,i=1,2,n试讨论A A的特征值的分布. 解解 由A A确定的3个圆盘分别为所以,当特征值为实数时， 315 -222 -63

3、2 n 相应的特征向量为 x1, x2, xn, 线性无关, Axi= ixi . 我们要求1 和 x1 . 乘幂法的基本思想是取初始向量v v(0),作迭代 v(k+1) =Av(k), k=0,1,2, 产生迭代序列v(k).显然有： v(k+1) =Av(k)= A2v(k-1)=Ak+1v(0) 由于 x1, x2, xn 线性无关, 从而 v(0) =a1x1+a2x2+anxn 利用 AxAxi= ixi 可知 Akxi= ikxi。 故有 v(k) = Akv(0)= Ak (a1x1+a2x2+anxn ) =a11kx1+a22kx2+annkxn (5.4) =a1Akx1

4、+a2Akx2+anAkxn (5.4)式可写成由于|i/1|r+1n ,这时,(5.4)式可写成)()(1111122111)(nknrkrrrkknraaaaaxxxxxv若a1,a2,ar不全为零, 则对充分大的k有 22111)(rrkkaaaxxxv由于a1x1+a2x2+arxr 是对应1的特征向量, 若仍记为x1 ,则有: v(k) 1kx1 ,故前面的结论仍然成立. 设1=-2,且 |1=|2|3 n ,这时,(5.4)式可写成)()() 1(1133122111)(nknkrkkknaaaaxxxxv则对充分大的k有 v v(k+1)1k+1(a1x x1-(-1)ka2x

5、x2) , v v(k+2)1k+2(a1x x1+(-1)ka2x x2) ) 1(22111)(xxvaakkk于是有nikiki, 2 , 1,/)()2(1vv x x1v v(k+1)+1v v(k) , x x2v v(k+1)-1v v(k)对于规范化的幂法，由于 v v(k+1)=AuAu(k), u u(k+1)=v v(k+1)/k+1,则 u u(k+2)=v v(k+2)/k+2=Au(k+1)/k+2 =AvAv(k+1)/k+1k+2=A A2u u(k)/k+1k+2于是有，211kk12 x x1k+1u u(k+1)+1 ku u(k) , x x2k+1u

6、u(k+1)-1 ku u(k)的按模最大特征值和相应的特征向量。例例4 用乘幂法求矩 解解 取初始向量u u(0)=(1,1,2)T ,计算可得K ku(k)012345678910111213113.5536284.6792043.4611244.6354653.4526554.6321163.4543154.6319293.4542914.6319203.4542884.631924(1,1,2)T(0.454545, 0.909091, 1)T(0.537222, 0.972116, 1)T(0.465201, 0.994041, 1)T(0.539392,

7、0.998269, 1)T(0.465721, 0.999627, 1)T(0.539487, 0.999892, 1)T(0.465890, 0.999975, 1)T(0.539495, 0.999993, 1)T(0.465893, 0.999999, 1)T(0.539495, 1, 1)T(0.465893, 1, 1)T(0.539495, 1, 1)T(0.465893, 1, 1)T1.2 加速技术加速技术由于)5 .5(),()max()max(121)1()(kkkkoAuv所以,乘幂法收敛速度取决于比值|2/1|,当|2/1|1时,收敛是很慢的. 1 1.Aitken A

8、itken 加速方法加速方法由(5.5)式可知 x x2=13u u(13)-112u u(12)=(-5.2963,-9.1852,-9.1852 )T.4, 4631924. 4454288. 3213121 x x1=13u u(13)+112u u(12)=( 9.6122,18.4490,18.4490 )T, 实际上, , A A的特征值为1=4,2=-4,3=1.可见,序列k线性收敛于1 .会达到加速收敛的目的.0lim12111kkk 构造Aitken序列kkkkkkk12212)( 如把Aitken加速方法用于例3,则有 ku(k)107.26.56.0033526.0016

9、756.000837(1,1,1)T(1,0.8,0.1)T(1,0.75,-0.111)T(1,0.730769,-0.188034)T.(1,0.714405,-0.249579)T(1,0.714346,-0.249790)T(1,0.714316,-0.249895)Tk01231011126.2666676.0000176.0000036.000000k 2 2.原点位移法原点位移法 取定数p, 作矩阵B B=A A-pI I, 则 B B的特征值为mi=i-p, (i=1,2,n),而且对应的特征向量相同.求1 求m1。则对B B应用乘幂法可达到加速收敛的目的。 解解 由于A A的

10、特征值为1=6,2=3,3=2,故取p=2.5,则B B的特征值为m1=3.5,m2=0.5,m3=-0.5,则如果选取p,使m1仍然是B B的按模最大特征值，且满足取初始向量u u(0)(0)=(1,1,1)T,由规范化计算公式:121212ppmm例例5 用原点位移法求例3中矩阵A A的按模最大的特征值和特征向量.5.00105.1050145.65.2IAB计算可得, 3 , 2 , 1,/)()(kkkkvu)max()(kkv)1()(kkBuvkku(k)01234567.53.7666623.5353963.5050023.5007183.500102(1,1,1,)T(1,0.

11、733333,-0.2)T(1,0.716814,-0.238938)T(1,0.714643,-0.249061)T(1,0.714337,-0.249777)T(1,0.714293,-0.249981)T(1,0.714287,-0.249995)T 反幂法是求矩阵按模最小的特征值和相应特征向量的方法.取,1=m1+p6+2.5=6.000102,x x1u u(6)=(1,0.714287,0.249995)T ,比例3中收敛的快.这是因为|2/1|=1/2,而|m2/m1| =1/7,故对B B应用乘幂法远比对A A应用乘幂法收敛的快.1.3 反幂法反幂法 设A是n阶非奇异矩阵, 其

12、特征值为 1 , 2 , ,n-1 ,n, |1| |2| |n-1| |n| 0对应的特征向量为x1,x2,xn, 则A-1的特征值为1211111nnn对应的特征向量为x1,x2,xn. 要想求n和xn只需对A-1应用乘幂法。任取初始向量u(0)0, 进行迭代计算 ,1,1,121n也可将上式改写成, 3 , 2 , 1,/)()(kkkkvu)max()(kkv)1(1)(kkuAv)1()(kkuAv)max()(kkv, 3 , 2 , 1,/)6 . 5()()(kkkkvu式(5.6) 称为反幂法. 显然有)max(/lim,1lim,1lim)(nnkkkknnkkxxu每一步

13、求v(k)需要求解线性方程组, 可采用LU分解法求解.反幂法还可与原点位移法结合求A A的更精确的特征值。设已知矩阵A A的特征值i的某个近似值 ，作原点位移,令B B=A A- I I, 对B B应用反幂法可求出 设已求得例3中矩阵A A的特征值的近似值16.003,和相应的特征向量x x1(1,0.714405,-0.249579)T, 试用带原点位移的反幂法求1和x x1的更精确的值.iip,取i)( ,221ijijiiiniii且必有例例6 解解 取p=6.003, 作矩阵B B=A A-6.003I I,则003. 4010997. 65014003.10B则B B的特征值为)/(

14、1iik,得到精度更高的i.取初始向量u u(0)=(1,0.714405,-0.249579)T,对B B用反幂法计算可得:可见收敛速度非常快,这是因为B B的3个特征值为m1=-4.003, m2=-3.003,m3=-0.003,|m3/m2|0.000999很小. Jacobi方法是求实对称矩阵全部特征值和特征向量的一种矩阵变换方法。T)250000. 0,714286. 0, 1 (,00000167. 6003. 61)1(1uT)250000. 0,714286. 0, 1 (,000000007. 6003. 61)2(2u2 Jacobi 方法方法 实对称矩阵A A具有下列性

15、质： （1）A A的特征值均为实数； （2）存在正交矩阵R R，使R RTARAR=diag(1,2,n),而 R的第i个列向量恰为i的特征向量; 直接求正交矩阵R是困难的 . Jacobi提出用一系列所谓平面旋转矩阵逐次将A约化为对角矩阵.平面解析几何中的平面坐标旋转变换表示平面上坐标轴旋转角的变换. （3）若记A A1=R=RTARAR,则A A1仍为对称矩阵. 2.1 平面旋转矩平面旋转矩阵阵 1122cossinsincosyxyx 在三维空间直角坐标系中,ox1y1平面绕着oz1轴旋转角的坐标变换为1112221000cossin0sincoszyxzyx 一般地, 在n维向量空间R

16、n中, 沿着xpyq平面旋转角的变换矩阵为行第行第qppq11cossin11sincos11)(R称Rpq()为平面旋转矩阵平面旋转矩阵. Rpq()具有下列性质: 设实对称矩阵A A=(aij)nn ,记B B=R RpqT()ARARpq()=(bij)nn则它们元素之间有如下关系: （1）R Rpq()为正交矩阵,即R Rpq-1()=R RpqT(); （2）如果A A为对称矩阵, 则R RpqT()ARARpq()也为对称矩阵, 且与A A有相同的特征值. （3）R RpqT()A A仅改变A A的第p行与第q行元素,ARARpq()仅改变A A的第p列与第q列元素.)4 . 5(

17、),(cossinsincos2cos2sin)(2sincossin2sinsincos212222qpjiabaabbaabbaaabbaaabaaabijijjqjpqjjqjqjppjjppqppqqqppqpqqqppqqpqqqpppp所以有从而2222qipiqipiaabb),(2222qpiaabbiqipiqip)5 . 5(22222222pqqqpppqqqppaaabbb22FFAB)6 . 5(,1212njiijnjiijab，即有(5.5)、(5.6)式可得222222pqjiijpqjiijaabb21221222pqniiipqniiiaabb 如果apq0

18、, 适当选取角, 使02cos2sin)(21pqppqqqppqaaabb只需角满足从而 如果取|apq|=若记)7 . 5(4|,22pqqqppaaactgjiijpqjiijjiijaaab22222niiipqniiiniiiaaab12212122|maxijjia于是jiijpqanna2211）（，则jiijjiijannb22) 1(21 (,)(2jiijaA 则上式可记为)8 . 5()() 1(21 ()(ABnn 由式(5.7),令t=tg,则t满足方程 t2+2t-1=0 经典Jacobi算法是对A A(0)=A A施行一系列平面旋转变换:为保证|/4,取绝对值较小

19、的根,有于是)9 . 5(0,10,)1|/(|)(2signt,)1 (cos212t)10. 5(cossint2.2 Jacobi 方法方法 A A(1)=R R1TA A(0)R R1 ,A A(2)=R R2TA A(1)R R2 , A A(k)=R RkTA A(k-1)R Rk ,每一步变换选择A A(k-1)=(aij(k-1)nn 的非对角线元素中绝对值最大者apq(k-1)(称为主元素)作为歼灭对象, 构造平面旋转矩阵R Rk=R Rpq(), 经变换得到A A(k)=(aij(k)nn ,且apq(k)=0,这时由(5.8)式有从而由此递推得到 当k充分大时,或者(A

20、A(k),或者)() 1(21 ()()1()(kknnAA)() 1(21 ()()(AAkknn)( , 0k),(lim21)(nkkdiagDA 另外,由于 A A(k)=R RkTA A(k-1)R Rk=R RkTR Rk-1TR R1TARAR1R R2R Rk=R RTARAR 是给定的精度要求,则A A的特征值可取为iaii(k),i=1,2,n.,|max)(kijjia的全部特征值. 解解 记 A A(0)=A,A,取p=1,q=2,apq(0)=a12(0)=2,于是有因此,R=RR=R1R R2R Rk 的列向量x xj (j=1,2,n)为A A的近似特征向量.例例

21、7 用Jacobi 方法计算对称矩阵612152224A25. 02)0(12)0(22)0(11aaa780776. 0)1|/(|)(,2signt788206. 0)1 (cos212t615412. 0cossin,t从而有所以 再取p=2,q=3,apq(1)=a23(1)=2.020190,类似地可得以下依次有1000788206. 0615412. 00615412. 0788206. 01000cossin0sincos)(1pqRR6020190. 2961. 0020190. 2561552. 60961. 00438448. 21)0(1)2(RARAT241166. 40724794. 00320386. 8631026. 0724794. 0631026. 0438448. 2)2(A496424. 4209614. 00209614. 0320386. 8595192. 00595192. 0183185. 2)3(A496424.