资源描述：
《有限差分法的Matlab程序.doc》由会员上传分享，免费在线阅读，更多相关内容在行业资料-天天文库。
1、有限差分法的Matlab程序（椭圆型方程）functionFD_PDE(fun,gun,a,b,c,d) %用有限差分法求解矩形域上的Poisson方程 tol=10^(-6); %误差界 N=1000; %最大迭代次数 n=20; %x轴方向的网格数 m=20; %y轴方向的网格数 h=(b-a)/n;%x轴方向的步长 l=(d-c)/m;%y轴方向的步长 fori=1:n-1 x(i)=a+i*h; end%定义网格点坐标 forj=1:m-1 y(j)=c+j*l; end%定义网格点坐标 u=zero
2、s(n-1,m-1);%对u赋初值 %下面定义几个参数 r=h^2/l^2; s=2*(1+r); k=1; %应用Gauss-Seidel法求解差分方程 whilek<=N %对靠近上边界的网格点进行处理 %对左上角的网格点进行处理 z=(-h^2*fun(x(1),y(m-1))+gun(a,y(m-1))+r*gun(x(1),d)+r*u(1,m-2)+u(2,m-1))/s; norm=abs(z-u(1,m-1)); u(1,m-1)=z;
3、 %对靠近上边界的除第一点和最后点外网格点进行处理 fori=2:n-2 z=(-h^2*fun(x(i),y(m-1))+r*gun(x(i),d)+r*u(i,m-2)+u(i+1,m-1)+u(i-1,m-1))/s; ifabs(u(i,m-1)-z)>norm; norm=abs(u(i,m-1)-z); end u(i,m-1)=z; end %对右上角的网格
4、点进行处理 z=(-h^2*fun(x(n-1),y(m-1))+gun(b,y(m-1))+r*gun(x(n-1),d)+r*u(n-1,m-2)+u(n-2,m-1))/s; ifabs(u(n-1,m-1)-z)>norm norm=abs(u(n-1,m-1)-z); end u(n-1,m-1)=z; %对不靠近上下边界的网格点进行处理 forj=m-2:-1:2 %对靠近左边界的网格点进行处理
5、 z=(-h^2*fun(x(1),y(j))+gun(a,y(j))+r*u(1,j+1)+r*u(1,j-1)+u(2,j))/s; ifabs(u(1,j)-z)>norm norm=abs(u(1,j)-z); end u(1,j)=z; %对不靠近左右边界的网格点进行处理 fori=2:n-2 z=(-h^2*fun(x(i),y(j))+u(i
6、-1,j)+r*u(i,j+1)+r*u(i,j-1)+u(i+1,j))/s; ifabs(u(i,j)-z)>norm norm=abs(u(i,j)-z); end u(i,j)=z; end %对靠近右边界的网格点进行处理 z=(-h^2*fun(x(n-1),y(j))+gun(b,y(j))+r*u(n-1,j+1)+r*u(n-1,j-
7、1)+u(n-2,j))/s; ifabs(u(n-1,j)-z)>norm norm=abs(u(n-1,j)-z); end u(n-1,j)=z; end %对靠近下边界的网格点进行处理 %对左下角的网格点进行处理 z=(-h^2*fun(x(1),y(1))+gun(a,y(1))+r
