当前位置:文档之家› 实验解线性方程组的基本迭代法实验

实验解线性方程组的基本迭代法实验

数值分析实验报告0 a 12K a1,n 1K a2,n 1UOM则有:第一步: Jacobi 迭代法a1n a2nM , 则有: A D L Uan 1,nAx bA A x D b 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 Ua 11a12 La1n a11Aa 21a 22 La2na22M MMMOan1a n2 Lanna11得到 Da22Oann由a 21 0MMOan 1,1 an 1,2L 0anna n1an2La n,na21LMMOan 1,1 an 1,2Lan1an2Lan,n 1a12K a1,n 1 a1n0 Ka 2,n 1a2nOM Man 1,n10令G (D L) U,则称G为Gauss-Seidel 迭代矩阵 f (D L) b由此可得 Gauss-Seidel 迭代的迭代格式如下:x (0) ,初始向量第三步SOR 迭代法w0ADLU1(D1 wL ((1 w)D wU ))(D 1wL)((1 w)D wU )www令Mw1(D wL), N1((1 w)D wU )则有:AMNwwAx bAMLWNM (M N )x b Mx Nx b x MNx M bNM,令Wf Mb带入 N 的值可有L W ((1 w)D wU)(D wL) 1((1 w)D wU) (D wL)f 1bw 1(D wL) 1b1(D wL) w称 L W 为 SOR 迭代矩阵,由此可得 SOR 迭代的迭代格式如下:x(0) ,初始向量二、算法程序Jacobi 迭代法的 M 文件: function [y,n]=Jacobi(A,b,x0,eps) %*************************************************%函数名称 Jacobi 雅克比迭代函数 %参数解释 A系数矩阵 % b常数项 %x0估计解向量x(k 1)Gx (k) f ,k 0,1,2,L(k 1)f,k 0,1,2,L误差范围%eps%返回值解向量%y迭代次数%n%函数功能实现线性方程组的Jacobi 迭代求解%*************************************************n=length(A);if nargin<3error(' 输入错误,最少要输入三个参数'); return;endif nargin==3eps=1e-6;endD=diag(diag(A));L=-tril(A,-1);U=-triu(A,1);M=D;N=L+U;B=M\N;f=M\b;if max(abs(eig(B)))>=1disp(' 谱半径大于等于1,迭代不收敛,无法进行'); return;endn=1;for i=1:1:nif sum(A(i,i)~=n)~=nerror(' 输入的A 矩阵的对角线元素不能为0'); return;endendy=B*x0+f;while norm(y-x0)>=eps&n<100x0=y;y=B*x0+f;n=n+1;endGauss-Seidel 迭代法的M 文件和function[y,n]=GaussSeidel(A,b,x0,eps) %************************************************* %函数名称GaussSeidel高斯赛德尔迭代函数系数矩阵%参数解释A常数项%b%x0估计解向量%eps误差范围%返回值解向量%y迭代次数%n%函数功能实现线性方程组的Jacobi 迭代求解n=1;for i=1:1:nif sum(A(i,i)~=n)~=n0');error(' 输入的A 矩阵的对角线元素不能为return;endendy=B*x0+f;while norm(y-x0)>=eps&n<100x0=y;y=B*x0+f;n=n+1;endSOR 迭代法的M 文件function [y,n]=SOR(A,b,x0,w,eps)%*************************************************%函数名称SOR松弛迭代函数系数矩阵%参数解释A松弛因子%w常数项%b估计解向量%x0误差范围%eps%返回值%y解向量迭代次数endend y=B*x0+f;while norm(y-x0)>=eps&n<100x0=y;y=B*x0+f;n=n+1;end三、数值计算1)首先编写如下程序实现输入大矩阵A: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;endendendfor i=1:1:20for j=1:1:20if abs(i-j)==2A(i,j)=-1/4;endend第一次给定初始向量为20 行一列的0,end右端面项向量b=20 行一列的1迭代误差要求0.005Jacobi 迭代法求解: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;endendendfor i=1:1:20for j=1:1:20if abs(i-j)==2A(i,j)=-1/4;endendend>> b=ones(20,1);>> x0=zeros(20,1);>> 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.4813 n =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)==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);>> 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.63280.57340.4816 n =7>>第二次给定初始向量为20 行一列的0 右端面项向量b=20 行一列的1.001 迭代误差要求0.005 Jacobi 迭代法求解:在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.50030.50030.50030.50030.50030.50030.50030.50020.49990.49780.48560.41467>>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)==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.48560.49780.49990.50030.50030.50030.50030.50030.50030.50030.50030.50030.50030.50030.50030.50000.49800.48580.4146n =5观察计算结果得到的序列可以看出其是收敛,在较少的迭代次数下即可的到满足误差要求的解。

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: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);>> 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.0831-5.4598-3.5651n =100第二次给定初始向量为20 行一右端列的0,面项向量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);>> 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.33330.33330.33330.33330.33330.33330.33330.33340.33310.33480.324621 / 220.387419分析结果:通过对迭代次数及其迭代结果的分析,我的得出的结论是松驰系数ω在SOR 迭代中起着相当重要的作用,不同的松驰系数ω,可能对迭代结果带来很大的影响,恰当的松驰系数ω可以加速收敛,得到较为良好的迭代结果,而不恰当的松驰系数ω选取,则可能会得导致无法获得理想的结果,甚至还可能影响到迭代的收敛性。

相关主题