传热学的数值解法

导热问题的数值求解方法

数值解法的基本思想是用空间和时间区域内有限个离散点（称为节点）上温度的近似值，代替物体内实际的连续温度分布，然后由导热方程和边界条件推导出各节点温度间的相互关系的代数方程组（称为离散方程），求解此方程组，得到节点上的温度值，此即物体中温度场的解。只要节点分布的足够稠密，数值解就有足够的精度。求解导热问题的数值方法有有限差分法及有限元法，近几年又发展了边界元法和有限分析法。数值方法适用于求解各种导热问题，不管物体的几何形状有多复杂，不管线性或非线性问题，都能使用。由于计算机的飞速发展，计算技术软件发展也很快，数值方法的的地位越来越重要。

1 数值求解的基本思路及稳态导热内节点离散方程的建立

一、 解法的基本思路

1、基本思路：数值解法的求解过程可用框图4-1表示。

由此可见：

1）物理模型简化成数学模型是基础；

2）建立节点离散方程是关键；

3）一般情况微分方程中，某一变量在某一坐标方向所

需边界条件的个数等于该变量在该坐标方向最高阶导

数的阶数。

二、稳态导热中位于计算区域内部的节点离散方程的

建立方法

1、基本方法

方法：①泰勒级数展开法；②热平衡法。

1）泰勒级数展开法

如图4-3所示，以节点(m,n)处的二阶偏导数为例，

对节点(m+1,n)及(m-1,n)分别写出函数t 对(m,n)点的泰

勒级数展开式：

对(m+1,n)：

+∂∂∆+∂∂∆+∂∂∆+∂∂∆+=+4

44333,222,,,12462x t x x t x x t x x t x t t n m n m n m n m （a ）

对（m-1,n ）：

+∂∂∆+∂∂∆-∂∂∆+∂∂∆-=-44

4333,222,,,12462x t x x t x x t x x t x

t t n m n m n m n m （b ）

（a ）+（b ）得： +∂∂∆+∂∂∆+=+-+444,222,,1,1122x t x x t x t t t n m n m n m n m 变形为n m x t

,22∂∂的表示式得：

n m x t

,22∂∂)(0222,1,,1x x t t t n

m n m n m ∆+∆+-=-+ 上式是用三个离散点上的值计算二阶导数n m x t ,2

2∂∂的严格表达式，其中：

)(02x ∆―― 称截断误差，误差量级为2x ∆

在数值计算时，用三个相邻节点上的值近似表示二阶导数的表达式即可，则相应的略去)(02x ∆。于是得：

n m t x t

,2

2∂∂2,1,,12x t t t n m n m n m ∆+-=-+ （4-1a ）

同理： 21

,,1,,222y t t t y t

n m n m n m n m ∆+-=∂∂-+ （4-1b ）

根据导热问题的控制方程(导热微分方程)02222=∂∂+∂∂y t x t 得：

2,1,,12x

t t t n

m n m n m ∆+--+0221,,1,=∆+-+-+y t t t n m n m n m (4-2)

若△x=△y 则有： )(411,1,,1,1,-+-++++=

n m n m n m n m n m t t t t t

2) 平衡法：

其本质是傅里叶导热定律和能量守恒定律的体现。对每个元体，可用傅里叶导热定律写出其能量守恒的表达式。如图4-3所示，元体在垂直纸面方向取单位长度，通过元体界面(w,e,n,s)所传导的热流量可以对有关的两个节点根据傅里叶定律写出：

从节点(m-1,n)通过界面w 传导到节点(m,n)的热流量： x t t y n

m n m w ∆-∆=Φ-,,1λ （a ）

同理：通过界面e,n,s 传导给节点（m,n ）的热流量：

x t t y n m n m e ∆-∆=Φ+,,1λ （b ）

y t t x n m n m n ∆-∆=Φ+,1,λ （c ）

y t t x n

m n m s ∆-∆=Φ-,1,λ （d ）

对元体(m,n).根据能量守恒定律可知：

0=Φ+Φ+Φ+Φs n e w （4-3）

其中，规定：导入元体（m,n ）的热流量为正；

导出元体（m,n ）的热流量为负。

说明：① 上述分析与推导是在笛卡儿坐标系中进行的；

② 热平衡法概念清晰，过程简捷；

③ 热平衡法与§2—2建立微分方程的思路与过程一致，但不同的是前者是有限大小的元体，后者是微元体。

2 稳态导热边界节点离散方程的建立及代数方程的求解

（1） 对于第一类边界条件的导热问题，所有内节点的离散方程组成一个封闭的代数方程组，即可求解；

（2） 第二类或第三类边界条件的导热问题，所有内节点的离散方程组成的代数方程组是不封闭的，因未知边界温度，因而应对位于该边界上的节点补充相应的代数方程，才能使方程组封闭，以便求解。

一、用热平衡法导出典型边界点上的离散方程

假设物体具有内热源 ⋅Φ (不必均匀分布)，而且边界上有向该元体传递的热流密度w q ：

1、位于平直边界上的节点

如图所示4-4边界节点(m,n)只能代表半个元体，若边界上有向该元体传递的热流密度为w q ，据能量守恒定律对该元体有：

y x t t n

m n m ∆∆--,,1λ2,1,x y t t n m n m ∆⋅∆-++λ2

,1,x y t t n m n m ∆⋅∆-+-λ02,=∆+Φ∆∆+⋅w n m yq y x （4-4a ）

若y x ∆=∆时，则：

)22(41,2,1,1,1,λλφw n m n m n m n m n m xq x t t t t ∆+∆+++=⋅

-+- （4-4b ）

2、外部角点

如图4-5所示，二维墙角计算区域中，

该节点外角点仅代表1/4个以y x ∆∆,为边长

的元体。假设边界上有向该元体传递的热流

密度为w q ，则据能量守恒定律得其热平衡

式为： 2

,,1y x t t n m n m ∆⋅∆--λ2,1,x y t t n m n m ∆⋅∆-+-λ024,=∆+∆+Φ∆∆+⋅w n m q y x y x （4-5a ）

若y x ∆=∆时，则：

)22(21,21,,1,λλφw n m n m n m n m xq x t t t ∆+∆++=⋅

-- （4-5b ）

3、内部角点：

如图4-5所示内部角点代表了3/4个以y x ∆∆,为边界长的元体。

同理得： )22322(61,2,11,1,,1,λλw n m n m n m n m n m n m xq x t t t t t ∆+Φ∆++++=⋅

+-+- （4-6b ） 二、代数方程的求解方法

1、直接解法：通过有限次运算获得精确解的方法，如：矩阵求解，高斯消元法。

2、迭代法：先对要计算的场作出假设（设定初场），在迭代计算中不断予以改进，直到计算前的假定值与计算结果相差小于允许值为止的方法，称迭代计算收敛。目前应用较多的是：

1）高斯——赛德尔迭代法：每次迭代计算，均是使用节点温度的最新值。

2）用雅可比迭代法：每次迭代计算，均用上一次迭代计算出的值。

3 非稳态导热问题的数值解法

由前可知：非稳态导热和稳态导热二者微分方程的区别在于控制方程中多了一个非稳态项，其中扩散项的离散方法与稳态导热一样。