數值分析冪法和反冪法.doc_第1頁
數值分析冪法和反冪法.doc_第2頁
數值分析冪法和反冪法.doc_第3頁
數值分析冪法和反冪法.doc_第4頁
數值分析冪法和反冪法.doc_第5頁
免費預覽已結束,剩余18頁可下載查看

下載本文檔

版權說明:本文檔由用戶提供并上傳,收益歸屬內容提供方,若內容存在侵權,請進行舉報或認領

文檔簡介

1、題冪法和反冪法求矩陣特征值目具隨機產生一對稱矩陣,對不同的原點位移和初值(至少取 3 個)分別使用冪體法求計算矩陣的主特征值及主特征向量,用反冪法求計算矩陣的按模最小特征內值及特征向量,并比較不同的原點位移和初值說明收斂。容1.認真讀題,了解問題的數學原形;要2.選擇合適問題求解的數值計算方法;求3.設計程序并進行計算;4.對結果進行解釋說明;對于冪法和反冪法求解矩陣特征值和特征向量的問題將從問題分析,算法設計和流程圖,理論依據,程序及結果進行闡述該問題。一問題的分析:求 n 階方陣 A 的特征值和特征向量,是實際計算中常常碰到的問題,如:采機械、結構或電磁振動中的固有值問題等。對于 n 階矩

2、陣 A,若存在數和用n 維向量 x 滿足Ax=x(1)方則稱 為矩陣 A 的特征值, x 為相應的特征向量。法由高等代數知識可知,特征值是代數方程及|I-A|=n +a1n 1 + +a n 1 +a n =0(2)結的根。從表面上看,矩陣特征值與特征向量的求解問題似乎很簡單,只需果求解方程( 2)的根,就能得到特征值,再解齊次方程組說(I-A ) x=0(3)明的解,就可得到相應的特征向量。上述方法對于 n 很小時是可以的。但當n 稍大時,計算工作量將以驚人的速度增大, 并且由于計算帶有誤差, 方程(2)未必是精確的特征方程,自然就不必說求解方程(2)與( 3)的困難了。冪法是一種計算矩陣主

3、特征值(矩陣按模最大的特征值)及對應特征向量的迭代方法,特別是用于大型稀疏矩陣。反冪法是計算海森伯格陣或三角陣的對應一個給定近似特征值的特征向量的有效方法之一。二算法設計及流程圖1、冪法算法(1)取初始向量 u ( 0 ) (例如取 u( 0) =(1,1,1) T ), 置精度要求,置 k=1.(2)計算v ( k ) =Au( k1) ,m k =max(v( k ) ), u (k ) = v ( k) / m k(3)若 | m k = mk 1| ,則停止計算( mk作為絕對值最大特征值1 ,u ( k) 作為相應的特征向量)否則置k=k+1,轉( 2)2、反冪法算法(1)取初始向量

4、 u ( 0 ) (例如取 u( 0) =(1,1,1) T ), 置精度要求,置 k=1.(2)對 A 作 LU分解,即 A=LU(3)解線性方程組Ly(k ) =u( k 1) ,Uv ( k ) =y (k)(4)計算mk =max(v(k ) ), u ( k ) = v ( k) / m k(5)若 |m k =mk 1 |,則停止計算( 1/m k作為絕對值最小特征值n ,u ( k) 作為相應的特征向量) ; 否則置 k=k+1,轉( 3).冪法流程圖:開始輸入 A;m,u,index=pow(A,1e-6)k=0;m1=v=A*uvmax,i=max(abs(v)m=v(i);

5、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*um1=m;k=k+1vmax,i=max(abs(v)m=v(i);u=v/mabs(m-m1)|2 | | n |則計算最大特征值與特征向量的迭代格式為v (k ) =Au( k 1),m k =max(v( k) ), u (k ) = v ( k ) / m k( 1)其中 max(v( k) ) 表示向量 v ( k)絕對值的最大分量。2、對于冪法的定理按式( 1

6、)計算出 mk 和 u (k ) 滿足lim mk =1 ,lim u ( k) =x1kkmax(x1 )(二)反冪法算法的理論依據及推導反冪法是用來計算絕對值最小的特征值忽然相應的特征向量的方法。是對冪法的修改,可以給出更快的收斂性。1、反冪法的迭代格式與收斂性質設 A 是非奇異矩陣,則零不是特征值,并設特征值為| 1 | | 2 | | n 1 | n | 則按 A 1 的特征值絕對值的大小排序,有|1 |1| | 1 |nn 11對 A 1 實行冪法,就可得 A 1 的絕對值最大的特征值 1/ n 和相應的特征向量,即 A的絕對值最小的特征值和相應的特征向量。由于用 A1代替 A 作冪

7、法計算,因此該方法稱為反冪法,反冪法的迭代格式為v( k ) = A 1 u ( k1) ,m k =max(v( k) ), u ( k ) = v ( k) / m k( 2)2、對于反冪法的定理按式( 2)計算出的 m和 u ( k ) 滿足:kxnlim mk = 1 ,lim u (k ) =knkmax(xn )在式( 2)中,需要用到 A 1,這給計算帶來很大的不方便,因此,把(2)式的第一式改為求解線性方程組A v(k ) = u ( k 1)(3)但由于在反冪法中,每一步迭代都需求解線性方程組(3)式,迭代做了大量的重復計算,為了節省工作量,可事先把矩陣A 作 LU分解,即

8、A=LU所以線性方程組( 3)改為Ly( k) =u( k 1) ,Uv ( k ) =y (k)四、算法程序設計代碼冪法程序,在 matlab 中建立一個 M 文件并保存。%pow.mfunction m,u,index,k=pow(A,u,ep,it_max)if nargin4it_max=1000;endif nargin3ep=1e-5;endn=length(A);index=0;k=0;m1=0;m0=0;I=eye(n);T=A-m0*I;while k=it_maxv=T*u;vmax,i=max(abs(v);m=v(i);u=v/m;if abs(m-m1)ep;inde

9、x=1;break;endm=m+m0;m1=m;k=k+1;end在 matlab 輸入面板,輸入A=rand ( 4); %產生一個 4 維隨機矩陣B=A+A ;u=1 1 1 1 ;% 設立初始向量m,u,index,k=pow(B , u,ep,it_max)% 最多 可省略 2 個參數程序結束。在 M 文件中可以通過改變m0 的值改變原點位移,從而達到原點位移加速。反冪法程序設計代碼:在 matlab 中建立一個 M 文件并保存。%pow_inv.mfunctionm,u,index,k=pow_inv(A,u,ep,it_max)if nargin4it_max=1000;endi

10、f nargin3ep=1e-5;endn=length(A);index=0;k=0;m1=0;m0=0;I=eye(n);T=A-m0*I;invT=inv(T);while k=it_maxv=invT*u;vmax,i=max(abs(v);m=v(i);u=v/m;if abs(m-m1)B=rand(4);A=B+B A =0.26750.57760.63441.31300.57761.15030.76410.13670.63440.76410.02570.41931.31300.13670.41931.2248 u=1 1 1 1; m,u,index,k=pow(A,u) m

11、=2.6813u =0.85760.69340.56231.0000index =1k =49修改 M0=1e-3m =2.6814u =0.85760.69340.56231.0000index =0k =1001修改 M0=0%此時為冪法m =2.6815u =0.85760.69350.56231.0000index =1k =10修改 U=1 2 3 4修改 M0=1e-4m =2.6813u =0.85760.69340.56231.0000index =1k =9修改 M0=1e-3m =2.6805u =0.85760.69340.56221.0000index =1k =7修改

12、 M0=0m =2.6814u =0.85760.69340.56231.0000index =1k =9修改 U=3 5 6 7修改 M0=1e-4m =2.6819u =0.85770.69370.56241.0000index =1k =7修改 M0=1e-3m =2.6814u =0.85760.69340.56231.0000index =0k =1001修改 M0=0m =2.6820u =0.85770.69370.56241.0000index =1k =7總結以上 ,冪法如下:Um0muindexk0.00012.68130.8576 0.6934 0.5623 1.0000

13、14911110.0012.68140.5876 0.6934 0.5623 1.00000100102.68150.8576 0.6935 0.5623 1.00001100.00012.68130.8576 0.6934 0.5623 1.00001912340.0012.68050.8576 0.6934 0.5622 1.00001702.68140.8576 0.6934 0.5623 1.0000190.00012.68190.8577 0.6937 0.5624 1.00001735670.0012.69140.8576 0.6934 0.5623 1.00000100102.6

14、920.8577 0.6937 0.5624 1.000017反冪法結果顯示:在m0 為 0 時M0=0.001U=1 1 1 1M0=0.1u=1 1 1 1M0=0 u=1 3 5 7M0=0.1u=1 3 5 7M0=0.5u=1 3 5 7M0=0u=2 3 4 5M0=0.1u=2 3 4 5M0=0.7u=2 3 4 5綜上,反冪法結果如下:um0muindexk11110.10.3847-0.8996 1.0000 0.2726 -0.23641150.0010.3847-0.8996 1.0000 0.2726 -0.236411600.3847-0.8996 1.0000 0.2726 -0.236411613570.50.3847-0.8995 1.0000 0.2726 -0.23641270.10.3847-0.8996 1.0000 0.2726 -0.236411700.3847-0.8996 1.0000 0.2726 -0.2364120234 50.70.7091-0.6962 -0.4497 0.2196 1.0000150.10.3847-0.8995 1.000

溫馨提示

  • 1. 本站所有資源如無特殊說明,都需要本地電腦安裝OFFICE2007和PDF閱讀器。圖紙軟件為CAD,CAXA,PROE,UG,SolidWorks等.壓縮文件請下載最新的WinRAR軟件解壓。
  • 2. 本站的文檔不包含任何第三方提供的附件圖紙等,如果需要附件,請聯系上傳者。文件的所有權益歸上傳用戶所有。
  • 3. 本站RAR壓縮包中若帶圖紙,網頁內容里面會有圖紙預覽,若沒有圖紙預覽就沒有圖紙。
  • 4. 未經權益所有人同意不得將文件中的內容挪作商業或盈利用途。
  • 5. 人人文庫網僅提供信息存儲空間,僅對用戶上傳內容的表現方式做保護處理,對用戶上傳分享的文檔內容本身不做任何修改或編輯,并不能對任何下載內容負責。
  • 6. 下載文件中如有侵權或不適當內容,請與我們聯系,我們立即糾正。
  • 7. 本站不保證下載資源的準確性、安全性和完整性, 同時也不承擔用戶因使用這些下載資源對自己和他人造成任何形式的傷害或損失。

評論

0/150

提交評論