抛物形扩散方程的有限差分法与数值实例.docx
- 文档编号:4583465
- 上传时间:2022-12-07
- 格式:DOCX
- 页数:11
- 大小:368.35KB
抛物形扩散方程的有限差分法与数值实例.docx
《抛物形扩散方程的有限差分法与数值实例.docx》由会员分享,可在线阅读,更多相关《抛物形扩散方程的有限差分法与数值实例.docx(11页珍藏版)》请在冰豆网上搜索。
抛物形扩散方程的有限差分法与数值实例
偏微分方程数值解
所在学院:
数学与统计学院
课题名称:
抛物形扩散方程的有限差分法及数值实例
抛物形扩散方程的有限差分法及数值实例
1.1抛物型扩散方程
抛物型偏微分方程是一类重要的偏微分方程。
考虑一维热传导方程:
2
u
a2
x
其中a是常数,f(x)是给定的连续函数。
按照初边值条件的不同给法,可将(1.1.1)的定解分为两类:
第一,初值问题(Cauchy问题):
求足够光滑的函数ux,t,满足方程(1.1.1)和初始条件:
ux,0x,x(1.1.2)
第二,初边值问题(也称混合问题):
求足够光滑的函数ux,t,满足方程(1.1.1)和初始条件:
ux,0x,0xl(1.1.3)
及边值条件
u0,tul,t0,0tT(1.1.4)
假定fX和x在相应的区域光滑,并且于0,0,l,0两点满足相容条件,则上述问题有唯一的充分光滑的解。
1.2抛物线扩散方程的求解
F面考虑如下热传导方程
(1.2.1)
u(0.t)u(L,t)0
u(x,0)(x)
其中,0xl,0tT,a(常数)是扩散系数。
平行直线xXjjh,j0,1,,N和ttkk,k0,1,,M将矩形域
G0xl;0tT分割成矩形网格。
其中Xj,tk表示网格节点;Gh表示网格点(位于开矩形G中的网格节点)的集合;Gh表示位于闭矩形G中的网格节点的集合;h表示Gh-Gh网格边界点的集合。
uk表示定义在网点xj,tk处的待求近似解,0
现在对方程进行差分近似:
(一)向前差分格式
0
UjjXj,
kk
U0=Un=0
(1.2.3)
计算后得:
1
UjrUj1(1
2r)Uk
kr
V1fj
(1.2.4)
其中,r
各,j0,1,,N1,k
h
0,1,
M1。
显然,
这是一个四点显示格式,
每一层各个节点上的值是通过一个方程组求
k1k
UjUj
h2
kkk
Uj12UjUj1a
(fjf(Xj))
解到的。
方程组如下:
若记
则显格式(1.2.4)可写成向量形式
(1.2.6)
U
其中
1.3抛物型热传导方程数值算例
2
uu
tx2
初值条件:
u(Xj,0)sinxj,j0,1,...,M
对时间和空间进行分割,令M=40,N=1600,通过Matlab计算得到该方程的解析解,数值解以及相对误差如下:
图
(1)解析解的图像
图
(2)数值解的图像
绝对误差
x1D
丫00X
图(3)M=40,N=1000的相对误差的图像
我们取部分精确解和数值解进行比较,结果如表⑴
x
t
数值解
精确解
相对误差
0.1
0.1
0.1151
0.1152
4
1.170010
0.2
0.2
0.0815
0.0816
1.6580104
0.3
0.3
0.0418
0.0419
1.2752104
0.4
0.4
0.0183
0.0184
5
7.445610
0.5
0.45
0.0117
0.0118
5.3755105
0.6
0.35
0.0300
0.0301
4
1.067410
0.7
0.25
0.0684
0.0686
1.7410104
0.8
0.15
0.5084
0.5085
5
7.589710
0.9
0.05
0.3654
0.3654
1.7407105
表
(1)数值解与精确解的比较
由表
(1)我们可以看出,精确解和数值解的绝对误差在104以,因此可以得出,
在分割M=40,N=1600下,该有限差分万法对万程(1.3.1)是收敛和稳定的。
下面,我们比较在不同的分割下对有限差分算法精度的影响。
在扩散系数a1不变的情况下,讲时间和空间进行更加细密的分割,取
M50,N10000,其中,M表示空间上的分割,N表示时间上的分割。
观察数
值解与精确解在节点Xj,tk处的绝对误差值,如下图所示:
绝对误差
x10
-1
321U.8
O.
图(4)M=50,N=10000的相对误差的图像
由图(3)和图(4),两者在节点处的误差收敛分别是在104和105以,因此,可以得出的结论是:
在收敛围,随着时间和空间的分割越细,节点数越多,精确解和解析解之间的绝对误差也越小,有限差分法的算法精度也越高。
最后,我们比较网比r1/2以及r1时扩散方程的收敛情况。
当网比r1时,此时我们取M=10,N=50,这时,方程的数值解与解析解还
有相对误差图如下:
«画呂<迺報呂0gHN-0LHIAI(9)®
e®旨<卑<旨ogHN-OLHIAI(g)®
03
O
CO
对误差
x10
图(7)M=10,N=50的绝对误差的图像
此时,我们观察绝对误差发现,扩散方程(1.3.1)时不收敛不稳定的。
而前面
1
我们已经知道,至丽格比为r丄时,方程是收敛稳定的。
所以,我们可以验证,
2
11
当网比r—时稳定,当r—时不稳定。
22
[参考文献]
[1]荣华,播•微分方程数值解法[M]..高等教育.2009.1.
[2]王曰朋.偏微分方程数值解[0L].
wenku.baidu./link?
url=-UAAXRI8_V547zVS6UwenT8rGKsbwNf9CuDmh2qmsy5K8eQ32PzhYazZ_sWiHfz_Pj3LA7ufHH4O3t8NIIUCnXNUiyYhssqmNBeArsrLwQG
[3]未知.偏微分方程的Matlab解法[OL].
wenku.baidu./link?
url=WBIR0q9nO2YsNg6UOLkCQ4Z5QOCefDKNdEglrNR2TJ8VlxaJ9lkCrLlaE
DIJ_OHCLehf19UJU_B7PRe1skQTYaaLcAa1UWcwlv7GYsWNtG7
[4]周品,何正风.MATLAB数值分析.[M]..机械工业.2009.1.
附录:
L=1;
M=40;
N=1600;
alfa=1;
lambda=0.5;%网格比
%**********************************************%
h=L/M;%空间步长
x=0:
h:
L;
x=x';
tao=lambda*hA2/alfa;%时间步长
tm=N*tao;%热传导的总时间tm
%tm=0.1;
t=0:
tao:
tm;
t=t';
%计算初值和边值
U=zeros(M+1,N+1);
U(:
1)=sin(pi*x);
U(1,:
)=0;
U(M+1,:
)=0;
%*************用差分法求出温度U,与杆长L,时间T的关系****************%
fork=1:
N
j=2;
whilej<=M
U(j,k+1)=lambda*U(j+1,k)+(1-2*lambda)*U(j,k)+lambda*U(j-1,k);
j=j+1;
end
end
length(U);%***************
fori=1:
N+1
X(:
i)=x;
endforj=1:
M+1
Y(j,:
)=t;
endmesh(X,Y,U);
legend('数值解');xlabel('X');ylabel('T');zlabel('U');
z=zeros(M+1,N+1);forj=1:
M+1
fork=1:
N+1z(j,k)=exp(-pi*pi*t(k))*sin(pi*x(j));end
end
%mesh(x,t,z')
legend('解析解');
xlabel('X');
ylabel('T');
zlabel('Z');
forj=1:
M+1
fork=1:
N+1
y(j,k)=abs(z(j,k)-U(j,k));
end
endmesh(x,t,y');
legend('绝对误差');xlabel('X');
ylabel('Y');
zlabel('error');
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 抛物形 扩散 方程 有限 差分法 数值 实例