二维问题的有限差分方法

西北农林科技大学实习报告

学院：理学院 专业年级：信计061 姓名：袁金龙 学号：15206012 课程：微分方程数值解 报告日期：2008-12-3

实习二、二维问题的有限差分方法

一）

实习问题：

二维经典初边值问题：

2

22

2,01(,0),01(0,)(1,)0,01x

u u te t t t u x x

x u t u t t ???=+<≤?????=<

二）算法描述：

由古典显格式利用线性组合的思想，可得古典隐格式：

1

111

11

2

k k

k k k j

j

j j

j k

j u u u u u a

f h

τ

+++++---+=+

上式可简化为：1

1

111(12)k k k k k

j j j j j j j j r u r u r u u f τ+++-+-++-=+

写成矩阵形式：

1111112122211111111

1

11200

0*00

12k k k

k k k k k k n n n n n k k k k k k

r r u u f r u u f u u f r r U F U F ττ+++-----+--++-?????

?????-?????

???????????=+?????????????????????????

?????-+?

?

=+=+G U U G G

L

O O O

O M M M O O L

由题知a=2

三）matlab 程序：

function []=erchubian(a,n,m) %n,m 分别表示空间与时间分成的分数 %r 表示网格比 %a 表示二阶导数的系数 h=1/n; g=1/m; r=(a*g)/(h*h);

for i=1:n

x(i)=h*i;

end

for i=1:m

t(i)=g*i;

end

for i=1:(n-1)

G(i,i)=1+2*r;

end

for i=1:(n-2)

G(i,i+1)=-r;

end

for i=2:n-1

G(i,i-1)=-r;

end

for i=2:n

for j=1:m

F(i,j)=f(x(i),t(j));

end

end

%求U(1:n,1)

for i=2:n

U(i,1)=x(i)*x(i);

end

for i=2:m

U(2:n,i)=inv(G)*U(2:n,i-1)+g*inv(G)*F(2:n,i-1); end

U

%绘制三维图

for i=1:n

for j=1:m

plot3(x(i),t(j),U(i,j),'H');

hold on;

end

end

四）图形显示的计算结果：

1