微分方程数值解法实验2.docx

1、微分方程数值解法实验2一、实验概述：【实验目的】熟练掌握求解方程第一边值问题的五点差分格式。【实验原理】1.取定沿X轴和Y轴的步长h1和h2，作两族与坐标轴平行的直线：，两族直线的交点称为网点或节点，记为，假定为正则内点。沿x,y方向分别用二阶中心差商代替,则得，由于差分方程中只出现u在（i,j）及其四个邻点上的值，故称为五点差分格式。2.迭代法设线性方程组 （1）的系数矩阵可逆且主对角元素均不为零，令，并将分解成 （2）从而（1）可以写成，令，其中 （3）以为迭代矩阵的迭代法 （4）称为迭代法，用向量的分量来表示，（4）为 （5）其中为初始向量。3.迭代法把矩阵分解成 （6）其中，分别为的主

2、对角元除外的下三角和上三角部分，于是方程组（1）可以写成，即，其中 （7）以为迭代矩阵构成的迭代法 （8）称为迭代法，用向量表示的形式为4.迭代【实验环境】MATLAB R2014a二、实验内容：【实验方案】1.求解边值问题取步长h=k=1/64,1/128,做五点差分格式。用jacobi迭代，Gauss-Seidel迭代和SOR迭代（取）求解差分方程，以前后两次重合到小数点后四位的迭代值作为解的近似。比较三种解法的迭代次数以及差分解与精确解的精度。2.在软件上编写代码；3.运行得出结论。【实验过程】（实验步骤、记录、数据、分析）1.分析问题。2.在上编写代码，调试。3.运行代码，发现问题后改

3、进代码。4.分别输入n值为64和128，得出结论。【实验结论】（结果）n=64时:1.Jacobi迭代迭代次数i =78112.Guass-Seidel迭代迭代次数i =42113.SOR超松弛法迭代次数i =121n=128时：1.Jacobi迭代迭代次数i =200012.Guass-Seidel迭代迭代次数i =156603.SOR超松弛法迭代次数i =242【实验小结】（收获体会）通过本次实验，掌握了求解方程第一边值问题的五点差分格式，比较了三种解法的迭代次数以及差分解与精确解的精度，更加熟悉的掌握了利用MATLAB求解数学问题的方法。三、指导教师评语及成绩：评 语评语等级优良中及格不

4、及格1.实验报告按时完成,字迹清楚,文字叙述流畅,逻辑性强2.实验方案设计合理3.实验过程（实验步骤详细,记录完整,数据合理,分析透彻）4实验结论正确. 成 绩： 指导教师签名： 批阅日期：附录1：源 程 序clc;n=input(请输入等分数 n=);h=1/n;A=zeros(1);A=sparse(A);b=zeros(n-1)2,1);b=sparse(b);u=zeros(n-1)2,1);u=sparse(u);I=speye(n-1)2);D=zeros(n-1)2,1);D=sparse(D);for i=1:(n-1) for j=1:(n-1) xi=i*h; yj=j*h

5、; b(j+(i-1)*(n-1),1)=2*pi2*exp(pi*(xi+yj)*(sin(pi*xi)*cos(pi*yj)+cos(pi*xi)*sin(pi*yj); u(j+(i-1)*(n-1),1)=exp(pi*(xi+yj)*sin(pi*xi)*sin(pi*yj); %精确解 endend %-建立系数矩阵A-j=0;for i=1:(n-1)2 if(i=n) j=j+1; A(i,j)=1/h2; A(j,i)=1/h2; endendfor i=1:(n-1)2A(i,i)=-4/h2;endj=0;for i=2:(n-1)2 j=j+1; if(mod(j,n-

6、1)=0) A(i,j)=1/h2; A(j,i)=1/h2; endendA;D=diag(diag(A);L=tril(A,-1);R=triu(A,1);L=sparse(L);R=sparse(R);Z=speye(1);x=1/n:1/n:(n-1)/n;y=1/n:1/n:(n-1)/n;X,Y=meshgrid(x,y);Z1=exp(pi*(X+Y).*sin(pi*X).*sin(pi*Y);subplot(2,2,1);mesh(X,Y,Z1); %精确解图像title(精确解图像);xlabel(x);ylabel(y);zlabel(u); %-Jacboi迭代使用-G

7、=-inv(D)*(L+R);H=inv(D)*b;%-Guass-Seidel迭代使用-M=inv(D+L);M=sparse(M);R1=-M*R;R1=sparse(R1);M=M*b;%-SOR超松弛法使用-Aw=zeros(1);Aw=sparse(Aw);B=zeros(1);B=sparse(B);B=I-inv(D)*A;q=eig(B);q=q(1);w=2/(1+sqrt(1-q2);P=inv(D-w*(-L);Aw=P*(w*(-R)+(1-w)*D);P=w*P*b; %-for k=1:3 switch k%-Jacobi迭代-case 1 U=speye(n-1)

8、2,1); U0=zeros(n-1)2,1); U0=sparse(U0); disp(1.Jacobi迭代); i=0; while norm(U-U0)1e-4 U0=U; U=G*U0+H; i=i+1; if(i20000) break end end disp(迭代次数); i for a=1:n-1 for b=1:n-1 Z(a,b)=U(n-1)*(a-1)+b,1); end end subplot(2,2,2); mesh(X,Y,Z); title(Jacobi迭代图像); xlabel(x); ylabel(y); zlabel(u); %-Guass-Seidel迭

9、代- case 2 U=speye(n-1)2,1); U0=zeros(n-1)2,1); U0=sparse(U0); disp(2.Guass-Seidel迭代); i=0; while norm(U-U0)1e-4 U0=U; U=R1*U0+M; i=i+1; if(i20000) break end end disp(迭代次数); i for a=1:n-1 for b=1:n-1 Z(a,b)=U(n-1)*(a-1)+b,1); end end subplot(2,2,3); mesh(X,Y,Z); title(Guass-Seidel迭代图像); xlabel(x); yl

10、abel(y); zlabel(u); %-SOR超松弛法-case 3 U=speye(n-1)2,1); U0=zeros(n-1)2,1); U0=sparse(U0); disp(3.SOR超松弛法); i=0; while norm(U-U0)1e-1 U0=U; U=Aw*U0+P; i=i+1; if(i20000) break end end disp(迭代次数); i for a=1:n-1 for b=1:n-1 Z(a,b)=U(n-1)*(a-1)+b,1); end end subplot(2,2,4); mesh(X,Y,Z); title(SOR超松弛法图像);

11、xlabel(x); ylabel(y); zlabel(u); endend附录2：实验报告填写说明 1实验项目名称：要求与实验教学大纲一致。2实验目的：目的要明确，要抓住重点，符合实验教学大纲要求。3实验原理：简要说明本实验项目所涉及的理论知识。4实验环境：实验用的软、硬件环境。5实验方案（思路、步骤和方法等）：这是实验报告极其重要的内容。概括整个实验过程。对于验证性实验，要写明依据何种原理、操作方法进行实验，要写明需要经过哪几个步骤来实现其操作。对于设计性和综合性实验，在上述内容基础上还应该画出流程图、设计思路和设计方法，再配以相应的文字说明。对于创新性实验，还应注明其创新点、特色。6实验过程（实验中涉及的记录、数据、分析）：写明具体实验方案的具体实施步骤，包括实验过程中的记录、数据和相应的分析。7实验结论（结果）：根据实验过程中得到的结果，做出结论。8实验小结：本次实验心得体会、思考和建议。9指导教师评语及成绩：指导教师依据学生的实际报告内容，给出本次实验报告的评价。