7特征值問題的迭代解法1

1、7特征值問題的迭代解法17.1投影算法17.2 Raylcigh-Ritz 算法 27.2.1 對(duì)稱矩陣27.3 Lanczos 算法 47.4 Arnold!算法 67.5 非對(duì)稱Lanczos算法 7參考文獻(xiàn)第七章特征值問題的迭代解法但矩陣規(guī)模很大時(shí)，計(jì)算其所有的特征值和特征向量是非常困難的而在實(shí)際應(yīng)用中, 我們通常只對(duì)其中的某些特征值和特征向量感興趣，因此也沒必要計(jì)算所有的特征值和特 征向量.本章討論計(jì)算部分特征值和特征值向昴的迭代解法.這些算法的存儲(chǔ)鼠要遠(yuǎn)小于 O(“2),運(yùn)算量也遠(yuǎn)小于O(n3).7.1投影算法最簡(jiǎn)單的特征值問題就是僅僅計(jì)算一個(gè)特征值，如計(jì)算模最大的特征值.這時(shí)我們可

2、 以使用幕迭代算法.扇迭代:計(jì)算最大特征值1 Given x<°>2. for until converge3 y(0 = Ax(,_1)4. x(0 = y(/)/|y(0|25. end幕迭代所產(chǎn)生的迭代向fit *0), x，N")生成一個(gè)Krylov子空間驗(yàn)(A,x(0) = span快°),Ax(°),Am-1x(0).在幕迭代中,我們?nèi)榻铺卣飨蛄?顯然，如果我們?cè)谥姓页?最住”的 近似特征向量，則收斂速度就可能會(huì)大大加快.下面我們討論如何在心=心(A,/®)中尋找最佳”的近似特征向量.設(shè)A G時(shí)匕 并設(shè)兒和厶是R”的

3、兩個(gè)也維子空間.投影算法就是在尋找A的近似特征對(duì)(入丘),満足 下面的Petrov-Galerkin條件find A eC and x G "Km such that Ax - Ax 丄厶.(711)這樣的算法我們稱為斜投影算法.如果我們?nèi)≯潭?，則上面的算法就是一個(gè)正交投影 算法，此時(shí)條件(7.1.1)稱為Galcrkin條件.設(shè)V0,和Wo，VVb.,IVm_i分別是冷和厶的一組標(biāo)準(zhǔn)正交基.令卩沏二vo,Vjvw_!, Wm =凹0，呦，閃加_.則對(duì)任意X e %、存在向量y G R"使得x = Vmy,7.2 Raylcigh-Ritz 算法2即丘可以由巾*線性表出根

4、據(jù)條件(7.1.1),我們可得AVmy-AVmy. w) = 0, i = 1,2,.tm(7.1.2)(7 丄 3)WAVmy = AWVmy.這是一個(gè)廣義特征值問題.如果我們?nèi)≯潭牟⒘頦m = Vmt則(7.1.2)就化為Tmy =血,其中Tm = V,AVm e RM".這意味著Qy)是矩陣Tm的一個(gè)特征對(duì).由于m通常比較小，因 此我們可以使用先前討論的方法(如QR迭代)來計(jì)算(“).這樣我們就可以計(jì)算岀A的一 個(gè)近似特征對(duì)(入対,其中"Vwy.7.2 Rayleigh-Ritz 算法事實(shí)上，我們可以在心中找岀m個(gè)最佳近似特征向量及相應(yīng)的最佳近似特征 值.這些近似特

5、征值和近似特征向量就是Ritz值和Ritz向量.定義71設(shè)驗(yàn)是Rm的一個(gè)加維子空間，它的一組標(biāo)準(zhǔn)正交基為v0,vi,.,vw_b并 令5 = vo.vj，.記Tm = VAVm,設(shè)(入刃是Tm的一組特征對(duì)，即Tmy二血且 Hylh = 1.則我們成A是A的一個(gè)Ritz值丘二Vmy是A的一個(gè)Ritz向量.Ravleigh-Ritz算法就是用Ritz值和Ritz向駅來近似A的特征值與特征向昴;.算法 7.1 (Rayleigh Ritz procedure)1. Compute an orthonormal basis of %: Vm = vq,Vj,.,2. Compute the eigen

6、values of the matrix Tm = VAVm9 i.e. 9 the Ritz values of A3. Select k desired ones: A.Xz.Ak where k < m.4. Compute the eigenvectors yi of Tm associated with 4/, / = L 29. 9k.5. Approximate the eigenvectors of A by the Ritz vectors x, = Vmyi, i=l,2,7.2.1對(duì)稱矩陣這里我們討論4對(duì)稱時(shí),Rayleigh Ritz算法的性質(zhì).設(shè) w Rg&qu

7、ot;是一個(gè)列正交矩陣,且使得V =必,Vu g 是一個(gè)正交矩陣.于是我們有VTAV = VmtVurAVm,Vu =?；疺AVUVAVm V；AVU7.2 Raylcigh-Ritz 算法3由于A對(duì)稱的，故Tm = VAVm e Rmxm也是對(duì)稱的此時(shí),Ritz值和Ritz向量具有下面的最 優(yōu)性質(zhì).定理7.1設(shè)A e Rnxn對(duì)稱，則對(duì)任意的對(duì)稱矩陣R e 有IHVm 一 曲lb > I叫一匕幾lb,即 AVm - VtnR2 在 R = Tm 處取最小值，此時(shí) AVm- VmRh = VAVm2-證明 Let R = Tm +Z where Z e RmXm is symmetric

8、. Note that both A and Tm arc symmetric and Tm = VAVm. we haveI叫-VW/?H1 = (AVm - VmR)T(AVm - VmR)=(AVm - Vm(Tm + Z)T(AVm - Vm(Tm + Z) =AVm - VmTm - VAVmZ + TmVlVmZ-ZVA Vm + ZV $Tm + ZT 叮 VmZ=AVm 一 VmTm + |Z|.It follows that AVm -is minimized when Z = 0 andAVm - VmTm2 = VVTAVm - VmTw|= |Va(VVm)|= |vJ

9、avw|2.注：定理7.1中的2范數(shù)可以改成任意的酉不變范數(shù)，如F范數(shù).定理7.2設(shè)A w Rnxn對(duì)稱，并設(shè)Tm = UW是臨=VAVW的特征值分解.設(shè)0 w Rnxm 是滿足span(0 = %的任意列正交矩陣,Dwirx/n是任意對(duì)角矩陣.我們有AQ-QD2>AVm-VmTm2,且當(dāng)Q = VmU, D二A時(shí)等式成立.證明 Because span(Q) = 7C = span(Vm)t it follows that there exist a matrix W 6 RM,XW such thatQ = 5 W. As Q is column orthogonal, we hav

10、e/ = er(2 = (VmW)TVmW = WTVVmW = WTW,7.3 Lanczos 算法4which implies that W is orthogonal. Let WDWT = +Z Then we haveIIAe= AVmW-VmWD=I叫-VmWDWT= l|Avm-rm|g + |Zii> I 叫-vwrm|?.If we choose W = U and D = A, then Z = WDWT - Tw = UNljT - Tm = 0 and, hence, the above equality holds.This result shows that

11、the minimum of AQ - QD2 over all byk column orthogonal matrices Q with span(Q) = % and over all diagonal matrices D is attained by Q 二 VmU and D = 2定理7.1和定理7.2表明，在|42 - QD2極小的意義下,Ritz值是特征值的最佳”近 似.所以我們用Ritz值作為特征值的近似是有道理的.7.3 Lanczos 算法設(shè)A 6 Rnxn是對(duì)稱矩陣.Lanczos就是利用Lanczos算法來計(jì)算 心 的基和Tm = V,AVmt然后計(jì)算A的Ritz值

12、和Ritz向量.算法 7.2 (Lanczos Algorithm)1. Choose a vector vq such that |voll = 1, and set 0o = 02. For j = 0,1,., Do3.Compute w = Avj -4.aj+i = (w. vj)5.w = w - aj+ivj6.0j+l = 11創(chuàng)127.If 0j+i = 0, then stop8.Vj+i = w/0丿+i9.Compute the eigenvalues and eigenvectors of T;10.Check the convergence11 EndDo在Lanc

13、zos算法7.2中,迭代m步后，向量vo, vi,., vm_i構(gòu)成子空間vo) = span vo, Avo,.,A111"1 vo),的一組基,并且有7.3 Lanczos 算法-5 -其中知二0,0，0,1卩釵”,設(shè)(“)是&的一個(gè)特征對(duì)，則有= VmTmy +0力(刃如.于是Ax-x=m(ey)vm,即|Ax - Xxi -隔2：刃|,其中x = Vmy.如果弘(圧刃|很小,則我們就認(rèn)為是A的某個(gè)特征值的很好的近似.事實(shí) 上,關(guān)于Ritz值入我們有下面的性質(zhì).引理7.1設(shè)A w R”x“是對(duì)稱矩陣，設(shè)r = Ax - Ax,其中x H 0.則min |八祗陣衣 o(A

14、)|x|2其中tr(A)表示A的譜，即所有特征值組成的集合.證明 Let A = UNljT be the eigendecomposition of A. Then we haver = (A Al)x = t/(A Al)UT x.Denote A = diag(4i9 Az ,4/i). For any vector z =乙 1,乙2八,Zn Rnt it holds thatni=l> Y min |久-碉J ds) 1=1=Iklll min A 一 1|2.A&r(A)Therefore, we have怖=I叭12 = IKA -測(cè) 21叫勰)IWb 勰片九whi

15、ch complete the proof.7.4 Arnoldi 算法-6-由引理71可知，存在A的某個(gè)特征值九使得His1爲(wèi)(龍y)“mll21%州2欝=皿.(7.3.1)7.4 Arnoldi 算法-#-注：在前面的討論中，我們都沒有考慮實(shí)際計(jì)算時(shí)可能的誤差在實(shí)際計(jì)算中，由于浮點(diǎn)運(yùn) 算的誤差，即使m很小（如m= 10或加二20）,也可能會(huì)導(dǎo)致向量v,失去正交性.這時(shí)我 們必須釆取一些補(bǔ)救措施，最簡(jiǎn)單的方法就是對(duì)它們重新來一次正交話，即在算法7.2的第 5步后加上一條語句Jw = w 工w,1=1這個(gè)過程就稱為帶全正交過程的Lanczos算法.顯然，這個(gè)過程是非常費(fèi)時(shí)的.另外一個(gè)可行的方法

16、就是選擇性正交7.4 Arnoldi 算法這里考慮非對(duì)稱情形，即計(jì)算非對(duì)稱矩陣A的特征值與Lanczos算法的思想相類似, 我們可以使用Arnoldi算法,即通過Arnoldi算法計(jì)算他的標(biāo)準(zhǔn)正交基巾,力，,和 上 Hessenberg 矩陣 Hm = VAVm,使得Ab = VmHm +and VAVm =弘.但此時(shí)Hm只是上Hessenberg,而不是對(duì)稱三對(duì)角但我們同樣可以通過計(jì)算Hm的特征值 和特征向量來得到A的Ritz值的Ritz向量，并用它們來近似A的特征值和特征向量.設(shè)（入刃是Hm的一個(gè)特征對(duì)，其中|y|2 = 1,則A(V%y) = mH my += AVmy +所以IlAx

17、一 Ax2 = II亦+1"2%MII2 = I亦+1"I kiLyh其中x = Vmy.若心+sl ey足夠小我就認(rèn)為（A.x）是A的某個(gè)特征對(duì)的近似.算法 7.3 (Arnoldi Algorithm)1 Choose a vector vq such that Hvolh = 12. For j = 0919. .9 Do3. Compute = Avy4. For / = 0,Do5 hij =7.4 Arnoldi 算法-#-7. EndDo8.5非對(duì)稱Lanczos算法78勿+1J = I|W>1|29. If h田j = 0 then stop勺+1 =

18、 Wj+Jhj+w11. Compute the eigenvalues and eigenvectors of Tj12. Check the convergence13. EndDo由于a是非對(duì)稱的其特征值可能是復(fù)的，或者是壞條件的，此時(shí)Lanczos算法的一些 最優(yōu)性質(zhì)就不再成立.盡快如此，目前還是存在一些有效Amoldi算法的實(shí)現(xiàn)方式,可參見 5, & 7.7.5 非對(duì)稱Lanczos算法非對(duì)稱Lanczos算法就是Lanczos算法在非對(duì)稱矩陣上的推廣,它是基于Lanczos雙正 交化過程.設(shè)vo和»vo是任意的非零向量，并設(shè)vo) = spanvo, Avo,.,

19、vo)和心wo) = span) wo, wo. 9 (A7)1 wq.Lanczos雙正交化過程就是計(jì)算v0)和*朋(A, w()的基vf|和w“,滿足(v,)和w“相 互正交，即所以其中 ym = vo,vb.,vm_i, Wm = wo,VV1，.,»VW_1).注意，通常(v,)«0或呵星。本身并不一定正交，故匕，和通常并不列正交.根據(jù)Lanczos雙正交化過理，我們可得ATWm = WmT +0”MM仁7.5非對(duì)稱Lanczos算法 10 其中em = 0所以= WlVmTm + ym(W；vm) = Tm.設(shè)；S是險(xiǎn)的特征值，其對(duì)應(yīng)的右特征向就和左特征向帚分別為

20、y和Z,且llylb = |剛2二1,即于是有my = ymTmy + ymmy = Vmy + ymeyvmi(WmZ)TA = (ATWmZ)T = (WmTz)T +0泌 吒=WmZ)T +0泌腳若(必刃I和(金z)l足夠小，我們就認(rèn)為；I是A的某個(gè)特征值的近似，而vmy和WM就 是相應(yīng)的右特征向量和左特征向量的近似.算法 7.4 (Nonsymmetric Lanczos Algorithm)1 Choose two vectors v© and wq such that切,wq) = 12 Set 0o = 0 and ya = 03. For j = 0,1,., Do4

21、.Compute a+i = (Avj.Wj)5.可+1 =人與 _ ajvj-0pj76.也汁 1 = Atwj - aWj - yjWj_i7.yj+i = l®+i,叼+i卩佗8.If yy+i = 0f then stop9.0 田=(Vj+i9Wj+l)/yj+l10.vj+i =巧+1/如11.旳+1 =町+1/0田12.Compute the eigenvalues and eigenvectors of Tj13.Check the convergence14.EndDo非對(duì)稱Lanczos算法的顯著優(yōu)點(diǎn)就是節(jié)省運(yùn)算駅缺點(diǎn)是更容易被中斷.參考文獻(xiàn)1 J O Aascn

22、, On the reduction of a symmetric matrix to tridiagonal form, BIT, 11(1971), 233-242.2 Nicholas J. Higham, Accuracy and Stability of Numerical Algorithms. Second Edition, SIAM, Philadelphia, 2002.3 R.A. Hom and C.R. Johnson, Matrix Analysis. Cambridge University Press, New York, 1985.4 R.A. Hom and C.R. Johnson, Topics in Matrix Analysis. Cambridge University Press, New York, 1991.5 R. Lehoucq, Analysis and Implementation of an Implicitly Restarted Arnoldi Iteration. Ph.D. thesis, Rice University,