热传导方程(有限差分法之热传导方程的求解)

  • A+
所属分类:知识

热传导方程
偏微分方程的数值求解在数值分析中占有重要地位,很多科学技术问题的数值计算都包括了偏微分方程的数值解问题。本文主要讨论有限差分方法在偏微分方程数值求解中的应用。有限差分方法的应用主要集中在依赖于时间的问题(双曲型和抛物型)。
有限差分基本思想就是把连续的定解区域用有限个离散点构成的网格来代替,这些离散点称为网格的节点;把连续定解区域上的连续变量的函数用在网格上定义的离散变量函数来近似;把原方程和定解条件中的微商用差商来近似,积分用积分和来近似。然后再利用插值方法便可以从离散解得到定解问题在整个区域上的近似解【1】。
本文以一维热传导方程为算例,简单介绍有限差分在偏微分方程数值求解中的应用。
用有限差分求解偏微分方程的时候必须把连续问题进行离散化。为此首先要进行对求解区域做网格剖分。
利用有限差分来离散偏微分方程时,方法是多种多样的。而且对于同一偏微分方程可以建立各不相同的差分方程。针对上述算例,本文采用对t和x都采用向前差商,建立古典显示格式的差分方程来代替原偏微分方程。
稳定性分析
利用有限差分格式进行计算时是按照时间层逐层推进的。如果考虑二层差分格式,那么计算第i+1层上的值时,需要用到第i层上计算出来的结果值。而计算第i层时的舍入误差(包括i=0的情况,不过此时是由于初始数据不精确而引起的)必然会影响到第i+1层的值。从而就要分析这种误差传播情况。希望误差的影响不至于越来越大,达到掩盖差分格式的解的面貌,这便是所谓的稳定性【2】。
进一步,本文用矩阵方法来分析上述差分方程的稳定性,首先把上述差分方程写成矩阵形式。
此时矩阵A是一个对称矩阵,所以矩阵A满足谱半径条件即可。矩阵A的特征值为
在这里,本文取热传导系数为1,时间步长取0.01,迭代次数为90次,空间长度取1。通过MATLAB编程求解可得:
由上图可知数值解与精确解【3】误差并不大,几乎看不出数值解与精确解之间的差别,其最大的相对误差为0.8512%,因此可认为热传导方程的数值解是稳定收敛的。

数值解 精确解 相对误差(%)
0 0 0
0.1585 0.15990 0.8512
0.2929 0.29088 0.7088
0.3827 0.38006 0.7070
0.4143 0.41137 0.7079
0.3827 0.38006 0.7070
0.2929 0.29088 0.7088
0.1585 0.15990 0.8512
0 0 0

附MATLAB代码
clc
clear all
dx=0.125;dt=0.001;%确定空间步长和时间步长
L=1;T=1;a=1;mu=dt/dx^2;L_bc=91;%确定求解区域、热传导系数和迭代次数
x=0:dx:L;t=0:dt:T;
n=length(t);m=length(x);
%给定初值
u=zeros(n,m);
u(1,:)=sin(pi.*x/L);
for i=2:L_bc;
    for j=1:m
        %给定边界条件
        if j==1||j==m
           u(i,j)=0;
        else
            %迭代求解
          u(i,j)=u(i-1,j)+a*mu*(u(i-1,j+1)-2*u(i-1,j)+u(i-1,j-1));
        end
    end
end
%用三次样条曲线插值函数,对数值解进行插值作图并与精确解对比
y0=u(91,:);
xx0=0:0.01:1;
yy0=interp1(x,y0,xx0,'cubic');
subplot(1,2,1)
plot(xx0,yy0)
axis([0,1,0,0.5]);
x1=0:0.125:1;
y1=[0,0.1599,0.29546,0.38603,0.41784,0.38603,0.29546,0.1599,0];
yy1=interp1(x1,y1,xx0,'cubic');
%计算相对误差
y11=[0,15.99,29.546,38.603,41.784,38.603,29.546,15.99,0];
y00=100*y0;
wc=100*((y11(2:8)-y00(2:8))./y11(2:8));
xxx=0.125:0.125:0.875;
subplot(1,2,2)
plotyy(xx0,yy1,xxx,wc)
参考文献
【1】谢焕田,吴艳.拉普拉斯有限差分法的MATLAB实现[J].四川理工学院学报,2008,21(3):1-2
【2】陆金甫,关治.偏微分方程数值解法[M].北京,清华大学出版社.1987.
【3】史策.热传导方程有限差分法的MATLAB实现[J].咸阳师范学院学报,2009,24(4):4-5

热传导方程相关文章

发表评论

:?: :razz: :sad: :evil: :!: :smile: :oops: :grin: :eek: :shock: :???: :cool: :lol: :mad: :twisted: :roll: :wink: :idea: :arrow: :neutral: :cry: :mrgreen: