資源描述：
《有限差分法的Matlab程序.doc》由會員上傳分享，免費在線閱讀，更多相關內(nèi)容在行業(yè)資料-天天文庫。
1、有限差分法的Matlab程序（橢圓型方程）functionFD_PDE(fun,gun,a,b,c,d)???%用有限差分法求解矩形域上的Poisson方程???tol=10^(-6);?%誤差界???N=1000;?%最大迭代次數(shù)???n=20;?%x軸方向的網(wǎng)格數(shù)???m=20;?%y軸方向的網(wǎng)格數(shù)???h=(b-a)/n;%x軸方向的步長???l=(d-c)/m;%y軸方向的步長???fori=1:n-1???????x(i)=a+i*h;???end%定義網(wǎng)格點坐標???forj=1:m-1???????y(j)=c+j*l;???end%定義網(wǎng)格點坐標???u=zero
2、s(n-1,m-1);%對u賦初值???%下面定義幾個參數(shù)???r=h^2/l^2;???s=2*(1+r);???k=1;???%應用Gauss-Seidel法求解差分方程???whilek<=N???????%對靠近上邊界的網(wǎng)格點進行處理???????????%對左上角的網(wǎng)格點進行處理???????????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、????%對靠近上邊界的除第一點和最后點外網(wǎng)格點進行處理???????????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???????????%對右上角的網(wǎng)格
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;???????%對不靠近上下邊界的網(wǎng)格點進行處理???????????forj=m-2:-1:2???????????????%對靠近左邊界的網(wǎng)格點進行處理?????
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;???????????????%對不靠近左右邊界的網(wǎng)格點進行處理???????????????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???????????????%對靠近右邊界的網(wǎng)格點進行處理???????????????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???????%對靠近下邊界的網(wǎng)格點進行處理???????????%對左下角的網(wǎng)格點進行處理???????????z=(-h^2*fun(x(1),y(1))+gun(a,y(1))+r
