特征值问题的计算方法PPT课件

1、1212det()() ()()pnnnpIA 其中12;()pijnnnnij 称 为 的代数重数（简称重数）；ini ()iimnrankIA 为 的几何重数。i iimn 2Def设 ，n nAC iimn 对于矩阵 的特征值 ，如果Ai ，则称该特征值 为 的一个半单特征值。Ai 若 的所有特征值都是半单的，则称 是非亏损的。AA是非亏损的等价条件是 有n个线性无关的特征向量AA第1页/共57页3Def设 ，,n nA BC 若存在矩阵 ，使得P1BP AP 则称 和 是相似的。AB相似矩阵有相同的特征值1AxxPAP PxPx BPxPx Axx 设寻求已知矩阵 的相似矩阵 ，要求：

2、矩阵 的特征值和特征向量容易计算ABB本章QR算法的基本思想：第2页/共57页8 1 1. .Th1(, )iir 设 ，有 r个互不相同的特征值 ，n nAC 其重数分别为 ，则一定存在非奇异矩阵11()(),();iiinniikiJdiag JJCir 使得（Jordan分解）1(, )in ir n nPC 其中112( (), (), ()rP APdiag JJJJ 11()iijiiJ ()jiJ 且除了 的排列次序外， 是唯一的。J称作 的Jordan标准型AJ第3页/共57页8 1 2. .Th设 ，则存在酉矩阵 ，使得：n nAC （Schur分解）其中 是上三角矩阵，且适

3、当选择 ，可使 的元素HUAUT n nUC TUT按任意指定的顺序排列。8 1 3. .Th设 ，令()n nijAaC （圆盘定理）/*Disc Theorem*/则1( ):;,iiiijj iG AzCzaain 12( )( )( )( )nAG AGAGA 第4页/共57页8 1 4. .Th设 为对称矩阵，则存在正交矩阵n nAR （谱分解定理）/*Spectral Decomposition*/其中 是 的n个特征值。n nQR 使得1(,)TnQ AQdiag 1,nA8 1 5. .Th设 为对称矩阵，且 的特征值为n nAR （极大极小定理）其中 表示 中所有k维子空间的

4、全体。则有12n A0maxminniTiTuu Auu u 10min maxnn iTTuu Auu u nk nR第5页/共57页设 为对称矩阵，其特征值分别为8 1 6. .Th,n nA BR （Weyl定理）则有1212;nn 21 2;, ,iiABin 说明：对称矩阵的特征值总是良态的。注意：实际问题中矩阵一般都是由计算或实验得到，本身必然存在误差，不妨假设 BAA 21 2;, ,iiAin 第6页/共57页2 幂法与反幂法/*Power Method and Reversed Power Method*/幂法是计算一个矩阵的模最大的特征值和对应的特征 向量的一种迭代方法（又

5、称为乘幂法）。 一、幂法的基本思想与算法假设 是可对角化的，即 存在如下分解： n nAC A1AXX 其中1(,)ndiag 1;,n nnXxxC 不妨假设12n 对于0nuC 01122;nniuxxxC第7页/共57页011nnkkkjjjjjjjA uA xx 11121() )njkkjjjxx 011211() )knjkjjkjA uxx 11()x k 01kkkA uu 说明：当k充分大时, 的一个近似特征向量为1 特征向量可以相差一个倍数第8页/共57页01kkkA uu 因为向量 中含有未知量 ，实际不能计算1 ku但我们关心的仅是 的方向，故作如下处理：0kkkA u

6、u 令其中 为 的模最大分量k 0kA u11121011121() )() )njkkjjkjnkjkkjjjxxA uxx 11()xkx 第9页/共57页 幂法迭代算法：1kkkAuu 1limlimlimkkkkkkAuu 1Axx 1kFor k=1,2,3,1kkyAu kky if1kkuu 输出 和kuk kkkyu 001,uu 设 和 均收敛，由算法知kuk 幂法可以计算矩阵的模最大的特征值和对应的特征向量1ku 第10页/共57页解：Step12 1013 1014A 例1 1：利用幂法求下列矩阵 的模最大的特征值及相应的特征向量.A01 1 1()Tu (取初始向量为

7、)10355()TyAu15 11131 15()Tyu Step2212311555()TyAu25 222231112525()Tyu 第11页/共57页Step3321 84 24 92( .)TyAu34 92. 3330 36590 85371( .)Tyu Step4431 58543 92684 8537( .)TyAu44 8537. 4440 32660 80901( .)Tyu 特征值及相应的特征向量精确值为:4 7321. 0 26790 73201( .)Tu 第12页/共57页 幂法的收敛性：8 2 1. .Th12p 设 有 p个互不相同的特征值满足：n nAC 且

8、模最大特征值 是半单的，如果初始向量 在的特征子空间上的投影不为零，则由幂法算法产生的1 ku向量序列 收敛到 的一个特征向量 ，且数值1 0u1 1x序列 收敛到 。 k 1 特征子空间： 0Vx Axx 第13页/共57页证明：设 有如下Jordan分解：A11(,)pAXdiag JJX iinniJC 是属于 的Jordan块构成的块上三角矩阵i 111nJI 是半单的特征值1 10yX u 令将 和 如下分块：yX12(,)TTTTpyyyy 12pnnn12(,)pXXXX 12pnnn1010(,)kkkpA uXdiag JJX u 111222kkkPppX J yX J y

9、X J y第14页/共57页111222kkkpppX yX J yX J y 21112211()()pkkkppJJX yXyXy 0111222kkkkPppA uX J yX J yX J y021122111()()kpkkppkJA uJX yXyXy1111110()/()kiiiJJ 01110lim()kkkA uX y 第15页/共57页记11111X yxX y AXXJ 11111AXX JX 11111AX yX y 111Axx 1kkkAuu 011kkkA u 0110kkkkA uA u 1111()kX yukX y 1limkkux 是属于 的一个特征向量

10、1 1x1kkkAuu 1()kk 第16页/共57页几点说明：定理8.2.1条件不满足时，幂法产生的向量序列 ku可能有若干个收敛于不同向量的子序列；幂法的收敛速度取决于 的大小；21: 021122111()()kpkkppkJA uJX yXyXy加速方法：适当选取 ，对 应用幂法AI 称之为原点平移法1Axx 1Axxxx 1()()AI xx 原点平移法不改变矩阵 的特征向量A第17页/共57页幂法可以计算第二个模最大特征值 2 常用的方法：降阶方法（收缩技巧）设已经计算出模最大特征值 及其特征向量 1 1x111Axx 对向量 ，采用复的Household变换计算酉矩阵1xP111

11、Pxe 1111HPAPePAx 1111Px 1 1e 10HPAPB 其中 是n-1阶方阵B2 为 的模最大特征值 B第18页/共57页二、反幂法的基本思想与算法反幂法是求一个矩阵的模最小的特征值和对应的特征 向量的一种迭代方法（又称为反迭代法）。 设 ，则Axx 11A xx 1A 对 应用幂法就可以求得矩阵 的模最小的特征值和相应的特征向量。A不妨假设 的特征值为11nn A则 的特征值为11nn 1A 1ii 第19页/共57页 反幂法算法：For k=1,2,3,1kkAyz kky if1kkzz 输出 和kzk kkkyz 001,zz limknkzx nnnAxx 若 和

12、均收敛，由幂法知kzk 1kz limknk 收敛速度取决于 的大小1:nn 反幂法每次迭代都需要求解方程组1kkAyz 第20页/共57页 带位移的反幂法：实际应用中，反幂法主要用于求特征向量。且用某种方法已经得到 的特征值 的近似值A1 1 对矩阵 采用反幂法迭代格式为：AI 1 记假设 的特征值满足120n A1()kkAI vz kkkvzv For k=1,2,3,因为方程组的系数矩阵Doolittle分解化为两个三角方是固定的，通常采用( (选主元) )AI 程组求解，从而减少工作量。第21页/共57页AILU 求解方程组 化为：1()kkAI vz 1kkkkLyzUvy 带位移

13、的反幂法迭代格式： For k=1,2,3,1kkLyz kkUvy kkkvzv 收敛速度取决于 的大小12 当 时，收敛速度会非常快1 设矩阵 存在Doolittle分解：AI 第22页/共57页解：210131014A 例2 2：用带位移的反幂法求矩阵12679. ）的近似特征向量。33 对应特征值 （精确值为A1 2679. AILU 其中1136591027310 1.L 073211000366210000011.U 01 1 1()Tz 10Lyz 11 00000 63661 9999.y Step1第23页/共57页反幂法具有一次“迭代性”11Uvy 16776 393849

14、60 0001815 8199.v Step221Lyv 21 00001 10796 0073.y 22Uvy 220327 4114880 725446 73.v 2221 00000 73200 2679.vzv 所求近似特征向量为：第24页/共57页3 Jacobi方法Jacobi法：计算实对称矩阵全部特征值和相应特征向量 基本思想对,n nTARAA Q存在正交矩阵 ，满足12(,)TnQ AQdiag 记12(,)nQq qq 则1 2;, ,iiiAqq in 寻找正交相似变换 ，将矩阵 约化为对角阵即可QA正交相似变换求法：通过Givens变换来实现第25页/共57页 经典Ja

15、cobi方法设,n nijijjiAaRaa 令1122222111( )()()nnniiijFiijj iE AAaa 非对角“范数”当 时， 趋于一个对角阵0( )E A A( , , )( , , )TGp qAG p q先来研究一下矩阵的元素和矩阵 的元素之间的关系。AGivens变换记为 ，下面通过Givens变换( , , )G p q 对矩阵 进行约化，使得0( )E A A第26页/共57页例如取424;,npqcos ;sincs 记11a12a24a13a14a13a22a12a23a23a33a34a14a24a34a44a10s000c000100s 0c10s000

16、c000100s 0c11a12a 13a14a13a 23a33a34a10s000c000100s 0c 11a 13a 13a 33a ( )ipipiqbcasa ;iqipiqbsaca (, )ip q 第27页/共57页11a 13a 13a 33a 222( )pppppqqqbc ascas a 222qqpppqqqbs ascac a 22()()pqpqppqqbcs asc aa 选取适当的 ，由Givens变换将矩阵的下三角元素 尽可能多的化为零：即非对角“范数”尽可能的小。如果 ，则取0pqa 10;cs否则，令2tan;qqpppqaastca 2210tt 2

17、20()()pqpqppqqbcs asc aa 首先由 确定 第28页/共57页2210tt 21sgn( )t 14()t 2211111cossectantc ;stc 其次确定旋转平面( , )p q2222211( )nniiiiFFiiEBBbAb由F-范数的正交不变性设 经过一次正交相似变换后变为矩阵AB22222211()nniiiippqqppqqiibaaabb第29页/共57页222222222( )( )()( )ppqqppqqpqEBEAaabbEAa 222222ppqqppqqpqaabba 注意到旋转平面 的选取方法：( , )p q1maxpqijij na

18、a 经典Jacobi方法的迭代格式：101 2 3( );, ,kTkijkkkAaG AGAA k 8 3 1. .Th对于经典Jacobi方法产生的矩阵序列 ， kA存在 的特征值的一个排列 满足A12,n 12lim(,)knkAdiag 证明见教材第30页/共57页 经典Jacobi方法的迭代算法给定矩阵(),ijijjiAaaa选取最佳旋转平面 ：( , )p q1maxpqijij naa For k=1,2,3,21sgn( );t 2qqpppqaaa 计算211;tc ;stc ipipiqacasa ;iqipiqasaca ,ifip q 计算222pppppqqqac

19、ascas a 222;qqpppqqqas ascac a 220()()pqpqppqqacs asc aa ( )E A 直到需比较 个元素12()n n 第31页/共57页习惯上称 次Jacobi迭代为一次“扫描”12()n n 循环Jacobi方法每一次Jacobi迭代不是去选择最佳旋转平面，而是直接按照某种预先指定的顺序来“扫描”自然顺序：1 21 312 32 42( , )( , ),( , ),( , );( , ),( , ),( , );p qnn 3 43 531( , ),( , ),( , );(, );nnn 按照自然顺序的循环Jacobi方法是渐进平方收敛的第3

20、2页/共57页4 QR 方 法 基本思想利用正交相似变换将一个给定的矩阵逐步约化为上三角矩阵或拟上三角矩阵的一种迭代方法 QR方法的迭代格式设0n nAAC 令111ARQ 011AQ R 对矩阵 进行QR分解0A122AQ R 再对矩阵 进行QR分解1A一、QR基本迭代方法QR方法是目前计算矩阵全部特征值的最有效的方法之一；具有收敛快、算法稳定等特点。第33页/共57页一般地有：1mmmAQ R 1 2;, ,mmmAR Qm 1HmmmmAQ AQ 矩阵序列 中每一个矩阵都与原矩阵 相似 mAAQR方法的迭代算法：1mmmAQ R mmmAR Q For m=1,2,3,直到 近似为上三角

21、阵mA由迭代格式同时还得到：2112HHHmmmAQQ Q AQQQ HmmmAQ AQ 记12mmQQQQ 11mmmAQR 代入11mmmmQ QRAQ 第34页/共57页11mmmmQ QRAQ 等式两端同时右乘21mRR R112121mmmmmmQ QRRR RAQ RR R 记21kkRRR R 11mmmmQRAQ R 21111111mmmmmmQRA QRA Q RA mmmAQ R 即111 1()mmA er q 其中 是 的第一列， 是 的相应元素1()mq11rmQmR可以看作是对矩阵 用 为1()mqA初始向量的幂法所得到的向量。1eQR方法与幂法的关系：第35页/

22、共57页 QR方法的收敛性8 4 1. .Th设 的特征值满足 ，且 的LU则由QR迭代算法产生的矩阵 的对角线以120n Ai ()mmijAa ()miiaY第i行是 对应于 的左特征向量；若 有 分解，AY下的元素趋于0，同时对角元素 趋于1(, )iin （QR方法的收敛性质）证明：令1;XY 1(,)ndiag 则有AXY mmmAXYXLU ()mmmXLU ()mmX IEU 其中mmmIEL 0limmmE 第36页/共57页01(, )iirin且()mmmAX IEU ()mmQR IEU 1()mmmAQ IRE RRU 1mIRE R 当m充分大时， 存在唯一QR分解：

23、1mmmIRE RQ R 01()(, )miirin且mmQ RI 当m充分大时limlimmmmmQRI()()mmmmAQQR RU QR （QR分解）记 的QR分解为：XQR X第37页/共57页为保证上述QR分解中上三角矩阵的对角元为正，令111(,)nnDdiag 11211;(,)nnnnuuDdiaguu 11221()()mmmmmmAQQ D DD DR RU mmQ R 由QR分解唯一性知：2112()HHHmHHmmmmmmAQ AQDDQ Q AQQ D D 代入11HAXYXXQR R Q 12112()HHmHmmmmADDQ R R Q D D 趋于上三角阵第3

24、8页/共57页实际应用中遇到的多数特征值问题都是关于实矩阵的，所以自然希望设计只涉及实数运算的QR迭代。 实QR迭代格式设1n nAAR 二、实Schur标准形kkkAQ R 1kkkAR Q For k=1,2,3,为正交矩阵为上三角阵kQkR8 4 2. .Th（实Schur分解）设n nAR ，则存在正交矩阵n nQR ，满足：11121222mmTmmRRRRRQ AQR 其中 为实数或具有一对复共轭特征值的2阶方阵iiR第39页/共57页iiiiiiR 220()iiiiIR 特征值为iij，其中 为虚单位j见文献1311121222mmmmRRRRRR矩阵称为 的Schur标准形A

25、定理8.4.2说明：只要求得矩阵的Schur标准形，就很容易求得矩阵的全部特征值。缺陷： 很难得到Q第40页/共57页H 011h21h0022h32h12h0023h33h13h043h21,nh 31,nh 11,nh 41,nh 11,nnh 2,nh3,nh1,nh4,nh1,nnh 22,nh 32,nh 12,nh 0001,n nh ,n nh042,nh 12,nnh 称下述形状的矩阵为上Hessenberg矩阵三、上Hessenberg化第41页/共57页设 ，n nAR 寻求非奇异矩阵 ，使得 为上Hessenberg 矩阵，然后再对 进行 迭代。1QAQ Q1QAQ QR

26、 基本思想和约化过程：ijAa 记矩阵Q下面采用Householde变换寻找Step1 选取Householde变换 ，使得1H111Hpe 其中1111aAA 令11100HH 11111111101000aH AHHAH 第42页/共57页1111 1111aHHH AH 11a p2A0 2111AH A H 其中Step2 选取Householde变换 ，使得2H221Hqe 下面对作 同样的处理，以此类推2A其中223AA 令22100HH 2222322101000H A HAHH 22232HH A H 第43页/共57页令22100HH 1121122122101000aH H

27、 AH HHpeAH 111222apeH AH 3232AH A H 其中 11a p3A0 0 q0 22222232H A HHH A H 记2112H H AH H 0 11h21h0 0 22h32h12h 12H H为正交阵第44页/共57页按照前述方法，经过n-2步后，可以得到：221122nnHH H AH HHH 其中H 011h21h0022h32h12h0023h33h13h043h21,nh 31,nh 11,nh 41,nh 11,nnh 2,nh3,nh1,nh4,nh1,nnh 22,nh 32,nh 12,nh 0001,n nh ,n nh042,nh 12,

28、nnh 第45页/共57页122nPH HH 记 ，则TAPHP 称分解式 为矩阵 的上Hessenberg分解TAPHP A 上Hessenberg分解算法（8.4.1）1 ,( (: , )vhouse A kn k For k=1:n-211(: ,: )() (: ,: )TA kn k nIvvA kn k n 1111( : ,: )( : ,: )()TAn knAn kn Ivv 然后对上Hessenberg矩阵进行QR迭代设Hyy TP APyy APyPy 第46页/共57页上Hessenberg矩阵的QR迭代算法：1mmmHQ R mmmHR Q For m=1,2,3,

29、直到 的次对角元素趋于零mH注意：此时的QR分解可采用Givens变换来实现；用Givens变换得到的RQ仍然为上Hessenberg矩阵；上Hessenberg分解不唯一。若上Hessenberg矩阵的次对角元素均不为零，则称之为不可约的上Hessenberg矩阵。7Def定理8.4.3说明：在不可约的条件下，上Hessenberg分解在相差一个正负号的意义下唯一。见文献6第47页/共57页 实用的QR迭代算法（仅计算特征值）Step1 由算法8.4.1计算上Hessenberg矩阵：()TAHP AP Step2（QR分解） For k=1:n-11 ( ), ( )( , ),(, )c

30、 k s kGivens H k k H kk 1111( )( )( :, : )( :, : )( )( )c ks kH k knH k kns kc k Step3（计算RQ）Step4 重复步骤Step23，直到下次对角元素趋于零1111( )( )( : , :)( : , :)( )( )c ks kHn k kHn k ks kc k 第48页/共57页四、三对角化（对称矩阵的上Hessenberg化）TP APH 设 为对称矩阵， 的上Hessenberg分解为n nAR A其中 为对称三对角矩阵。HStep1 选取Householde变换 ，使得1H111Hpe 其中11111TaAA 令11100HH 111111111101000TaH AHHHA 11111 1111TaHHH AH 2111AH A H 其中 11ap2A0p0第49页/共57页111H A H主要工作量在于计算1THIvv 1111()()TTH A HIvvA Ivv 21111TTTTAAvvvv Avv Avv 21111()()TTTTAAv vv Avvv Avv 令112;()TuAvuv u v 1111()TTTTH A HAuvvuv v u v 1TTAvv 11122()()TTTTAuv uv vv uv uv 减少运算量第50页/