数值分析第一次作业matlab实验报告

几种线性方程组迭代算法的MATLAB实现和性能比较用有限差分方法（五点差分格式）求解正方形域上的Poisson方程边值问题用MATLAB语言编写算法程序求解线性方程组fAu的算法程序，采用下列方法，比较计算结果和算法性能，对计算结果给出讨论。一、算法实现解：由差分格式可得:2,1,1,,1,10,1,,0,14(1,)0(0,)ijijijijijijjNjiiNuuuuuhfijNuuuuijN，写成矩阵形式：Au=f11112222......,,...........NNNNAIvbIAvbAvbIIAvb其中：4114....,.....114NNiiiiAAR其中：11,12,1,121,22,2,21,2,,2211,12,1,121,22,2,221,2,,(,,...,),(,,...,),......,(,,...,)(,,...,),(,,...,),......,(,,...,)TTNNTNNNNNTTNNTNNNNNvuuuvuuuvuuubhfffbhfffbhfff(1)用Jacobi迭代法求解线性方程组fAufunction[u,k,er,t]=xsgs(n)%Jacobi迭代法0)1,()0,(),1(),0(1,0,2),(2222xuxuyuyuyxyxfyuxu　　%U表示方程组的解；h表示步长；A表示迭代矩阵；k表示迭代次数；n表示非边界点数%f表示线性方程组A*U=f的右端矩阵f；e表示允许误差界；er表示迭代误差%t表示计算时间tic;b(2:n+1,2:n+1)=(n+1)^(-2)*2;u=zeros(n+2,n+2);e=10^(-9);fork=1:1000%迭代求解er=0;ub=u;forj=2:n+1fori=2:n+1u(i,j)=(ub(i-1,j)+ub(i+1,j)+ub(i,j-1)+ub(i,j+1)+b(i,j))/4;er=er+abs(u(i,j)-ub(i,j));%估计当前误差endender=er/n^2;ifere,break;end%判断是否达到计算精度，如果达到则退出循环endtoc;t=toc;end计算结果：u=0000000000000.02560.04130.05080.05600.05770.05600.05080.04130.0256000.04130.06860.08590.09550.09860.09550.08590.06860.0413000.05080.08590.10880.12160.12580.12160.10880.08590.0508000.05600.09550.12160.13640.14120.13640.12160.09550.0560000.05770.09860.12580.14120.14620.14120.12580.09860.0577000.05600.09550.12160.13640.14120.13640.12160.09550.0560000.05080.08590.10880.12160.12580.12160.10880.08590.0508000.04130.06860.08590.09550.09860.09550.08590.06860.0413000.02560.04130.05080.05600.05770.05600.05080.04130.0256000000000000k=304er=9.7771e-10t=0.0140(2)用块Jacobi迭代法求解线性方程组fAufunction[u,k,er,t]=xsbgs(n)%块Jacobi迭代法%u表示方程组的解h表示步长；k表示迭代次数；n表示非边界点数；%f表示线性方程组A*u=f的右端矩阵f；q表示n+2维向量；a表示方程组系数矩阵的下对角线元素%b表示方程组系数矩阵的主对角线元素；c表示方程组系数矩阵的上对角线元素；d表示追赶法所求方程的右端向量%e表示允许误差界；er表示迭代误差；l表示系数矩阵A所分解成的下三角阵L中的下对角线元素l(i)；z表示系数矩阵A所分解成的上三角阵U中的主对角线元素z(i)tic;f=2*1/(n+1)^2*ones(n+2,n+2);a=-1*ones(1,n);b=4*ones(1,n);c=-1*ones(1,n);u=zeros(n+2,n+2);e=10^(-9);fork=1:1000%迭代求解er=0;ub=u;forj=2:n+1d(1:n)=f(2:n+1,j)+ub(2:n+1,j-1)+ub(2:n+1,j+1);x=zg(a,b,c,d);%用追赶法求解u(2:n+1,j)=x';er=er+norm(ub(:,j)-u(:,j),1);ender=er/n^2;ifere,break;end%判断是否达到计算精度，如果达到则退出循环endtoc;t=toc;end计算结果：u=0000000000000.02560.04130.05080.05600.05770.05600.05080.04130.0256000.04130.06860.08590.09550.09860.09550.08590.06860.0413000.05080.08590.10880.12160.12580.12160.10880.08590.0508000.05600.09550.12160.13640.14120.13640.12160.09550.0560000.05770.09860.12580.14120.14620.14120.12580.09860.0577000.05600.09550.12160.13640.14120.13640.12160.09550.0560000.05080.08590.10880.12160.12580.12160.10880.08590.0508000.04130.06860.08590.09550.09860.09550.08590.06860.0413000.02560.04130.05080.05600.05770.05600.05080.04130.0256000000000000k=163er=9.5903e-10t=0.4118(3)用共轭斜量法求解线性方程组fAufunction[x,k,t]=cg(n)%共轭向量法求解线性方程组a=zeros(n^2);fori=1:n^2%储存A矩阵a(i,i)=4;ifi+nn^2+1a(i,i+n)=-1;a(i+n,i)=a(i,i+n);endifmod(i,n)~=0a(i,i+1)=-1;a(i+1,i)=a(i,i+1);endendfori=1:n^2%储存b和初始解b(i,1)=0.02;x(i,1)=0;endk=0;r=b-a*x;p=r;q0=r'*r;tic;whileq01e-9k=k+1;w=a*p;t=q0/(p'*w);%求迭代步长x=x+t*p;r=r-t*w;q=r'*r;s=q/q0;p=r+s*p;%新的搜索方向q0=q;endtoc;t=toc;x=reshape(x,9,9);endx=0.02560.04130.05080.05600.05770.05600.05080.04130.02560.04130.06860.08590.09550.09860.09550.08590.06860.04130.05080.08590.10880.12160.12580.12160.10880.08590.05080.05600.09550.12160.13640.14120.13640.12160.09550.05600.05770.09860.12580.14120.14620.14120.12580.09860.05770.05600.09550.12160.13640.14120.13640.12160.09550.05600.05080.08590.10880.12160.12580.12160.10880.08590.05080.04130.06860.08590.09550.09860.09550.08590.06860.04130.02560.04130.05080.05600.05770.05600.05080.04130.0256k=12;t=0.0035二、结果分析和总结1、在相同的精度要求下，块Jacobi迭代法迭代次数与点块Jacobi迭代法相比明显要少。这是由于求解的方程组规模比点迭代法小得多，所以速度较快，迭代次数较少。2、在相同的精度要求下，共轭向量法与两种Jacobi迭代法相比更加迅速，迭代次数大大减少，从算法时间复杂度和空间复杂度来看都是非常简单的算法。