矩阵特征值与特征向量的计算PPT学习教案

会计学1矩阵特征值与特征向量的计算矩阵特征值与特征向量的计算3.1 乘幂法及其变体3.1.1 乘幂法定理 设A Rnn有完全特征向量系，若1, 2, n为A的n个特征值且满足n21 对任取初始向量x(0) Rn，对乘幂公式)1()(kkAxx确定的迭代序列xk，有下述结论： （1）当 时，对i = 1, 2, , n211)()1(lim kikikxx收敛速度取决于 的程度，r 越小收敛越快，r 1收敛慢，112r且x(k)（当k充分大时）为相应于1的特征向量的近似值。第1页/共69页321（2）当 时a）若1 = 2，则主特征值1及相应特征向量的求法同（1）；(2)21( )limkikkixx收敛速度取决于 的程度。向量 、113r)(1)1(kkxx c）若 ，则连续迭代两次，计算出x(k+1)，x(k+2)，210)()1()2( kjkjkjqxpxx)(1)1(kkxx 分别为主特征值1、2相应的特征向量的近似值。然后对j = 1, 2, , n 解方程b）若1 = -2，对i = 1, 2, , n第2页/共69页求出 、 后，由公式pq2122 pqip 2222 pqip 解出主特征值1、2。此时收敛速度取决于 的程度。113r)(2)1(kkxx)(1)1(kkxx向量 、 分别为相应于1，2的特征向量的近似值。第3页/共69页令max(x)表示向量x分量中绝对值最大者。即如果有某i0，使iniixx 1max0则 max (x) = xi对任取初始向量x(0)，记)max()0()0()0(xxy 则)0()1(Ayx 一般地，若已知x(k)，称公式 ), 1, 0()max()()1()()()(kAyxxxykkkkk第4页/共69页定理 设ARnn具有完全特征向量系，1, 2, , n为An21则对任初始向量x(0)，由规范化的乘幂法公式确定的1)()max(lim kkx(1)（2）y(k)为相应于主特征值1的特征向量近似值的n个特征值，且满足向量序列y(k)，x(k)满足第5页/共69页3.1.2 反幂法若 A 有| 1 | | 2 | | n |，则 A1 有11111 nnA1 的主特征根 A的绝对值最小的特征根)()1(kkxAx )(1)1(kkxAx 如何计算解线性方程组对应同样一组特征向量。设ARnn可逆，则无零特征值，由)0( xxAx 有 xxA 11 第6页/共69页规范化反幂法公式为 ), 1, 0()max()()1()()()(kyAxxxykkkkk为节省工作量，可先对A进行LU分解，再解三角形方程( )(1)(0,1,)kkLxyUxxk第7页/共69页3.1.3 乘幂法的加速（原点位移法） njjkjjkkkvAxx111)1()( ，希望 | 2 / 1 | 越小越好。取0（常数），用矩阵B = A - 0I 来代替A进行乘幂迭代。0 ii (i = 1, 2, , n)iiiiiivvAvvIABv)()(000 设i (i = 1, 2, , n)为矩阵B 的特征值，则B与A特征值之间应有关系式：第8页/共69页关于矩阵B的乘幂公式为 )0(0)0()()(xIAxBxkkk 为加快收敛速度，适当选择参数0，使00210()maxjj n 达到最小值。 0101 1210()knjkjjjvv jkjnjjkvv12111 第9页/共69页当i (i = 1, 2, , n)为实数，且12 n时，取)(212*0n 则为 (0) 的极小值点。这时122122122*01*02221212121 nnnn第10页/共69页若已知矩阵A的特征值i的一个近似值0 ，要求i和对记B = A - 0I，对任取初始向量x(0)Rn， ), 1, 0()max()()1()()()(kyBxxxykkkkk应的特征向量，可考虑利用原点移位加速的反幂法。第11页/共69页3.2 子空间迭代法斯密特（Schmidt）正交化过程： 设1，2，3 为R3上的三个线性无关的向量，令 ，则1为单位长度的向量，再令2111222211222,),( 可以验证（1, 2）= 0，即1与2正交。若令22311333),(),( 则0),(),(2313 第12页/共69页2333 即与1, 2正交，将其单位化为于是向量组1, 2, 3构成R3上一组标准正交基，且232322131221321321),(),(),(,QR其中Q = 1, 2, 3为正交矩阵，R是上三角阵。第13页/共69页对n维向量空间，设1, , n为Rn上n个线性无关的向量，类似有11 2111 11222),( 2222 22311333),(),( 2333 jjnnjnn ),(11 2nnn 第14页/共69页 232322322113122111),(),(),(),(),(),(,nnnnnn 即QR Q为正交阵，R 为上三角阵将n个线性无关向量变换为n个两两正交向量的方法称为 斯密特正交化方法。斯密特正交化过程将可逆阵A分解为正交阵与上三角阵的乘积。第15页/共69页子空间迭代法： 先取p个线性无关的向量 构成一个n*p的初始矩阵 nxxx,21,210nxxxY 再用可逆矩阵A左乘上式，得到01AYZ 进一步用schimidt正交化的方法将Z1分解为111RYZ 其中Y1的各列两两正交，R1为可逆上三角阵。110RYAY 第16页/共69页若已知第k步迭代矩阵为Yk，则子空间迭代法的公式为 11111kkkkkRZYAYZ第17页/共69页3.3 雅克比 (Jacobi) 旋转法 3.3.1 预备知识1）若B是上（或下）三角阵或对角阵，则B的主对角元素即是B的特征值。2）若矩阵P满足PTP = I，则称P为正交矩阵。显然PT = P-1，且P1, P2, 是正交阵时，其乘积P = P1P2Pk仍为正交矩阵。3) 实对称矩阵的特征值均为实数，且存在标准正交的特征向量系。第18页/共69页4）称矩阵 列列列列jiPij111111 cossinsincos为平面旋转矩阵 第19页/共69页2雅克比方法设矩阵ARnn是对称矩阵，记A0 = A，对A作一系列旋转相似变换TCP APP是一个正交阵，我们称它是（i, j）平面上的旋转矩阵 PTAP只改变A的第i行、j行、i列、j列的元素；A和C的元素仅在第i行（列）和第j行（列）不同，第20页/共69页( , )lklkcal ki j其余元素不变2222cossinsin2sincossin21sin2cos22iiiijjijjjiijjijijjiiijjijcaaacaaaccaaacossin,cossinikkiikjkjkkjjkikccaaki jccaa其中A与C的元素有下面的关系第21页/共69页22,1,1( )( )nnijiji ji jijAaAa 22( )( )( )( )22 ijijCACAac2222(, )iljliljlccaali j引入记号易知22222222iijjijiijjijcccaaa 第22页/共69页我们选取Pij，使得 ，因此需使 满足0ijjicc2tan2ijiijjaaa将 限制在下列范围内44 如果 0iijjaa0ija 4 0ija 4 第23页/共69页直接从三角函数关系式计算sin 和cos，记2()iijjijiijjyaaxa sign aa 则tan2xy 当 时，有下面三角恒等式：4tan222212cos1cos212yxy第24页/共69页于是 2221cos2yxy采用下面公式计算 sin tan22sin22sin cos2cos2xxy 第25页/共69页第26页/共69页 首先，在A=A1=(als(1)的非对角元素中选取绝对值最大元素，记为1111111111211122223 3 9=()0,( )( )( ),( . . ),ijTlsijj iijPPAAP A PAaaa 对对确确定定的的，用用公公式式确确定定出出 ，从从而而产产生生平平面面旋旋转转矩矩阵阵约约化化 为为且且的的非非对对角角元元素素( )( ),()1111110ijijl samax a 第27页/共69页然后，再在A2的非对角元素中选取绝对值最大元素()(),(),22222,22222322221112,12112,(3.3.9),1,2,kkijTTTkijTTTkkkijPPAijj iAP A PP P A P PPPkAPP P A P PP 确确定定出出其其位位置置，按按公公式式确确定定 角角，从从而而产产生生平平面面旋旋转转矩矩阵阵并并将将中中位位置置的的元元素素化化为为零零 得得到到如如此此做做下下去去 可可得得到到一一系系列列平平面面旋旋转转矩矩阵阵使使得得第28页/共69页(1)2( )2(1)2111111(3.3.8)()()2()2()()2()22()()()1()(1)(1)21()0 ()(1)kkkkkijijkijkkkkkAAaaAaAAAAn nn nAkn n 由由关关系系式式知知并并有有这表明：当k趋于无穷时，序列收敛且其极限矩阵是一个对角阵。第29页/共69页 Jacobi方法中方法中 ,每变换一次都要在非对角元素中每变换一次都要在非对角元素中扫描选取绝对值最大元素扫描选取绝对值最大元素,这难免增加很大的这难免增加很大的工作量工作量.为减少工作量和提高运行效率为减少工作量和提高运行效率,下面介下面介绍一种改进的方法绍一种改进的方法,即即Jacobi过关法。过关法。 具体做法具体做法: :第30页/共69页0112131232(1)111212( ),( ),.nnnnijijrrAAAnaaaaaaaanAn 1、1、计计算算矩矩阵阵 的的非非对对角角元元素素平平方方和和2、2、预预取取阀阀值值对对 的的非非对对角角元元素素进进行行扫扫描描对对的的元元素素进进行行旋旋转转点点还还 将将化化为为零零 其其余余元元素素视视为为过过关关不不作作变变换换 当当所所有有的的非非对对角角元元素素绝绝对对值值都都小小于于后后, ,缩缩小小阀阀值值,如,如取取重重复复前前面面的的步步骤骤.按.按此此方方法法,阀,阀值值不不断断缩缩小小 直直到到满满足足之之后后停停止止计计算算这这种种方方法法称称为为Jacobi过Jacobi过关关法法. .第31页/共69页第32页/共69页22(, ).cossin0 (, )1/,sinsjkkjjkkjjkikikkjjkikccki jccaaki jaaacoa 选选取取 使使得得某某个个变变为为零零只只要要取取, , , ,2,3,12,4,12, ,13,4,23,5,23, ,21, ,2, ,1,0(, ),;,;,2,3,1,1,2,.Ti j ki j ki j kjkkjnnnn ni j iPPCPAPccki jPPPPPPPinPjiinAAC 记记上上述述旋旋转转矩矩阵阵为为则则相相似似变变换换使使元元素素取取或或对对用用一一次次对对 进进行行正正交交相相似似是是变变换换, ,便便可可将将 化化为为三三对对角角矩矩阵阵第33页/共69页1、旋转变换设 为实对称阵()ijn nAa1cossinsincos1ijiPjij行行列列 利用旋转矩阵Pij，复习TijijCP AP对A作旋转变换，即得第34页/共69页则A与C的元素之间有一定关系：复习TijijCP AP1sin2cos22ijjiiijjijccaaacossin(, )jkkjjkikccaaki j Jacobi旋转法0ijjicc选取 ，使得第35页/共69页cossin(, )jkkjjkikccaaki j0jkkjcc选取 ，使得记得到的旋转阵为Pi,j,k，则相似变换, , ,Ti j ki j kCPAP0(, ).jkkjccki j使得C中的元素取22,ikjkdaacos,sin,ikjkadad 若d 0，则取否则，取0.第36页/共69页*用旋转变换化实对称矩阵为三对角阵的具体步骤：A*0*0*00*0*0*000*0*0*0*T2,3,12,3,1PPT2,4,12,4,1PPT2, ,12, ,1nnPPStep1 约化第一行第一列：1A 1ATT3, ,23,4,23,4,23, ,2nnPPPPStep2 约化第二行第二列：2A *000*000*00*00*第37页/共69页Step k 约化第k 行第k列：TT1, ,1,2,11,2,1, ,kkn kkkkkkkkkn kAPPAPPStep n-2 约化第n-2行第n-2列：T21, ,21, ,2nnn nnn nAPP*0000*0000*00*000*0000*0000*00*0*0000*3nA用旋转变换化实对称矩阵为三对角阵的具体步骤：第38页/共69页例3.4利用旋转变换将下述矩阵化为三对角阵1212221 1.11112111A解：首先，约化第一行第一列：这时k1。1、取i2,j=3,k=1.则22ikjkdaacos0.8944,ikad22215,sin0.4472.jkad 第39页/共69页2,3,1100000.89440.44720.00.44720.894400001P12,3,12,3,112.2361022.2361111.3416.0110.447221.34160.44721TAPAP2、取i2,j=4,k=1.则22ikjkdaacos0.7454,ikad222.236123,sin0.6667.jkad 第40页/共69页2,4,1100000.745400.6667.001000.666700.7454P22,4,11 2,4,1130032.33330.44720.1491.00.44722100.149110.3333TAPAP第41页/共69页取i3,j=4,k=2.则22ikjkdaacos0.9487,ikad 22( .4472)0.14910.4714,sin0.3163.jkad 2130032.33330.44720.1491.00.44722100.149110.3333A其次，约化第二行第二列：这时k2。第42页/共69页3,4,210000100.000.94680.3163000.31630.9468P33,4,223,4,2130032.33330.47140.00.47141.16631.5001.50.5TAPA P第43页/共69页则T2HIuu称为Householder矩阵，简称H矩阵.定义2、Householder变换设 为单位向量，即 ，nuRT1u u 第44页/共69页1. H矩阵是对称的正交矩阵.性质2. Householder变换也称为镜面反射变换.3. 对任意 ，存在Householder矩阵H，,nx yR使 的充要条件是 .yHx22yx第45页/共69页容易验证 H 矩阵是对称的且正交的矩阵，即2,TTHH H HHI由于 uuuuIHuT)2(且对如何与u 直交的向量v, 都有 0Tu v vvuuIHvT2因此，对任意 ，可设 ，则其H变换为nRxvuxvuHvHuHx第46页/共69页若将 中所用与向量u正交的方向视为一个镜面，有上述公式看到，H 变换不改变向量在镜面上的投影，并将向量沿法向量的投影改变为反方向等长度的向量，因此H变换也称为镜面反射变换nR由H变换的性质不难知道，对任意非零量 ，如果 , 则必存在H矩阵，使得 nRyx,22yxyHx 事实上，当取 时，即可验证由(3.4.2)式所定义的矩阵满足 的要求 2yxyxuyHx 第47页/共69页()121211/221,(,).(1),(,() ,0,0) ,TnTrrnii rHAaa aanraca aasign assarn 利利用用变变换换将将 化化为为三三对对角角矩矩阵阵的的过过程程 可可逐逐列列 行行 的的方方式式进进行行. .设设为为矩矩阵阵的的某某一一列列向向量量 如如想想将将此此向向量量的的后后个个分分量量化化为为零零 即即将将 变变成成其其中中1221122112,(0,0,() ,) ,21()rrnTrrrnrra aaacacHHacwacasign as aas ssign aaw 假假设设不不全全为为零零 由由于于且且故故存存在在矩矩阵阵 满满足足若若令令第48页/共69页(11200() ,) ,rrTn rTrrrnIHIwwwasign as aa 其其中中可取H矩阵为：第49页/共69页例3.5利用H变换将下述矩阵化为相似的三对角阵1212221 1.11112111A解：首先，约化第一行第一列：即取r1。2222123s (0,5,1,2)Tw 222115w第50页/共69页11000212033311420315152211031515H 111130077133151571181010157575110170157575AH AH 则：第51页/共69页其次，约化第二行第二列：即取r2。2271215153s 721(0,0,)(0,0,0.9381, 0.0667)15315TTw 2222.2612w第52页/共69页210000100000.98990.1416000.14160.9899H 2212130032.33330.47140.000100.47141.16611.500100.00011.50010.5005AH A H 则：第53页/共69页1122210,1,2,1nnniabbabAbbabin 考考虑虑对对称称的的三三对对角角矩矩阵阵这这里里约约定定第54页/共69页00112112( ),( )1.( )1( )( )()( )( )2,3,kkkkkkCIkppppapapbpkn 记记特特征征矩矩阵阵的的左左上上角角的的 阶阶主主子子式式为为并并约约定定利利用用行行列列式式的的展展开开公公式式,可,可得得到到递递推推公公式式: :特征多项式序列：(1) 序列中任意两个相邻的多项式在 内无公共根；(,) (2)设 ，则 与 反号。 0()0jpx 10()jpx 10()jpx 且 与 的根相互交错。 ( )0kp 1( )0kp 第55页/共69页定义：固定 ，称序列 中相邻两数符号相同的个数为同号数，记为 。01( ),( ),( )nppp( )S例：矩阵210121012的多项式序列为第56页/共69页012233( )1( )2( )(2)1( )(2)2(2)pppp 对于 同号数 的值间下表。 1,0,1,2,3,4, ( )S012233-+0-+0-0+-0+43210-10( )p1( )p2( )p3( )p( )S第57页/共69页)3.2( )( )0 ,nSp 定定理理同同号号数数等等于于特特征征方方程程在在区区间间上上的的根根的的个个数数矩阵特征值的计算方法-对分法(1) 由 ，可见 。( )pCC ipC 于是，在区间 上便可求出矩阵C的各个特征值。,ppCC 假设矩阵C的第i个特征值 ，由定理3.2， 。取区间 的中点 并计算 。00,ia b 00(), ()S ai S bi00,a b0001()2cab0()S c第58页/共69页(3) 若 则 ，记0(),S ci 00,ic b 1010,;ac bb否则， ，记 00,ia c 1010,.aa bc (4) 如此继续下去，经过n次对分过程，得到包含 的区间: 区间的长度依次是 。i0011,nna ba bab001(),0,1,2kbakn (5) 当n充分大时，区间的长度足够小，此时可取区间的中点 作为 的近似值。i第59页/共69页012233-+0-+0-0+-0+43210-10( )p1( )p2( )p3( )p( )S确定特征值的区间：第60页/共69页新的多项式序列：=11012112112112121( )( )( )( )()( )( )( )( )( ),( )( )0( )( ) 0( ) 0kkkkkkkkkkkkkkkkqppapapbpqppbaqqqaqq 0( )1.q 规定第61页/共69页计算同号数 的方法：( )S 01( ),( ),( )nqqq通过计算 中非负项的个数。第62页/共69页3.4.3 三对角矩阵特征向量的计算三对角矩阵特征向量的计算(1) 用反幂法求矩阵 的按模最小的特征值 和相应的特征向量y，则y是矩阵C的相应于 的特征向量，并可得到特征值 的更精确的近似值CI 00. (2) 由于 ,故矩阵C和矩阵A具有相同的特征值，矩阵A相应于的特征向量TCQ AQ .xQy 第63页/共69页若是用H变换将A化为C矩阵时，则有2221122nnnTCAHH H AH HHQ AQ 其中 为正交矩阵。 122nQH HH 若是用旋转变换将A化为C矩阵时，则有2,3,12,4,12, ,13,4,23,5,23, ,21, ,2nnnn nQPPPPPPP 第64页/共69页QR方法是一种通过逐次正交分解计算特征值的方法。 设A是一个n阶矩阵，按Schmidt正交化，将A1A正交分解为一个正交矩阵Q1和另一个上三角阵R1的乘积： A1=Q1R1交换因式矩阵Q1和R1的次序，得到新矩阵 A2=R1Q1因Q1为正交阵，故 Q1T =Q1-1，于是 A2= Q1T A1Q1这表明矩阵A2和 A1正交相似 第65页/共69页用A2代替 A1，重复上述步骤可得A3.一般地，若已知矩阵Ak ，且有分解 Ak=QkRk (3.5.1)其中Qk为正交阵，Rk为上三角阵，将Qk ,Rk换序之后，得到 Ak1= QkRk=QkT AkQk (3.5.2)通常，将公式(3.5.1) ，(3.5.2)称为矩阵A的QR算法 第66页/共69页.1221,n nkkkkkkkARQRAkAAAAQQ QQRRR RAAkQRAQR 定定理理33 对33 对任任意意由由算算法法产产生生的的矩矩阵阵序序列列具具有有以以下下基基本本特特性性：（1）1） 对对每每个个与与 正正交交相相似似，从从而而与与 具具有有相相同同的的特特征征值值。（2）2）若若记记则则 （ 的的 次次幂幂）有有分分解解式式第67页/共69页112(),0*,*n nijnnkknQRAaRAAAQRAAk 1 11212定定理理3.4(方3.4(方法法的的收收敛敛性性) 设) 设并并设设 的的特特征征值值 ，为为实实数数且且满满足足 则则有有的的算算法法产产生生的的矩矩阵阵序序列列本本质质收收敛敛于于一一个个上上三三角角矩矩阵阵 本 本质质收收敛敛 第68页/共69页