稳态热传导问题有限元法.docx
《稳态热传导问题有限元法.docx》由会员分享，可在线阅读，更多相关《稳态热传导问题有限元法.docx（20页珍藏版）》请在冰豆网上搜索。
稳态热传导问题有限元法
6.稳态热传导问题的有限元法
本章的内容如下：
6.1热传导方程与换热边界
6.2稳态温度场分析的一般有限元列式
6.3三角形单元的有限元列式
6.4温度场分析举例
6.1热传导方程与换热边界
在分析工程问题时，经常要了解工件内部的温度分布情况，例如发动机的工作温度、金
属工件在热处理过程中的温度变化、流体温度分布等。
物体内部的温度分布取决于物体内部
的热量交换，以及物体与外部介质之间的热量交换，一般认为是与时间相关的。
物体内部的
热交换采用以下的热传导方程（Fourier方程）来描述，
TTTT
cxyzQ（6-1）
txxyyzz
式中为密度，kg/m3；c为比热容，J/（kgK）；x,y,z为导热系数，wmk；T为温度，C；t为时间，s；Q为内热源密度，w/m3。
对于各向同性材料，不同方向上的导热系数相同，热传导方程可写为以下形式，
T2T2T2T
c22—Q（6-2）
txyz
除了热传导方程，计算物体内部的温度分布，还需要指定初始条件和边界条件。
初始条
件是指物体最初的温度分布情况，
Tt0T0x,y,z（6-3）
边界条件是指物体外表面与周围环境的热交换情况。
在传热学中一般把边界条件分为三
类。
1）给定物体边界上的温度，称为第一类边界条件。
物体表面上的温度或温度函数为已知，
TsTs
或TsTs（x,y,z,t）（6-4）
2）给定物体边界上的热量输入或输出，称为第二类边界条件。
已知物体表面上热流密度，
（x
Tnxx
T
yny
y
T
znz）
z
s
qs
T
T
T、
或
（x
nx
yny
znz）
s
qs（x,y,z,t）
（6-5）
x
y
z
3）给定对流换热条件，称为第三类边界条件。
4）
物体与其相接触的流体介质之间的对流换热系数和介质的温度为已知。
如果边界上的换热条件不随时间变化，物体内部的热源也不随时间变化，在经过一定时
间的热交换后，物体内各点温度也将不随时间变化，即
丄0t
这类问题称为稳态（Steadystate）热传导问题。
稳态热传导问题并不是温度场不随时间我们不关心物体内部的温度场如何从初始状态过渡
（Transient）
（6-7）
Laplace方程,
考虑物体不包含内热源的情况，各向同性材料中的温度场满足
（6-9）
2t2t2t0
~2220
xyz
在分析稳态热传导问题时，不需要考虑物体的初始温度分布对最后的稳定温度场的影响，因此不必考虑温度场的初始条件，而只需考虑换热边界条件。
计算稳态温度场实际上是
求解偏微分方程的边值问题。
温度场是标量场，将物体离散成有限单元后，每个单元结点上只有一个温度未知数，比弹性力学问题要简单。
进行温度场计算时有限单元的形函数与弹性力学问题计算时的完全一致，单元内部的温度分布用单元的形函数，由单元结点上的温度来
确定。
由于实际工程问题中的换热边界条件比较复杂，在许多场合下也很难进行测量，如何
定义正确的换热边界条件是温度场计算的一个难点。
6.2稳态温度场分析的一般有限元列式
在前面我们已经介绍了有限元方法可以用来分析场问题，稳态温度场计算是一个典型的
场问题。
我们可以采用虚功方程建立弹性力学问题分析的有限元格式，推导出的单元刚度矩
阵有明确的力学含义。
在这里，介绍如何用加权余量法（WeightedResidualMethod）建立稳
态温度场分析的有限元列式。
微分方程的边值问题，可以一般地表示为未知函数u满足微分方程组，
未知函数u还满足边界条件，
B（u）
B1（u）
B2（u）0（在边界上）
（6-11）
如果未知函数u是上述边值问题的精确解，则在域中的任一点上u都满足微分方程
（6-10），在边界的任一点上都满足边界条件（6-11）。
对于复杂的工程问题，这样的精确解
往往很难找到，需要设法寻找近似解。
所选取的近似解是一族带有待定参数的已知函数，般表示为
n
（6-12）
uuNiaiNa
i1
其中ai为待定系数，Ni为已知函数，被称为试探函数。
试探函数要取自完全的函数序列，是线性独立的。
由于试探函数是完全的函数序列，任一函数都可以用这个序列来表示。
采用这种形式的近似解不能精确地满足微分方程和边界条件，所产生的误差就称为余
量。
微分方程（6-10）的余量为，
RA（Na）（6-13）
边界条件（6-11）的余量为，
RB（Na）（6-14）
选择一族已知的函数，使余量的加权积分为零，强迫近似解所产生的余量在某种平均意义上等于零，
TT—
WjRdWjRd0（6-15）
Wj和Wj称为权函数，通过公式（6-15）可以选择待定的参数ai。
这种采用使余量的加权积分为零来求得微分方程近似解的方法称为加权余量法。
对权函
数的不同选择就得到了不同的加权余量法，常用的方法包括配点法、子域法、最小二乘法、力矩法和伽辽金法（Galerkinmethod）。
在很多情况下，采用Galerkin法得到的方程组的系
数矩阵是对称的，在这里也采用Galerkin法建立稳态温度场分析的一般有限元列式。
在
Galerkin法中，直接采用试探函数序列作为权函数，取WjNj，WjNj。
下面用求解二阶常微分方程为例，说明Galerkin法（参见，王勖成编著“有限元法基本
原理和数值方法”的1.2.3节）。
d2u
dx2
例，求解二阶常微分方程
（0x1）
边界条件：
当x0时，u0;当x1时，u0。
取两项近似解：
Nix（1x）
2
N2x（1x）
〜
uN1a1N2a2a1x（1x）a2x（1x）
W1N1,W2N2
由公式（6-15）可以得到两个加权积分方程,
积分后可以得到一个二元一次方程组，解得,
a10.1924,a20.1707
近似解为，〜x（1x）（0.19240.1707x）
sinx
该方程的精确解为，ux
sin1
近似解与精确解的结果比较见表6-1，
表6-1近似解与精确解比较
x=0.25
x=0.5
x=0.75
sinxux
sin1
0.04401
0.06975
0.06006
〜x（1x）（0.19240.1707x）
0.04408
0.06944
0.06008
假定单元的形函数为，
[N][N12…心]
单元结点的温度为，
{T}e[T1T2...Tn]T
单元内部的温度分布为,
T[N]{T}e
.维问
以二维问题为例，说明用Galerkin法建立稳态温度场的一般有限元格式的过程。
题的稳态热传导方程为,
（6-16a）
第一类换热边界为
Ts
s
第二类换热边界条件为,
（6-16b）
TT
xnxynyqs
xy
第三类边界条件为，
（6-16c）
T—nxx
Tnyh（TfTs）
y
（6-16d）
e
T
T
W1[-
（x
）（
y
）Q]d
0
x
xy
y
由分部积分得，
T、
W1,
T、
T、
（w1
x）
：
（x-
）
W1（x■
）
x
x
x
x
x
x
T、
W1,
T、
T、
（w1
y）
1（y
）
W1（y
）
y
y
y
y
y
y
在一个单元内的加权积分公式为,
（6-17）
应用Green定理，一个单元内的加权积分公式写为,
x（
-（
y
Tnxx
（6-18）
Tny）d0
y
采用Galerkin方法，选择权函数为,
W!
Ni
将单元内的温度分布函数和换热边界条件代入
（6-18）式，单元的加权积分公式为,
6[旦（x
x
[N]）旦（
xy
y3）]{T}edy
NiQd
e
2Mqsd
（6-19）
Mh[N]{T}ed
MhTfd0
换热边界条件代入后,
在（6-19）式内相应出现了第二类换热边界项
e
3Niqsd，第三
3
类换热边界项Nih[N]{T}ed
3
e
3NihTfd，但没有出现与第一类换热边界对应的
3
项。
这是因为，采用M作为权函数，第一类换热边界被自动满足。
写成矩阵形式有,
型）（叫
xy
et
2[N]Tqsd
』）]{T}edy
Te
h[N][N]{T}d
3[N]ThTfd0
n个结点的温度Ti。
按有限元格式将
公式（6-20）是n个联立的线性方程组，可以确定
（6-20）表示为，
[K]e{T}e{P}e
（6-21）
其中矩阵[K]e为单元的导热矩阵或称为温度刚度矩阵,称为单元的温度载荷向量或热载荷向量（热矩阵[K]e和温度载荷向量
{T}e为单元的结点温度向量，{P}e
Thermalloadvector）。
对于某个特定单元，单元导
{P}e的元素分别为，
NiQd
在整个物体上的加权积分方程是单元积分方程的和,
根据单元结点的局部编号与整体编号的关系，直接求和得到整体刚度矩阵，整体方程组为,
[K]{T}{P}
6.3三角形单元的有限元列式
（心亠）
图6-1三角形单元
回顾第三章的内容可以发现，与计算弹性力学平面问题时所采用的方法一样,
场问题计算中所采用的三角形单元可以使用相同的形函数，
壬（aibixCiy）
2A
±（ajbjXCjy）2A（ambmXCmy）
2A
二维温度
Ni
Nj
Nm
ai
aj
Xjym
Xmyj
bi
yj
ym
Ci
Xm
Xj
XmYi
XiYm
bj
ym
yi
Cj
Xi
Xm
Xiyj
Xjyi
bm
yi
yj
Cm
Xj
Xi
Xi
yi
Xj
yj
T
Xm
ym
2A
单元内的温度分布用结点上的温度值表示为,
T[NiNj
Nm]Tj
Tm
（6-25）
在三角形单兀上,
采用
Galerkin
法可得，
A[N]TL（
AX
xT）
X
（y_y
T
）Q]dA0y
（6-26）
Ti
-）
X
[N]tt
假定单元内的导热系数为常数,
[N]TT
X
dA
XX
1
4A
bibj[bibm
bj
bm]
TiTjTm
dA
4A
b2
bibj
bibm
bibj
b2
bjbm
bibm
bjbm
bm
Ti
Tj
Tm
（6-27）
[N]t
-）
y
（[N]y
-）
y
[N]tt
y
[N]TT
y
dAy
1
4A
C
Cj[Ci
cm
Cj
Cm]
Ti
Tj
Tm
dA
4ACiCj
CiCm
GC
2
Cj
CjCm
CiCm
CjCm
2
Cm
Ti
Tj
Tm
（6-28）
单元的刚度矩阵为,
[K]e/
4
bi2
bibj
bibj
2
Ci
CiCj
CiCm
bibj
bj2
bjbm
y
GCj
2
Cj
CjCmc；
bibm
bjbm
b；
4
CiCm
CjCm
X
Ni
—1
TZ-
Qa
QA
[N]TQdA
A
NjdA
1
3
Nm
1
由Green公式可得
A[（[N]T
XT）
（[N]
TyT]dA
AX
X
y
y
■-（[N]Tx
Tn
[N]Ty-
T
ny）dS
s
X
y
显然，单元的导热矩阵是对称的。
如果单元的内部热源为常数，由内部热源产生的温度载荷项为,
方便起见，把换热边界统一表示为第三类换热边界，
（6-29）
（6-30）
A[—（[N]TX丄）一（[N]Ty丄]dA
Axxyy（6-30）
TTTe
：
h[N]（TfTs）dS；h[N]TfdS[h[N][N]{T}dS
如果在单元边上存在热交换，各条边上的边界换热条件在单元刚度矩阵中生成的附加项为，
hlj
2
1
0
[K]e
6
1
2
0
（6-31）
0
0
0
0
0
0
[K]e
hljm
0
2
1
（6-32）
6
0
1
2
2
0
1
[K]e
hlmi
0
0
0
（6-33）
6
1
0
2
由边界换热条件生成的温度载荷向量为,
1
hTflj
{P}e
2
1
（6-34）
0
0
{P}e
hTfljm
2
1
（6-35）
1
1
e
hTflmi
{P}
2
0
（6-36）
1
6.4温度场分析举例
正方形截面的烟囱如图6-2所示，烟囱由混凝土建造，边长为60cm，通道的边长为20cm,
混凝土的导热系数为k1.4W/（mK）。
假定烟囱内表面的温度为100°C,烟囱外表面暴露
在空气中，空气的温度为30C,换热系数为h20W/（m2K）。
计算烟囱截面内的稳态温
度场。
（参见，FiniteElementMethodTheoryandApplicationwithANSYS,p279）
图6-5热流量分布
稳态温度场分布与物体的初始状态无关，那么是否与材料的导热系数相关？
我们把烟囱
的模型做些修改，假定烟囱壁由两层材料构成。
内层材料为混凝土，外表面的截面尺寸为30cm30cm，烟囱通道的尺寸不变，仍为20cm20cm。
外层材料的导热系数为k0.1W/（mK），外部表面的截面尺寸不变，内部表面的截面尺寸为30cm30cm。
换
热边界条件不变，双层烟囱的有限元模型如图6-6所示。
图6-8双层烟囱的热流量分布
对比两种不同结构烟囱的温度分布和热流量分布，稳态温度场分布与材料的换热系数相
关。
双层烟囱外层材料的导热系数比较小，接近保温材料，热量很快传进内层烟囱，但向外部环境传热慢了很多，所以内层烟囱的温度很高。
比较热流量分布，保温材料的却能够有效的阻止热量散失。
北方城市的供暖管道都包有一层隔热材料就是基于这个道理。
升级会员 