五点差分法(matlab)解椭圆型偏微分方程

用差分法解椭圆型偏微分方程

U(0,y)=si n(pi*y),U(2,y)=eA2si n( pi*y);

先自己去看一下关于五点差分法的理论书籍

Matlab 程序：

unction [p e u x y k]=wudianchafenfa(h,m,n,kmax,ep)

% g-s 迭代法解五点差分法问题

%kmax 为最大迭代次数

%m,n 为x,y 方向的网格数，例如(2-0 ) /0.01=200;

%e 为误差，p 为精确解

syms temp ;

u=zeros(n+1,m+1);

x=0+(0:m)*h;

y=0+(0:n)*h;

for (i=1:n+1)

u(i,1)=sin(pi*y(i)); u(i,m+1)=exp(1)*exp(1)*sin(pi*y(i)); end

for (i=1:n)

for ( j=1:m)

f(i,j)=(pi*pi-1)*exp(x(j))*sin(pi*y(i));

end

-(Uxx+Uyy)=(pi*pi-1)eAxsin(pi*y)

0

0=

end

t=zeros(n-1,m-1);

for (k=1:kmax)

for (i=2:n)

for ( j=2:m)

temp=h*h*f(i,j)/4+(u(i,j+1)+u(i,j-1)+u(i+1,j)+u(i-

1,j))/4;

t(i,j)=(temp-u(i,j))*(temp-u(i,j));

u(i,j)=temp;

end

end

t(i,j)=sqrt(t(i,j));

if (k>kmax)

break ;

end

if (max(max(t))

break ;

end

end

for (i=1:n+1)

for ( j=1:m+1)

p(i,j)=exp(x( j))*sin(pi*y(i));

e(i,j)=abs(u(i,j)-exp(x( j))*sin(pi*y(i)));

end

End

在命令窗口中输入：

[p e u x y k]=wudianchafenfa(0.1,20,10,10000,1e-6)

k=147

surf(x,y,u) ;

xlabel( ‘ x ');ylabel( ‘ y ' );zlabel( ‘u ');

Title( ‘五点差分法解椭圆型偏微分方程例 1' )

就可以得到下图

五点差分法解椭圆型偏懺分方程例

1

0 0

surf(x,y,p)

[p e u x y k]=wudia nchafe nfa(0.05,40,20,10000,1e-6)

[p e u x y k]=wudianchafenfa(0.025,80,40,10000,1e-6)

步长为1/3的曲面误差

0.06 "

□ Cl

为什么分得越小，误差会变大呢?

我们试试运行：

[p e u x y k]=wudia nchafe nfa(0.025,80

,40,10000,1e-8) K=2164

surf(x,y,e)

误差变小了吧还可以试试

[p e u x y k]=wudia nchafe nfa(0.025,80 K=3355

误差又大了一点,40,10000,1e-10)

再试试 [p e u x y k]=wudia nchafe nfa(0.025,80,40,10000,1e-11)

□ Cl k=3952

误差趋于稳定

总结：

最终的误差曲面

与网格数有关，也与设定的迭代前后两次差值（ ep, 看程序）有关；固定网格数，随着设定的迭代前后两次差值变小，误差由大比变小，中间有一个最小值，随着又增大一点，最后趋于稳定。

也许可以去研究一下那个误差最小的地方或者研究趋于稳定时的临界值。