椭圆形方的差分法.docVIP

数学与计算科学学院 实 验 报 告 实验项目名称 椭圆方程的差分格式 所属课程名称 微分方程数值解 实 验 类 型 综合型 实 验 日 期 2016-4-23 班 级 信计1301 学 号 201353100101 姓 名 向 溶 成 绩 一、实验概述： 【实验目的】 1、建立内点的五点差分格式。(取h=1) 2、建立包括边界点在内的五点差分格式方程组。 3、用雅克比迭代迭代法求解方程组。 4、计算结果（保留至小数点后4位）。 5、由计算结果，写出结论 【实验原理】 椭圆形方程的差分格式： 正方形区域中的Ladplace方程Dirichlet边值问题的差分模拟： 【实验环境】 Matlab R2010a 二、实验内容： 【实验方案】 在Matlab环境下分别用五点差分法编写解椭圆方程的算法，然后调试运行，表述结果。 【实验过程】（实验步骤、记录、数据、分析） 【实验结论（结果）y = 0 0.2000 0.4000 0.6000 0.8000 1.0000 1.2000 1.4000 1.6000 1.8000 2.0000 2.2000 2.4000 2.6000 2.8000 3.0000 3.2000 3.4000 3.6000 3.8000 4.0000 【实验小结结论Helmtz.m： function [u,x,y]= helmtz(f,g,bx0,bxf,by0,byf,D,M,N,tol,maxiter) x0=D(1); xf=D(2); y0=D(3); yf=D(4); dx=(xf-x0)/M; x=x0+[0:M]*dx; dy=(yf-y0)/N; y=y0+[0:N]*dy; %边界条件 for m=1:N+1 u(m,[1,M+1])=[feval(bx0,y(m)),feval(bxf,y(m))]; end for n=1:M+1 u([1,N+1],n)=[feval(by0,x(n)),feval(byf,x(n))]; end %边界的平均值作为初始值 bvaver=sum([sum(u(2:N,[1,M+1])),sum(u([1,N+1],2:M))]); u(2:N,2:M)=bvaver/(2*(M+N-2)); for i=1:N for j=1:M F(i,j)=feval(f,x(j),y(i)); G(i,j)=feval(g,x(j),y(i)); end end dx2=dx*dx; dy2=dy*dy; dxy2=2*(dx2+dy2); rx=dx2/dxy2; ry=dy2/dxy2; rxy=rx*dy2; for itr=1:maxiter for j=2:M for i=2:N u(i,j)=ry*(u(i,j+1)+u(i,j-1))+rx*(u(i+1,j)+u(i-1,j))+rxy*(G(i,j)*u(i,j)-F(i,j)); end end if itr1max(max(abs(u-u0)))tol break; end u0=u; end ex11.m: f=inline(0,x,y); g=inline(0,x,y); x0=0;xf=4;M=20; y0=0;yf=4;N=20; bx0=inline(exp(y)-cos(y),y); bxf=inline(exp(y)*cos(4)-exp(4)*cos(y),y); by0=inline(cos(x)-exp(x),x); byf=inline(exp(4)*cos(x)-exp(x)*cos(4),x); D=[x0 xf y0 yf]; maxiter=500; tol=1e-6; [u,x,y]= helmtz(f,g,bx0,bxf,by0,byf,D,M,N,tol,maxiter) clf,mesh(x,y,u),axis([0 4 0 4 -60 60