冪法和反冪法求矩陣特征值課程

20/20題目冪法和反冪法求矩陣特征值具體內容隨機產生一對稱矩陣,對不同的原點位移和初值(至少取3個）分別使用冪法求計算矩陣的主特征值及主特征向量，用反冪法求計算矩陣的按模最小特征值及特征向量,并比較不同的原點位移和初值說明收斂。要求1。認真讀題，了解問題的數學原形；２．選擇合適問題求解的數值計算方法；３.設計程序并進行計算；4．對結果進行解釋說明;采用方法及結果說明對于冪法和反冪法求解矩陣特征值和特征向量的問題將從問題分析，算法設計和流程圖，理論依據，程序及結果進行闡述該問題.一．問題的分析：求ｎ階方陣A的特征值和特征向量，是實際計算中常常碰到的問題,如：機械、結構或電磁振動中的固有值問題等。對于ｎ階矩陣A,若存在數和n維向量x滿足Ax=x(1）則稱為矩陣A的特征值,x為相應的特征向量。由高等代數知識可知，特征值是代數方程｜Ｉ-A｜=+a+…＋a+a=０（２）的根。從表面上看,矩陣特征值與特征向量的求解問題似乎很簡單，只需求解方程(2)的根,就能得到特征值，再解齊次方程組（I-Ａ）x=0（3）的解，就可得到相應的特征向量。上述方法對于ｎ很小時是可以的.但當n稍大時，計算工作量將以驚人的速度增大，并且由于計算帶有誤差,方程(2）未必是精確的特征方程，自然就不必說求解方程（2)與（3）的困難了.冪法是一種計算矩陣主特征值（矩陣按模最大的特征值）及對應特征向量的迭代方法，特別是用于大型稀疏矩陣。反冪法是計算海森伯格陣或三角陣的對應一個給定近似特征值的特征向量的有效方法之一。二。算法設計及流程圖1、冪法算法（1)取初始向量u（例如?。?(1,1,…1))，置精度要求，置k=1.（2）計算ｖ=Au,m=max（v)，ｕ=ｖ/ｍ(3)若｜m=m|〈，則停止計算(m作為絕對值最大特征值,ｕ作為相應的特征向量)否則置k=k+１,轉(２）2、反冪法算法(1）取初始向量ｕ(例如取u＝（1，1，…1）)，置精度要求，置k=1。（2）對A作LU分解，即A=LU(３）解線性方程組Ly=ｕ，Uv＝y(tǒng)(４）計算ｍ=maｘ（v)，u=v/m（5）若|ｍ=ｍ|〈,則停止計算（1/m作為絕對值最小特征值，u作為相應的特征向量）；否則置k＝k+１，轉(3).冪法流程圖：開始開始輸入A；[m,u,index]=pow(A,1e-6)k=0;m1=0v=A*u[vmax,i]=max(abs(v))m=v(i);u=v/mabs(m-m1)<1e-6index=1;break;輸出：m，u，index結束m1=m;k=k+1反冪法流程圖開始開始輸入A；[m,u,index]=pow_inv(A,1e-6)k=0;m1=0v=invA*u[vmax,i]=max(abs(v))m=v(i);u=v/mabs(m-m1)<1e-6index=1;break;輸出：m，u，index結束m1=m;k=k+1輸入A；[m,u,index]=pow(A,1e-6)三、算法的理論依據及其推導（一）冪法算法的理論依據及推導冪法是用來確定矩陣的主特征值的一種迭代方法，也即,絕對值最大的特征值。稍微修改該方法，也可以用來確定其他特征值。冪法的一個很有用的特性是它不僅可以生成特征值，而且可以生成相應的特征向量。實際上，冪法經常用來求通過其他方法確定的特征值的特征向量。1、冪法的迭代格式與收斂性質設n階矩陣A的特征值,，…,是按絕對值大小編號的,ｘ（i＝1，2，…,n）為對應的特征向量，且為單根，即|｜>｜｜≥…≥｜｜則計算最大特征值與特征向量的迭代格式為v=Ａu,m=mａx（v），ｕ=v/m（1）其中ｍax(ｖ）表示向量v絕對值的最大分量.2、對于冪法的定理按式（1）計算出ｍ和u滿足m=，u＝（二）反冪法算法的理論依據及推導反冪法是用來計算絕對值最小的特征值忽然相應的特征向量的方法。是對冪法的修改，可以給出更快的收斂性.1、反冪法的迭代格式與收斂性質設A是非奇異矩陣，則零不是特征值,并設特征值為｜|≥｜|≥…≥｜｜>|｜則按Ａ的特征值絕對值的大小排序，有｜|＞||≥…≥|｜對A實行冪法，就可得A的絕對值最大的特征值1／和相應的特征向量，即A的絕對值最小的特征值和相應的特征向量。由于用A代替A作冪法計算,因此該方法稱為反冪法,反冪法的迭代格式為v=Aｕ,ｍ＝maｘ(ｖ),u=ｖ/m(２）2、對于反冪法的定理按式（2）計算出的m和ｕ滿足:m＝，u=在式(2）中,需要用到A,這給計算帶來很大的不方便，因此，把（２)式的第一式改為求解線性方程組Av=u(3）但由于在反冪法中，每一步迭代都需求解線性方程組（3)式,迭代做了大量的重復計算,為了節(jié)省工作量，可事先把矩陣A作LＵ分解,即A＝LU所以線性方程組(3）改為Ｌy=ｕ，Uv=y四、算法程序設計代碼冪法程序,在ｍatlａb中建立一個M文件并保存。%pｏｗ。mfuｎctiｏn[m，u,indeｘ，k]=pow(A，u，ep，it_ｍax)ifｎａrgｉn<4iｔ＿mａx=1000；enｄiｆnargiｎ<3ep=１ｅ—5；endn=lengｔh（Ａ）；iｎｄｅx=0;k=0；m1=0;ｍ０=０;I=eyｅ（ｎ）;Ｔ＝A—m０＊I;whilｅｋ<=it_ｍaｘv=Ｔ*u；[vmax,i］＝ｍax(ａbs(v)）；m＝ｖ(i);u＝v/m；ifabs(m—m１）<ep;ｉndex=１；bｒeak；ｅndm＝ｍ+m0；m1＝m；ｋ=k+1；enｄ在mａt(yī)lａｂ輸入面板，輸入A=rａnｄ(4）；%產生一個4維隨機矩陣B=A+A’;u＝［111１］’;％設立初始向量［m，u，indeｘ，k］=pow（B，ｕ,eｐ，iｔ＿ｍax)%最多可省略2個參數程序結束。在M文件中可以通過改變m0的值改變原點位移，從而達到原點位移加速.反冪法程序設計代碼:在matlab中建立一個M文件并保存。%ｐｏw_inｖ.mｆunction［m，ｕ，ｉnｄex,k］=pow_ｉnv（A，u,eｐ,iｔ_mａx）ｉfnargｉn＜4it＿ｍax=1000；endifnarｇｉｎ<3ep＝1e-5；enｄn=ｌｅngtｈ（A)；indｅｘ＝0；ｋ＝0;m1＝０；m０=０;I＝ｅye(ｎ);T=A-m０＊I；iｎｖＴ=iｎｖ(T)；whilek<=it_ｍaxｖ=invT＊u；［ｖｍaｘ,i］=maｘ(aｂs（v））;m=v(i）；ｕ=v/ｍ;ifaｂs（ｍ－m1)＜epindex=1;ｂｒeａk;ｅndm１=m；k=ｋ+1；eｎdm=1/ｍ;ｍ=m＋ｍ0；在mａｔlab輸入面板,輸入Ａ=rａnd（4）；％產生一個4維隨機矩陣B=A+A’；u=［１１11］’；％設立初始向量［m，u,ｉndｅx,k]=pｏw＿inv（B，u，ep,ｉt＿maｘ）％最多可省略2個參數程序結束.在Ｍ文件中可以通過改變m０的值改變原點位移，從而達到原點位移加速?！窘Y果顯示】%在M０=1e—4>〉B=ranｄ（4）；＞〉A=B＋B’A=0。２６750.57760.６34４1.３１３00.５7761。15030.7６410.1３670．６3440．764１0．02５70.４１931．３1300。13６7０.41９31.２248＞>ｕ=［1111］＇；>〉［m,u，inｄｅx,ｋ］=pow(A,u)m＝２.6813ｕ=０。857６0.6９340.５6231.0０00inｄex=1k＝４９修改M０＝１e-3m=2．６8１４ｕ=0．8576０.6９３40．5623１.０000indｅｘ＝０k=１0０1修改Ｍ０=０％此時為冪法m=2.６815ｕ＝0．8５７60．69350。562３1．0000ｉｎdｅx=１k=10修改U=[1２34］修改M0=1e-4m＝２。6８13u=0。85760。69340.５6231．0０0０ｉnｄex＝１ｋ=9修改Ｍ０=１e-3m＝2．6805u=0。85７60。69３40.56221.00００ｉndｅx=１k=7修改M０=0ｍ=2．６８14u＝0．8５７６0.6９３40。5623１。０000indｅx=1k=9修改U＝［35６7］修改M0=1e-4m＝2.6819u=0．857７0。６9370.56241．０000ｉnｄex＝1k＝7修改M0＝１e—3m=2.６814ｕ=０。85760．69３40.562３１．000０indｅx=0k＝１0０1修改M0=0ｍ=２.6８２０u=０。8５770.６937０.56241.0000ｉｎdex=1k=7總結以上，冪法如下：Um０ｍuindeｘk[１111]0．00012.681３［0.８576０。6９３4０．56２31。0000］14９０.0012。6814[0.５8760。6９3４0。５6２3１．0000］０1０0１０２。68１5［0.８5760。6９３50。56231．０0０0]１10［1234］0．00012。６81３［0。857６0.６9３40．5623１.0０00]１90．0012．680５［0。85760。６９340。56221。０0０0］1702．6814[０。８５760。69340.56２31。0000]１9［3567］0。０0０１2.68１9［0。85770．6937０.５６２41.0000]1７０。0012．6９１４[０。85760．69340.５６231．０000］01０0102．692[0.8５７70.６９3７０．５6２41。０0０0]17反冪法結果顯示：在ｍ０為0時Ｍ0=0.0０１Ｕ=[11１1］Ｍ0=0。１u＝［1111］M0=0u=［1357]M0=0．1u＝[１35７］Ｍ0=0.5u＝［1357］Ｍ０＝0u=［２345]M０＝0。１u=[2345］M０=0．7ｕ=［2345]綜上，反冪法結果如下：um0mｕiｎdｅxk［１111］0．10.3847[-0．8９９６1．00000.27２6－0。２364］1150．0010．３847[-0。8９9６1．00000.２7２６—0.2364］1１60０.3８47［-０。８９9６1.000０0。272６-０.2364]１１6［13５7］0．50.3847［—0．899５１.00００0.2726—0.2364］1270。1０。38４７［—0.89９61。０0０00.2７26—0。23６4］117０0。3847［—0.８９961.00000。2７2６－０。2364］120［23４5］0。70．７0９1[－０。6962-0.４４９7０。2