實驗解線性方程組的基本迭代法實驗

1、數(shù)值分析實驗報告姓名忘川學號1205025106系別數(shù)學系班級12 級主講教師指導教師實驗日期2014/6/25專業(yè)信息與計算科學專 業(yè)課程名稱數(shù)值分析同組實驗者無一、實驗名稱：實驗八、解線性方程組的基本迭代法實驗二、實驗目的：1深入理解 Jacobi 迭代法、 Gauss-Seidel 迭代法和 SOR 迭代法；2通過對 Jacobi 迭代法、 Gauss-Seidel迭代法和 SOR 迭代法的程序設計，提高程序設計的能力； 3應用編寫的程序解決具體問題，掌握三種基本迭代法的使用，通過結(jié)果的分析了解每一種迭代法的 特點。三、實驗內(nèi)容及要求：1根據(jù) Matlab 語言特點，描述 Jacobi

2、迭代法、 Gauss-Seidel迭代法和 SOR 迭代法。2編寫 Jacobi迭代法、 Gauss-Seidel迭代法和 SOR迭代法的 M 文件。20 203給定 A R20 20 為五對角矩陣11324111322411113422411113422411134221 1 342(1) 選取不同的初始向量 x(0) 及右端面項向量 b，給定迭代誤差要求，分別用編寫的Jacobi 迭代法和Gauss-Seidel 迭代法程序求解，觀察得到的序列是否收斂？若收斂，通過迭代次數(shù)分析計算結(jié)果并得出你的結(jié)論。(2)用編寫的 SOR 迭代法程序，對于 (1)所選取的初始向量 x(0) 及右端面項向量

3、 b 進行求解，松馳系數(shù) 取 12的不同值，在 x(k) x(k 1) 10 5 時停止迭代，通過迭代次數(shù)分析計算結(jié)果并得出你的結(jié)論。四、實驗步驟(或記錄)一、算法描述5 / 22a11a12 La1na11Aa21a22 La2na22MMMMOan1an2 Lanna11得到Da22Oann由0a210MMOan 1,1an 1,2L0annan1an2Lan,n0a210LMMOan 1,1an 1,2Lan1an2L0an,n 100a12Ka1,n 1a1n0Ka2,n 1a2nOMM0an 1,n1000 a12Ka1,n 10Ka2,n 1UOM0則有：第一步： Jacobi 迭

4、代法a1na2nM , 則有： A D L Uan 1,n0Ax bAAx Db L U (D L U)x b Dx (L U)x b x D (L U)x D b令 J D (L U) 則稱 J 為雅克比迭代矩陣 f D b由此可得雅克比迭代的迭代格式如下：x(0) , 初始向量x(k 1) Jx(k) f ,k 0,1,2,L第二步Gauss-Seidel 迭代法Ax b(D L U )x b (D L)x Ux b x (D L) Ux (D L) b A D L U令 G (D L) U ，則稱G為Gauss-Seidel 迭代矩陣 f (D L) b由此可得 Gauss-Seidel

5、 迭代的迭代格式如下：x(0) ,初始向量x(k 1) Gx(k) f ,k0,1,2,L第三步SOR 迭代法w0ADLU 1 (D1wL (1 w)D wU )(D1 wL)(1 w)D wU )www令Mw1 (DwL),N1(1 w)D wU )則有：AMNwwAx bAMLWNM(MN)x b Mx Nx b x MNx MbNM,令W fMb帶入N 的值可有LW (1 w)D wU) (D wL) 1(1 w)D wU) (D wL)f 1 bw 1(D wL) 1b1 (D wL)w稱 LW 為 SOR 迭代矩陣，由此可得 SOR 迭代的迭代格式如下：x(0) ,初始向量(k 1)

6、f,k 0,1,2,L二、算法程序Jacobi 迭代法的 M 文件： function y,n=Jacobi(A,b,x0,eps)%*%函數(shù)名稱Jacobi雅克比迭代函數(shù)%參數(shù)解釋A系數(shù)矩陣%b常數(shù)項%x0估計解向量%eps誤差范圍%返回值%y解向量%n迭代次數(shù)%函數(shù)功能實現(xiàn)線性方程組的 Jacobi 迭代求解%*n=length(A);if nargin=1disp( 譜半徑大于等于 1，迭代不收斂，無法進行 ); return;endn=1;for i=1:1:nif sum(A(i,i)=n)=nerror( 輸入的 A 矩陣的對角線元素不能為 0); return;endendy=B

7、*x0+f;while norm(y-x0)=eps&n100x0=y;y=B*x0+f;n=n+1;endGauss-Seidel 迭代法的 M 文件和 function y,n=GaussSeidel(A,b,x0,eps) %*%函數(shù)名稱GaussSeidel高斯賽德爾迭代函數(shù)%參數(shù)解釋A系數(shù)矩陣%b常數(shù)項%x0估計解向量%eps誤差范圍%返回值%y解向量%n迭代次數(shù)%函數(shù)功能實現(xiàn)線性方程組的 Jacobi 迭代求解%*n=length(A);if nargin=1disp( 譜半徑大于等于 1，迭代不收斂，無法進行 );return;n=1;for i=1:1:nif sum(A(i,

8、i)=n)=n0);error( 輸入的 A 矩陣的對角線元素不能為 return;endendy=B*x0+f;while norm(y-x0)=eps&n100x0=y;y=B*x0+f;n=n+1;endSOR 迭代法的 M 文件function y,n=SOR(A,b,x0,w,eps)%*%函數(shù)名稱SOR松弛迭代函數(shù)%參數(shù)解釋A系數(shù)矩陣%w松弛因子%b常數(shù)項%x0估計解向量%eps誤差范圍%返回值%y解向量迭代次數(shù)%函數(shù)功能 實現(xiàn)線性方程組的 Jacobi 迭代求解%*n=length(A);if nargin=eps&n b=ones(20,1); x0=zeros(20,1);

9、eps=0.005; y,n=Jacobi(A,b,x0,eps) y =0.48130.57290.63210.65130.66000.66320.66460.66510.66530.66530.66530.66530.66510.66460.66320.66000.65130.63210.57290.4813n =9在 Command Window 中輸入： Gauss-Seidel 迭代法求解： 在 Command Window 中輸入： A=zeros(20,20);for i=1:1:20A(i,i)=3;endfor i=1:1:20for j=1:1:20if abs(i-j)=

10、1A(i,j)=-1/2;endendendfor i=1:1:20for j=1:1:20if abs(i-j)=2A(i,j)=-1/4;endendend13 / 22 b=ones(20,1); x0=zeros(20,1); eps=0.005; y,n=GaussSeidel(A,b,x0,eps)y =0.48140.57320.63250.65180.66060.66400.66540.66600.66620.66630.66630.66630.66610.66560.66420.66090.65210.632813 / 220.57340.4816 n =7第二次給定初始向量

11、 為 20 行一列的 0 右端面項向量 b=20 行一列的 1.001 迭代誤差要求 0.005Jacobi 迭代法求解：在 Command Window 中輸入A=zeros(20,20);for i=1:1:20A(i,i)=3;endfor i=1:1:20for j=1:1:20if abs(i-j)=1A(i,j)=-1/2;endendend b=1.001*ones(20,1); x0=zeros(20,1); eps=0.005; y,n=Jacobi(A,b,x0,eps)y =0.41460.48560.49780.49990.50020.50030.50030.50030

12、.50030.50030.50030.50030.50030.50030.50030.50020.49990.49780.48560.41467Gauss-Seidel 迭代法求解： 在 Command Window 中輸入 A=zeros(20,20);for i=1:1:20A(i,i)=3;endfor i=1:1:20for j=1:1:20if abs(i-j)=1A(i,j)=-1/2;endendend b=1.001*ones(20,1); x0=zeros(20,1); eps=0.005; y,n=GaussSeidel(A,b,x0,eps)y =0.41450.4856

13、0.49780.49990.50030.50030.50030.50030.50030.50030.50030.50030.50030.50030.50030.50030.50000.49800.48580.4146n =5在較少的迭代次數(shù)下即可的到滿足誤差要求的解。觀察計算結(jié)果得到的序列可以看出其是收斂，2）第一次給定初始向量 為 20 行一列的 0，右端面項向量 b=20 行一列的 1 迭代誤差要求 0.00005 松弛因子為 1.5 在 Command Window 中輸入A=zeros(20,20);for i=1:1:20A(i,i)=3;endfor i=1:1:20for j=1

14、:1:20if abs(i-j)=1A(i,j)=-1/2;endendendfor i=1:1:20for j=1:1:20if abs(i-j)=2A(i,j)=-1/4;endendend b=ones(20,1);19 / 22 x0=zeros(20,1); w=1.5; eps=1e-5; y,n=SOR(A,b,x0,w,eps)y =1.0e+012 *-0.5082-0.9690-1.5400-2.1738-2.8767-3.6356-4.4375-5.2635-6.0901-6.8885-7.6243-8.2578-8.7437-9.0319-9.0675-8.7940-7

15、.0831-5.4598-3.5651n =100列的 0，第二次給定初始向量 為 20 行一 右端面項向量 b=20 行一列的 1 迭代誤差要求 0.00005 松弛因子為 1.2在 Command Window 中輸入A=zeros(20,20);for i=1:1:20A(i,i)=3;endfor i=1:1:20for j=1:1:20if abs(i-j)=1A(i,j)=-1/2;endendend b=ones(20,1);23 / 22 x0=zeros(20,1); w=1.2; eps=1e-5; y,n=SOR(A,b,x0,w,eps)y =0.27920.32460.33190.33310.33330.33330.33330.33330.333