姓名-学号-中国海洋大学偏微分方程的数值解法第七次作业

偏微分方程的数值解法上机习题七题目：求解边值问题：取步长作五点差分格式。用Jacobi迭代，Gauss-Seidel迭代和SOR迭代，（取）求解差分方程，以前后两次重合到小数点后四位的迭代值作为解的近似。比较三种解法的迭代次数以及差分解与精确解的精度。求解：给定步长，建立方程组系数矩阵和常数项。计算出迭代法的迭代矩阵、常数项及最佳松弛因子进行迭代计算，输出迭代次数求解精确解并作出精确解和迭代解对应图像作出误差图像程序：clc;clear;n=input('请输入等分数n=');h=1/n;%---------------------建立常数项b--------------------------fori=1:(n-1)forj=1:(n-1)b(j+(i-1)*(n-1),1)=2*pi^2*exp(pi*(i*h+j*h))*(sin(pi*i*h)*cos(pi*j*h)+cos(pi*i*h)*sin(pi*j*h));endend%---------------------建立系数矩阵A--------------------------j=0;fori=1:(n-1)^2ifi>=nj=j+1;A(i,j)=1/h^2;A(j,i)=1/h^2;endA(i,i)=-4/h^2;endj=0;fori=2:(n-1)^2j=j+1;ifmod(j,n-1)~=0A(j,i)=1/h^2;endend%--------------迭代法矩阵----------------D=diag(diag(A));L=tril(A,-1);R=triu(A,1);N=20000;%迭代最大次数%------------Jacboi迭代法迭代矩阵及常数项-----------------G=-inv(D)*(L+R);H=inv(D)*b;%------------Guass-Seidel迭代法迭代矩阵及常数项------------M=inv(D+L);R1=-M*R;M=M*b;%-----------SOR超松弛迭代法迭代矩阵、常数项及最佳松弛因子-------------I=speye((n-1)^2);B=I-inv(D)*A;q=eig(B);q=q(1);w=2/(1+sqrt(1-q^2));P=inv(D-w*(-L));Aw=P*(w*(-R)+(1-w)*D);P=w*P*b;%------------------------迭代法计算并作图-----------------------------str={'1.Jacobi迭代的迭代次数是','2.Guass-Seidel迭代的迭代次数是','3.SOR超松弛法的迭代次数是'};fork=1:3U=speye((n-1)^2,1);U0=zeros((n-1)^2,1);U0=sparse(U0);time=0;whilenorm(U-U0)>1e-4U0=U;switchkcase1U=G*U0+H;case2U=R1*U0+M;case3U=Aw*U0+P;endtime=time+1;ifi>Nbreakendenddisp([char(str{k}),num2str(time)]);fora=1:n-1forb=1:n-1switchkcase1Z_1(a,b)=U((n-1)*(a-1)+b,1);case2Z_2(a,b)=U((n-1)*(a-1)+b,1);case3Z_3(a,b)=U((n-1)*(a-1)+b,1);endendendend%-----------------------作方程解图需要---------------------------------x=1/n:1/n:(n-1)/n;y=x;[X,Y]=meshgrid(x,y);figure(1);%-------------------------精确解计算并作图-----------------------------Z=exp(pi*(X+Y)).*sin(pi*X).*sin(pi*Y);subplot(2,2,1);mesh(X,Y,Z);title('精确解图像');xlabel('x');ylabel('y');zlabel('u');%-----------------------迭代法图像--------------------------------subplot(2,2,2);mesh(X,Y,Z_1);title('Jacobi迭代图像');xlabel('x');ylabel('y');zlabel('u');subplot(2,2,3);mesh(X,Y,Z_2);title('Guass-Seidel迭代图像');xlabel('x');ylabel('y');zlabel('u');subplot(2,2,4);mesh(X,Y,Z_3);title('SOR超松弛法图像');xlabel('x');ylabel('y');zlabel('u');%------------------------迭代法误差图像----------------------------figure(2);subplot(2,2,1);mesh(X,Y,Z_1-Z);title('Jacobi迭代误差图像');xlabel('x');ylabel('y');zlabel('u');subplot(2,2,2);mesh(X,Y,Z_2-Z);title('Guass-Seidel迭代误差图像');xlabel('x');ylabel('y');zlabel('u');subplot(2,2,3);mesh(X,Y,Z_3-Z);title('SOR超松弛误差图像');xlabel('x');ylabel('y');zlabel('u');运行结果步长为1/64迭代次数请输入等分数n=641.Jacobi迭代的迭代次数是382.Guass-Seidel迭代的迭代次数是283.SOR超松弛法的迭代次数是26解的图像误差图像步长为1/32迭代次数请输入等分数n=321.Jacobi迭代的迭代次数是392.Guass-Seidel迭代的迭代次数是293.SOR超松弛法的迭代次数是27解的图像误差图像步长为1/16迭代次数请输入等分数n=161.Jacobi迭代的迭代次数是382.Guass-Seidel迭代的迭代次数是283.SOR超松弛法的迭代次数是26解的图像误差图像结果分析从迭代次数来看，Jacobi迭代法最高，Guass-Seidel迭代法次之，SOR超松弛法最少。步长JacobiGuass-SeidelSOR1/643828261/323929271/16382826从上面三种迭代法的三维图来看，当步长比较大