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

偏微分方程的数值解法上机习题七题目：求解边值问题：22sincoscossin,0,10,10,,.xyuexyxyxyGuxyG，，取步长1/64,1/128,hk作五点差分格式。用Jacobi迭代，Gauss-Seidel迭代和SOR迭代，（取opt）求解差分方程，以前后两次重合到小数点后四位的迭代值作为解的近似。比较三种解法的迭代次数以及差分解,1/64,1/128huxyh与精确解,sinsinxyuxyexy的精度。求解：1给定步长，建立方程组系数矩阵和常数项。2计算出迭代法的迭代矩阵、常数项及最佳松弛因子3进行迭代计算，输出迭代次数4求解精确解并作出精确解和迭代解对应图像5作出误差图像程序：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;ifiNbreakendenddisp([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）步长为1/64a)迭代次数请输入等分数n=641.Jacobi迭代的迭代次数是382.Guass-Seidel迭代的迭代次数是283.SOR超松弛法的迭代次数是26b)解的图像00.20.40.60.8100.510204060x精确解图像yu00.20.40.60.8100.51-0.500.51xJacobi迭代图像yu00.20.40.60.8100.51-0.500.51xGuass-Seidel迭代图像yu00.20.40.60.8100.51-0.500.51xSOR超松弛法图像yuc)误差图像00.20.40.60.8100.51-60-40-200xJacobi迭代误差图像yu00.20.40.60.8100.51-60-40-200xGuass-Seidel迭代误差图像yu00.20.40.60.8100.51-60-40-200xSOR超松弛误差图像yu（2）步长为1/32a)迭代次数请输入等分数n=321.Jacobi迭代的迭代次数是392.Guass-Seidel迭代的迭代次数是293.SOR超松弛法的迭代次数是27b)解的图像00.20.40.60.8100.510204060x精确解图像yu00.20.40.60.8100.51-2024xJacobi迭代图像yu00.20.40.60.8100.51-2024xGuass-Seidel迭代图像yu00.20.40.60.8100.51-2024xSOR超松弛法图像yuc)误差图像00.20.40.60.8100.51-60-40-200xJacobi迭代误差图像yu00.20.40.60.8100.51-60-40-200xGuass-Seidel迭代误差图像yu00.20.40.60.8100.51-60-40-200xSOR超松弛误差图像yu（3）步长为1/16a)迭代次数请输入等分数n=161.Jacobi迭代的迭代次数是382.Guass-Seidel迭代的迭代次数是283.SOR超松弛法的迭代次数是26b)解的图像00.20.40.60.8100.510204060x精确解图像yu00.20.40.60.8100.51-50510xJacobi迭代图像yu00.20.40.60.8100.51-50510xGuass-Seidel迭代图像yu00.20.40.60.8100.51-50510xSOR超松弛法图像yuc)误差图像00.20.40.60.8100.51-60-40-200xJacobi迭代误差图像yu00.20.40.60.8100.51-60-40-200xGuass-Seidel迭代误差图像yu00.20.40.60.8100.51-60-40-200xSOR超松弛误差图像yu结果分析1从迭代次数来看，Jacobi迭代法最高，Guass-Seidel迭代法次之，SOR超松弛法最少。步长JacobiGuass-SeidelSOR1/643828261/323929271/163828262从上面三种迭代法的三维图来看，当步长比较大时，4张图像的光滑度和精确度略微存在差异，随着步长的减小，等分区间增多，网格越来越密，光滑度和精确度也越来越高，基本与精确解的图像一致。3由误差图可以看出，三种迭代法误差的趋势是一致的，当步长大时，误差图像基本一致，随着步长的减小，迭代法的运算结果越来越精确。4由于步长为1/128时，矩阵过大，运行缓慢，故改为探究步长为1/16、1/32、1/64三组。