★★★★★有限元法的讲解

第四章求解导热问题的有限单元法

第节概述

第节泛函变分原理

第节有限单元法

第节概述

粗略地讲：有限元法是获得微分方程近似解的一种方法，是一种适合计算机来求解的数值计算方法。（元素特性方程和总体合成方程的建立可以采用直接法，变分法，加权余数法和能量平衡法等四种方法之一，所以粗略地说有限元法是获得微分方程近似解的一种方法也有道理）

比较严格的定义：有限单元法是求解泛函变分问题的一种近似方法。

那么这两种说法有什么联系，或者说是共同之处呢？

变分和微分是对未知函数的不同描述，同一连续介质问题往往都可以找到微分和变分的等价表达方式。变分和微分几乎是同时发展起来的两个数学分支，其目的是相同的，都是求解未知函数，但是方法上有很大差别。

在已知边界条件的情况下，求微分方程的精确解析虽然已有完整的理论，但是真正能解出的只有极少数的几种简单情况，因为在很多情况下，微分方程并不存在初等函数解析解。（对于各种各样的映射，初等函数的表达能力实在太有限了，初等函数包括：冥函数、指数函数、对数、三角函数，以及它们的四则运算等。）由于寻求微分方程的初等函数解析解有困难，所以我们在前一章讲述了微分方程的近似解法，即差分法。

泛函变分原理虽然也可以用解析法（即积分）求得未知函数，但是因为有很多被积函数根本无法找到初等原函数，也就不能积分，尤其是对于二维和三维问题，解析法更加困难。所以我们也要寻求泛函变分的近似解法。泛函变分的近似解法包括里兹法和有限元法（里兹法是有限元法的前身），这两种方法的原理完全相同，即：构造一个近似的初等函数，用近似的初等函数去逼近未知函数。因为任何未知函数都可以找到它的近似初等函数（如：包含待定系数的多项式或三角函数），所以从根本上克服了解析法（无法找到初等原函数）的局限性—牺牲极小的理论计算精度，却换回了对大量复杂二维和三维工程问题的适用性。

微分方程的近似解法：差分法

泛函变分的近似解法：里兹法，有限元法

第节泛函变分原理

一、泛函的概念（借助讲解）

二、变分的概念

借助普通函数微分的概念，用类比法讲解

三、泛函的极值条件

借助普通函数的极值条件，用类比法讲解

四、里兹法（补充内容，但是很重要）

泛函变分的近似解法

一、泛函的概念

通过教材§

泛函的概念：函数的函数

泛函与普通函数的区别就在于：函数的自变量是数；而泛函的自变量则是函数，泛函的定义域由具有一定条件的一组函数组成。泛函是一个函数集到一个数集的映射；普通函数则是一个数集到另一个数集的映射。

泛函的表达式：J=J(y)=J[y(x)] J=J(T)=J[T(x ，y)]

泛函的一般式：dx y y x F x y J x x ?=2

1)',,()]([

从物理意义上讲，暂时你也可以把泛函理解成熵，自由能等。对于泛函的具体数值我们并不是特别关心，而更关注它何时取得极值，即取什么样的自变量函数，泛函有极值。

二、变分的概念

三、泛函的极值条件（欧拉方程）

泛函dx y y x F x y J x x ?=2

1)',,()]([的极值条件等价于'0y y d

F F dx

-

= 欧拉方程给出了泛函极值条件与微分方程的关系！ 利用欧拉方程解教材§

在变分问题中，使泛函J(y)为极小值的条件，除0J δ=外，还应有20J δ>（二

阶变分）

四、里兹法（泛函变分的近似解法）

（变分原理在求解微分方程中的应用）连续介质问题经常有着不同的但是等价的表达公式－－微分表达公式和变分表达公式，从例§

求解变分精确解的过程中需要进行各种积分运算，而许多情况下被积函数根本无法找到与相应的初等函数形式的原函数，这说明通过求原函数来计算积分有它的局限性，甚至于可以说这种形式的变分运算根本无法体现出它的运算较微分解方程有什么优越性。 变分法的优越性体现在：我们可以找到一种适用于求得以变分形式表达问题的近似解的简便方法，这种方法叫里兹法，是有限元法的前身 。 例：用里兹法解微分方程：01''=++y y

边界条件0,0

1,0

x y x y ==??==?

解：构造泛函1

22011

[()][(')]22J y x y y y dx =--?，在满足边界条件0,01,0x y x y ==??==?

情况下，

该泛函的极值条件与微分方程01''=++y y 同解。

利用欧拉方程可以很容易的证明两者同解（1--=y F y ，''y F y =） 令近似函数（或称为试探函数）

其中123,,a a a 为待定系数，因为近似函数必须满足边界条件，所以我们构造了这样一个函数。

将近似函数代入泛函，则此时泛函J [()]y x 已经“变质”了，它不再是函数()y x 的

函数（泛函）。而实际上是关于未知数123,,,......n a a a a 的多元函数J(123,,,......n a a a a )

普通多元函数取极值的条件：

0k

J

a ?=? k=1,2,3,…n 从这n 个代数方程中显然可以求得123,,,......n a a a a 等n 个未知数 这里，我们令 21()y a x x =-

＝1

22342110[(144)(2)()]a x x a x x x x x dx -+--+--?

＝1

234111110

[(41)(31)2]a a x a x a x a x dx -++++-?

＝111111111(41)(31)2325a a a a a -

++++-

＝113

610

a -+＝0 得：159a =，所以：)(9

5

2x x y -=

检验：

近似解：1()0.1389

2

y=

由解析法得精确解：

1cos1.0

cos sin1

sin1.0

y x x

-

=+-

述评：

i）一般而言近似函数（试探函数）的项数愈多，达到的精度愈高。

ii）这个方法只是从一族假定中给我们最好的解，因此，非常明显，近似解的精确度和试探函数的选择有关。

iii）我们要求试探函数应定义在整个求解区域上，而且它们至少要满足一些边界条件，通常是要满足全部边界条件。

iv）通常试探函数是由次数连续增大的多项式构成，但在有些情况下，采用其它类型的函数可能是有好处的。

第节有限单元法

一、里兹法的不足（有限元法与里兹法的异同点）

相同点：

两者实质是相同，数学基础都是泛函变分，求解方法都是以初等函数（多项式）去近似未知函数，利用普通多元函数的极值条件来求解多项式中的待定系数。

不同点（里兹法的缺点）：

ⅰ）里兹法将近似函数定义在整个定义域上，而构造近似函数时，又要求它满足所有的边界条件，所以说定义域的边界只能是简单的多边形或多面体（它只能应用在形状相当简单的求解区域上），而不能是复杂形状的边界，故适用范围有限；

有限元素法将近似函数定义在一个单元内（给出未知函数的分片近似函数），虽然对单元同样存在几何上的限制，但是由于形状简单的单元可以被集合起来表示非常复杂的几何形状，因此，有限单元法比里兹法的用途要广泛得多。

有限单元法的近似函数的自变量（即待定系数）是节点上的场变量，应变量当然是单元内区域上各点的场变量。

ⅱ）对里兹法，一般来说，试探函数的项数越多精度越高；对有限元法，一般可选用线性（一次）多项式试探函数，通过缩小单元格尺寸提高计算精度。ⅲ）一般情况下，里兹法要求解出试探函数表达式；而有限元法则一般只需求得离散节点上的未知函数值。

二、有限元法的解题思路（即一般步骤）

1．有限元方法的一般步骤（摘录自其它教材）：

（1）连续介质的离散化

把连续介质划分成很多元素（单元）。即使在一个求解区域中亦可应用不同形状的元素。

（2）指定节点选择插值函数

场变量可能是一个标量，或是一个更高阶的张量。通常选择多项式作为场变量的插值函数，因为它们易于积分和微分，多项式次数的选择取决于每个元素上指定的节点数目，每个节点上未知数的性质和数目以及加在节点

上和元素边界上的某些连续性的要求。节点上场变量的大小以及它的导数的大小很可能是未知的。 （3） 求出元素特性

应用直接法，变分法，加权余数法和能量平衡法等四种方法之一，建立表示各个元素特性的矩阵方程。

（4） 集合元素特性以求得系统方程组（总体合成）

将表示元素性态的矩阵方程组加以合并，形成表示整个求解区域或系统性态的矩阵方程组

系统的矩阵方程组和一个单独元素的方程组具有相同的形式，只是系统方程组包括更多的项，因为它们包括所有的节点。 总体集合成的基础是：对于每个共有某个互连节点的元素而言，在这个节点上的场变量值是相同的

在准备求解系统方程组之前，必须引入边界条件对系统方程组加以修正 （5） 求解系统方程组

求解系统方程组，可求得场变量的未知节点值。 （6） 按照需求要进行附加计算。

利用系统方程组的解来计算其它的重要参数。 有限元素法的应用范围：

① 平衡问题（或不依赖于时间的稳态问题）。

② 定常状态的特征值问题，如固体或流体振动的固有频率。

③ 传播问题（或非稳态问题）。将时间量纲加入前两类问题而产生的新问题。 2．有限元法的解题思路（我的总结）★★★

① 构造泛函J D [T D （x,y ）]= {单元体泛函J e [T e (x,y)]与整个定义域泛函形式相同},使泛函的极值条件与需要求解的微分方程组等价（如傅立叶导热微分方程和边界条件构成的微分方程组）。

② 网格剖分：如果将求解区域D （以二维平面为例）划分为E 个单元（如三角单元）和n 个节点。（注意单元和节点编号规则）

③ 单元分析：在每个单元内部应用温度近似函数123e e e e T a a x a y =++（123,,e e e a a a 是

待定系数），通过解方程可以将待定系数转换为节点温度T i e ,T j e ,T m e （待定系

数的转换）。这样处理以后，温度近似函数写作：T e =N i e T i e +N j e T j e +N m e T m e , 一般称之为温度插值函数，其中N i e ,N j e ,N m e 称为形函数。将温度插值函数

T e =N i e T i e +N j e T j e +N m e T m e 代入单元体泛函J e [T e (x,y)]，单元体泛函J e [T e (x,y)]实际上变质为普通多元函数J e (T e i ,T j e ,T m e )。{即：由关于未知函数T e (x,y)的泛函J e [T e (x,y)]转化成关于待定系数的普通多元函数J e (T e i ,T j e ,T m e )} 然后求解普通多元函数J e

(T

e

i ,T j e

,T m e

)对节点温度T i e ,T j e

,T m e 的偏导数

e

i e T J ??,e j e T J ??,e m

e

T J ??。（因为里兹法只有一个单元，而有限元素法包括很多单元，所以它们均不等于零，只有总体合成后才等于零。）

④ 对全部E 个单元都进行类似③的运算

⑤ 总体合成：为使泛函J e [T e (x,y)]取得极值，要求 其中k=1,2,3…,n

上式包含n 个代数方程，解之可得n 个节点的温度。

三、温度场泛函（与傅里叶导热微分方程等价的泛函表达式）

一般来说，泛函可根据物理学原理（如熵的最大化，能量最低原量等）或者由微分方程推导而来，具体的推导过程已超出了课程的范围。

对于一个二维稳态温度场，函数族(,)T x y 所构成的泛函: 式中D 表示温度场所对应的积分区域。

对于一个二维非稳态温度场，函数族(,,)T x y t 所构成的泛函，理论上: 但实际上，目前非稳态温度场的处理并不采用上式，而是先把t （或T

t

??）当作常数（即先把非稳态温度场作为稳态温度场处理），然后将

T

t

??用一价差商代替。 泛函J 的形式与原微分方程及其边界条件密切相关，原微分方程和边界条件决定泛函J 的形式。（原微分方程决定了泛函在区域内部的积分式，边界条件决定了泛函在区域边界上的积分式）

1．一般形式傅立叶方程所描述的温度场的泛函 微分方程的一般形式

其中：'q 表示单位体积内热源的发热强度；

k x ，k y ，k z 表示不同坐标方向上的导热系数，且导热系数随温度不同亦即随坐标位置（x,y,z ）不同而变化。

边界条件的一般形式可表示为

其中：h 为综合换热系数； T ∞是环境温度；q 为单位边界表面上的热流强

度；cos ,cos ,cos αβγ 分别为边界曲面S 外法线的方向余弦。

与上述二式对应的泛函

2．无内热源三维稳态温度场,但x y z k k k ≠≠时，泛函

由于q ’=0,0T

t ?=?,故

3．无内热源三维稳态温度场的泛函，且x y z k =k =k k = (并且是一常数) 1) 若同时具备第三类边界条件，即（对流和辐射换热边界条件）

()0T h T T n k ∞?+-=?，

T

n

??为边界外法线方向上的温度梯度 其中:,,h k T ∞是常数

2）若具备第二类边界条件（边界上热流已知）

0T k q n

?+=?,则此时()0h T T ∞-= 3）若具备第一类边界条件（表面温度已知）

(),,s

T

f x y z = 则()0h T T ∞-=,0q =

因为在这种边界条件下,边界上温度是固定的，边界上的积分为常数，其变分为零，因而在泛函J 中没有与s T 有关的项。 4．无内热源二维稳态温度场，且x y k =k =k (并且是一常数) 1）若具备第三类边界条件：

2222[(,)]D T T h h J T x y dxdy T T T ds x y k k ∞Γ??????????

=++-?? ? ? ?????????????

????（改） 教材p117写作

两者差一常数2

k

,完全等价。

2）若具备第二类边界条件： 3）若具备第一类边界条件：

四、网格剖分

网络剖分亦称为时空离散化，就是将时间和空间分割成若干有限的小单元。 （一）时间离散化

对于稳态导热问题，只涉及空间的剖分。而对于非稳态导热问题，目前对这类方程的泛函变分问题尚未很好解决，通常的处理方法是在空间域内用有限单元网格划分，而在时间域内则用有限差分网格划分（即：有限元法处理空间变量，有限差分法处理时间变量）。具体而言，有两种处理方式（两种处理方法得到的结果完全相同）：

(1) 令时间变量暂时固定，即先考虑在某一瞬间对泛函变分，然后再考虑t 的变化，把T

t ??用差分展开为1

n n T T t

--?（亦可采用其它差分格式）。 (2) 先把T t ??用差分展开为1

n n T T

t

--?（亦可采用其它差分格式），然后进

行泛函变分运算，而在变分运算过程中，把1n T -与t ?均作为常量处理（n T T =，即n T 作为变量）。

（二）空间离散化

有限元法（Finite Element Method ，简称FEM ），它的最大特点是单元形状和疏密程度可以任意变化，因而对具有复杂几何形状和条件的物体极为适用。 1．网格剖分的规则

从理论上说，有限元法可以采用多种形状的多边形或多面体网格单元，但是通常以三角形和四边形单元用得最多。平面温度场计算一般取三角形网络，三角形的形状和疏密程度可任意变化，对复杂形状的物体尤其能显示其优点。 网格剖分的规则：（结合课本P115图讲解规则的应用）

（1）一个单元中只能包含一种材料。

（2）不能把一个三角形的顶点取在另一个三角形的边上。

（3）把求解区域划分为内部单元和边界单元，规定边界单元只能有一条边落在边界 上（注意：对曲线边界，实际是两个节点落在边界上）。求解区域的边

界为曲线时，剖分时用直线代替，并取为三角形单元的一边。

（4）三角形的三条边的长度不宜相过大，图为计算精度受单元最长边长与最短边长之比值的控制。

2．单元和节点编号规则

为了减少失误，提高计算效率，便于计算机实现。单元和节点编号有一定规则：

（1）单元编号规则（结合课本P115图讲解规则的应用）

单元分为两大类：内部单元，边界单元

单元编号规则：

先编内部单元后编边界单元，一、二、三类边界单元依次编排。（原因

是单元体泛函由依序由简单到复杂，有助于简化总体合成后矩阵方程中

的系数矩阵，从而提高算法效率。其原因仍然有待进一步深入思考）

记作：①②③④⑤…

（2）节点编号规则

每一个节点分别有两个序号：局部序号，全局序号

节点全局序号编号规则：（结合课本P136图讲解规则的应用）

依序编排，力求使邻近节点的编号尽可能接近，尤其是同一单元中三个

节点的编号相差不宜太大。（原因：合理地编排节点全局序号，可以使

得总体合成后矩阵方程中的系数矩阵有规则地分布在主对角线附近相

对狭小的宽度内，构成所谓带状矩阵。这种对称正定带状矩阵对利用消

去法求解极为有利，有助于大大提高算法效率。）

记作：1，2，3，4…

节点局部序号编号规则：

(Ⅰ)内部单元：i，j，m按逆时针方向依次编排，起始位置可以是任意

的；

(Ⅱ)边界单元：对于三角形单元，节点i与边界相对（即内节点），i，

j，m逆时针方向依次编排；对于其他多边形单元，保证

j，m落在边界上。

记作：i，j，m，…

五、单元分析

1. 温度插值函数

温度插值函数实际上就是里兹法试探函数（近似函数）的变形。

插值函数一般由不同次幂的多项式组成，但也有采用三角函数或指数函数组成的乘积表示，但最常用的多项式插值函数。有限元插值函数分为两大类，一类只要求插值多项式本身在插值点取已知值，称为拉格朗日(Lagrange)多项式插值；另一种不仅要求插值多项式本身，还要求它的导数值在插值点取已知值，称为哈密特(Hermite)多项式插值。单元坐标有笛卡尔直角坐标系和无因次自然坐标，有对称和不对称等。常采用的无因次坐标是一种局部坐标系，它的定义取决于单元的几何形状，一维看作长度比，二维看作面积比，三维看作体积比。在二维有限元中，三角形单元应用的最早，近来四边形等参元的应用也越来越广。对于二维三角形和四边形电源单元，常采用的插值函数为有Lagrange插值直角坐

标系中的线性插值函数及二阶或更高阶插值函数、面积坐标系中的线性插值函数、二阶或更高阶插值函数等。 单元内部的温度插值函数

设三角形单元e 上的温度e T 是x ，y 的线性函数，即

其中：123,,e e e a a a 为待定系数

将节点坐标及温度代入得

写成矩阵形式：

利用矩阵求逆的方法，可以解出待定系数123,,e e e a a a

式中：?为三角形单元的面积，其值为

（为了使面积不为负值，故要求节点的局部序号i,j,m 必须按逆时针方向排列）

将<2>式代入<1>式，得

1()()()2e e e e e e e e e e e e

i i i i j j j j m m m m a b x c y T a b x c y T a b x c y T ??=++++++++???

<3> 上式通常可以简写成

式中 []{}1()21()

2,1()2e e e e

i i i i e e e e e i j j j j e e j e e e e e

m m m m m e e e e

i j m N a b x c y T N a b x c y T T T N a b x c y N N N N ?=++???

???=++???

?=?

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

,,e e e

i j m N N N 称为形状因子，形状函数或简称形函数。

单元边界上的温度插值函数

〈4〉式适用于整个单元区域，对于单元边界当然也是适用的，但是，根据线性插值的概念，既然单元边界jm 上两节点的温度分别为e j T 和e m T ，那么直线jm 上任意一点的温度e T 应在e j T 和e m T 之间呈线性变化，而与e i T 无关，这样在边界上就可以构造一个更加简单的插值函数

式中g 为一参变量，[0,1]g ∈；g=0对应于节点j ，g=1对应于节点m 。 边界弧长（积分变量）s 与参变量g 的关系

显然jm 的边长e i S 是一个已知数

而曲线积分中的边界弧长变量s 与e i S 间的关系可利用g 联系起来 2. 单元变分运算（即普通多元函数偏导数计算） 二维稳态温度场单元变分计算

（实际上可理解为偏导计算，因为用温度插值函数代泛函以后，泛函实际上已经成为一个普通多元函数）

2221[(,)]()22e e e e e e b jm e k T T J T x y dxdy h T hT T ds x y ∞??????????=++-?? ? ? ?????????????

??? <1> （我不成熟的理解：泛函表达式的第二部分原来也是对整个单元的面积分，但是被积函数只有在边界上才有非零值，所以转换成了边界的线积分）

泛函表达式的第二部分,对整个求解区域而言，边界一般是封闭的，记作“??”；然而对某个单元而言，边界一般是不封闭的，记作“jm

?

”。

将三角形单元内部的插值公式以及边界上的插值公式代入泛函，泛函

(),e e

b J T x y ????实际上变成了普通多元函数

求e b e i J T ??，e b e j J T ??，e

b

e m

J T ??可以采用两种方法：

(1) 先将插值公式代入<1>式，然后再求e b e i J T ??，e b e j J T ??，e

b

e m

J T ??

Q 2221[(,)]()22e e e e e e b jm e k T T J T x y dxdy h T hT T ds x y ∞??????????=++-?? ? ? ?????????????

??? <1> 单元内部插值函数： 单元边界插值函数： 积分变量的转换：

接下去，重点看下列推导过程的后半部分。 同理：

写成矩阵形式

[]{}{}e b

e i e

e

e

e e ii ij im

i i e e e e e

e e e e b ji jj jm j j e j e e e e e mi mj

mm

m m e b e m J T k k k T p J k k k T p K T P T k k k T p J T ?????

??????????????????

??=-=-????

?????????????????

????????????

<2> 其中

(2) 先推导变分公式，再将插值函数代入

令e T u x

?=?，则()()e e i i f u f u u

T u T ???=???g 同理可求得e b e j J T ??，e

b

e

m

J T ??（后续推导略）

内部单元变分计算

内部单元的泛函

22

[(,)]2e e e e

i e k T T J T x y dxdy x y ????????=+?? ? ???????????

?? <3> 将三角形单元内部的插值公式以及边界上的插值公式代入泛函，泛函

(),e e

i J T x y ????实际上变成了普通多元函数

求e i e i J T ??，e i e j J T ??，e

i e m

J T ??可以采用两种方法：

(1) 先将插值公式代入<3>式，然后再求e i e i J T ??，e i e j J T ??，e

i e m

J T ??

Q 22[(,)]2e e e e i e k T T J T x y dxdy x y ??????

??=+?? ? ????????????? <3> 单元内部插值函数： 同理 同理：

写成矩阵形式 其中

从形式上说，内部单元的变分结果与边界单元基本相同；从本质上说，仅需

令边界单元变分结果中的0e i S =，即可得内部单元的变分结果。 (2) 先推导变分公式，再将插值函数代入

令e T u x

?=?，则()()e

e i i

f u f u u

T u T ???=???g 同理可求得e b e j J T ??，e

b

e

m

J T ??（后续推导略） 二维非稳态温度场单元变分计算，无内热源

二维非稳态温度场泛函变分问题尚未很好解决，目前采用的处理方法有两种：

①先用有限元法处理空间域，后用有限差分法处理时间域 ②先用有限差分法处理时间域，后用有限元法处理空间域

我们采用前一种处理方法。

二维非稳态温度场泛函（先固定时间变量，即先把T

t

??看作常数）：

无内热源二维非稳态温度场第三类边界单元泛函 无内热源二维稳态温度场第三类边界单元泛函

两者相比较，仅有e e

T c T t

ρ??项是新的，并且我们在整个变分运算（实际上是普通函数积分运算）过程中，把e

T t

??作为常数处理。

因此

可以证明 得到：

同理可得： 写成矩阵形式

注意：此处的

[]e N 并非插值函数中的形状函数

其中

由上可见，非稳态温度场的变分结果与稳态温度场基本相同，仅增加了一项

[]e

e

T N t ????????

而已。

内部单元变分计算 （略）

六、总体合成

总体合成的原理

（结合教材p129图讲解）

三角形单元的温度插值函数为

将三角形单元的温度插值函数代入单元体泛函，则泛函(),e e

J T x y ????实际上变质

为普通多元函数

同理，若整个求解域包含n 个节点，则定义在整个求解域上的温度插值函数可以写作

将定义在整个求解域上的温度插值函数代入关于整个求解域的泛函，则泛函

(),D J T x y ????实际上变质为普通多元函数

所以泛函(),D J T x y ????（或者说普通多元函数()12,,,D n J T T T L ）取极值的条件为

0D

l

J T ?=?，1,2,,l n =L ；即节点的全局序号 由于定义整个求解域的泛函是各单元体泛函之和；而且，定义在整个求解域上的温度插值函数，其待定系数,1,2,,l T l n =L 是各单元体的顶点（节点），所以

10D e

E

e l l

J J T T =??==??∑，1,2,,l n =L ；即节点的全局序号

注意：

①因为我们只要求定义整个求解域的泛函(),D J T x y ????取极值，并没有要求单元

体泛函也取极值，所以0e l

J T ????→?不要求；

②0e

e

l e

l

l

J T J T J T ?????????→????≠???不要求=0,节点l 不在单元e 上0,节点l 落在单元e 上 所以总体合成的主要任务就是找出节点l 究竟落在哪些单元体e 上，节点的全局序号l 如何与各相关单元体e 的局部编号i,j,m 一一对应。 无内热源二维稳态温度场的总体合成

现以图中的节点3为例来说明总体合成过程 总体合成结果

无内热源二维稳态温度场三角形边界单元变分结果 类推，无内热源二维稳态温度场总体合成矩阵方程

系数矩阵（包括温度刚度矩阵及其他系数矩阵）的建立方法 规律：

编号小的单元，影响总体刚度矩阵中左上角部分；

任意节点i 和节点j 若不在一个单元（i ，j 表示节点的全局序号），则总体刚度矩阵中的元素k ij =0；

与节点i 相邻的节点数量越多，则总体刚度矩阵i 行中非零元素个数越多。 相邻节点的全局序号越接近越好（实际上严格来说，并不特别要求构成一个单元的所有节点的全局序号很接近，只不过两种要求产生的效果是相同的）。 无内热源二维非稳态温度场的总体合成

无内热源二维非稳态温度场三角形边界单元变分结果 类推，无内热源二维非稳态温度场总体合成矩阵方程 记作：

对于某一特定时刻，上式可写成 如果边界条件恒定，则

将t

T t ???

?????用后向差商展开 例如：111t t t

t T T T t t

-??-=

??

无内热源二维非稳态温度场总体合成矩阵方程，后向差商展开形式

七、加权余量法