第1章有限单元法基本程式 对于平面应变问题，只要把上述各式中的E换成E/1-),把y换成/(1-) 即可。 1.3.3单元刚度矩阵 从平衡结构中取出的单元也是平衡的。显然，在与其相邻单元的公共边界面 上作用有分布力。此外，单元上就不再有其他力了（因为在有限单元法中，作用于 单元上的全部荷载都要按照等效的原则转化为作用于节点上的等效节点荷载)。 对于平衡的单元来讲，上述边界上的分布力显然是平衡的外力系，并与单元 的内部应力相对应。为方便起见，用等效的单元节点力来代替这些分布力。图 1.3所示为T,单元边界上的分布力和节点力。这样，单元节点力就与单元内部 应力相对应。节点力可理解为节点对单元的作用力。显然，对于单元而言，节点 力是外力；而对于整体结构来说，节点力则是内力。 1--U 图1.3单元节点力 2 根据虚位移原理，可建立单元节点力与节点位移之间的关系即单元刚度方 程，从而得单元刚度矩阵。虚位移就是任意的、微小的可能位移，虚位移原理可 表述如下：如果在虚位移发生之前物体处于平衡状态，那么在虚位移发生时，外 力所做虚功等于物体的虚应变能。现假定单元®发生虚位移，根据式(1.9)有 ou Noa (a) 其中，6a为单元节点虚位移向量。按式(1.13)，单元的虚应变为 8=Boa 6) 将虚位移原理应用于该单元，则单元节点力所做虚功等于单元的虚应变能，即 BaTF=BeTado (e) 其中，为单元体积域：F为单元节点力向量。例如，对于T3单元有 Fi= ∫U (i=1,2,3) F3)
对于平面应变问题，只要把上述各式中的Ｅ 换成Ｅ／（１－ν２），把ν换成ν／（１－ν） 即可。 １３３ 单元刚度矩阵 从平衡结构中取出的单元也是平衡的。显然，在与其相邻单元的公共边界面 上作用有分布力。此外，单元上就不再有其他力了（因为在有限单元法中，作用于 单元上的全部荷载都要按照等效的原则转化为作用于节点上的等效节点荷载）。 对于平衡的单元来讲，上述边界上的分布力显然是平衡的外力系，并与单元 的内部应力相对应。为方便起见，用等效的单元节点力来代替这些分布力。图 １３所示为Ｔ３单元边界上的分布力和节点力。这样，单元节点力就与单元内部 应力相对应。节点力可理解为节点对单元的作用力。显然，对于单元而言，节点 力是外力；而对于整体结构来说，节点力则是内力。 图１３ 单元节点力 根据虚位移原理，可建立单元节点力与节点位移之间的关系即单元刚度方 程，从而得单元刚度矩阵。虚位移就是任意的、微小的可能位移，虚位移原理可 表述如下：如果在虚位移发生之前物体处于平衡状态，那么在虚位移发生时，外 力所做虚功等于物体的虚应变能。现假定单元ｅ发生虚位移，根据式（１９）有 δｕ＝Ｎδａｅ （ａ） 其中，δａｅ 为单元节点虚位移向量。按式（１１３），单元的虚应变为 δε＝Ｂδａｅ （ｂ） 将虚位移原理应用于该单元，则单元节点力所做虚功等于单元的虚应变能，即 δａｅＴＦｅ ＝∫Ω ｅδεＴσｄΩ （ｃ） 其中，Ωｅ 为单元体积域；Ｆｅ为单元节点力向量。例如，对于Ｔ３单元有 Ｆｅ＝ Ｆ１ Ｆ２ Ｆ 烅 烄 烆 烍 烌 ３烎 ， Ｆｉ＝ Ｕｉ Ｖ 烅 烄 烆 烍 烌 ｉ烎 （ｉ＝１，２，３） ８ 第１章 有限单元法基本程式
1.3单元分析 将式a)和b)代入式(c),得 BDBadd 由于虚位移是任意的，等式两边与6T相乘的矩阵应当相等，从而得到单元刚 度方程 r=∫aB'dn 或Fe=Ka 1.17) 其中，K称为单元刚度矩阵，其表达式为 K-BDBdo 1.18) 在平面问题中，被积函数与之坐标无关，故有 BDBrdedyds DBdrdy 其中，A~为单元面积域：t为单元厚度。如果单元的材料是均质的，则矩阵D中 的元素是常量：如果单元厚度也是常量，再注意到ddy=A,则有 [k1k12k13 K=BTDBtA=k21 k22 k23 Lk31 k32 k33 对于平面应力问题，子矩阵为 E ko=BIDBA=4-)A i=1,2,3:j=1,2,3) 单元刚度矩阵是6×6阶矩阵。为研究其性质，将单元刚度方程1.17)写成 下述形式 k11 k12k13 k14 k16 k21 k22 k23 k24k25 k26 k32 k33 k34 k35 k西 u2 (d) k41 k42 k43 k44 k45 v2 U: k51 k52k53k54k55 us V; k61 k62 k6k64k65 k66v3
将式（ａ）和（ｂ）代入式（ｃ），得 δａｅＴＦｅ ＝δａｅ ∫Ｔ Ω ｅＢＴσｄΩ ＝δａｅ ∫Ｔ Ω ｅＢＴＤＢａｅｄΩ 由于虚位移是任意的，等式两边与δａｅＴ相乘的矩阵应当相等，从而得到单元刚 度方程 Ｆｅ ＝∫Ω ｅＢＴσｄΩ 或 Ｆｅ ＝ Ｋｅａｅ （１１７） 其中，Ｋｅ 称为单元刚度矩阵，其表达式为 Ｋｅ ＝∫Ω ｅＢＴＤＢｄΩ （１１８） 在平面问题中，被积函数与ｚ坐标无关，故有 Ｋｅ ＝ Ω ｅＢＴＤＢｔｄｘｄｙｄｚ＝Ａ ｅＢＴＤＢｔｄｘｄｙ 其中，Ａｅ 为单元面积域；ｔ为单元厚度。如果单元的材料是均质的，则矩阵Ｄ 中 的元素是常量；如果单元厚度也是常量，再注意到Ａ ｅｄｘｄｙ＝Ａ ，则有 Ｋｅ＝ＢＴＤＢｔＡ＝ ｋ１１ ｋ１２ ｋ１３ ｋ２１ ｋ２２ ｋ２３ ｋ３１ ｋ３２ ｋ 熿 燀 燄 ３３燅 对于平面应力问题，子矩阵为 ｋｉｊ＝ＢＴ ｉＤＢｊｔＡ＝ Ｅｔ ４（１－ν２）Ａ ｂｉｂｊ＋１－ν ２ ｃｉｃｊ νｂｉｃｊ＋１－ν ２ ｃｉｂｊ νｃｉｂｊ＋１－ν ２ ｂｉｃｊ ｃｉｃｊ＋１－ν ２ ｂｉｂ 熿 燀 燄 ｊ燅 （ｉ＝１，２，３； ｊ＝１，２，３） 单元刚度矩阵是６×６阶矩阵。为研究其性质，将单元刚度方程（１１７）写成 下述形式 Ｕ１ Ｖ１ Ｕ２ Ｖ２ Ｕ３ Ｖ 烅 烄 烆 烍 烌 ３烎 ＝ ｋ１１ ｋ１２ ｋ１３ ｋ１４ ｋ１５ ｋ１６ ｋ２１ ｋ２２ ｋ２３ ｋ２４ ｋ２５ ｋ２６ ｋ３１ ｋ３２ ｋ３３ ｋ３４ ｋ３５ ｋ３６ ｋ４１ ｋ４２ ｋ４３ ｋ４４ ｋ４５ ｋ４６ ｋ５１ ｋ５２ ｋ５３ ｋ５４ ｋ５５ ｋ５６ ｋ６１ ｋ６２ ｋ６３ ｋ６４ ｋ６５ ｋ 熿 燀 燄 ６６燅 ｕ１ ｖ１ ｕ２ ｖ２ ｕ３ ｖ 烅 烄 烆 烍 烌 ３烎 （ｄ） １３ 单 元 分 析 ９
10 第1章有限单元法基本程式 1)k,的物理意义 k,即单元节点位移向量中第)个自由度发生单位位移而其他位移分量为零 时，在第i个自由度方向引起的节点力。 2)K的对称性 根据功的互等定理，1=1时产生的节点力V,等于1=1时产生的节点 力U,即k21=k2。更一般地，有k=k。故K是对称矩阵。事实上，从式 1.18)很容易看出K的对称性。 3)K的奇异性 当单元所有节点位移分量均发生单位位移时，由式()得 U1=k1+k12+k13+k14+k15+k16 V3=k61+k62+k63+k64+k65+k6 由于上述位移属于单元刚体移动，故所有节点力分量应均为零。于是，单元刚度 矩阵中任何一行元素的代数和等于零。由对称性可知，任何一列元素的代数和 也等于零。可见，单元刚度矩阵的行列式为零，即是奇异的。 1.3.4等效节点荷载 离散体系的平衡分析是在节点上进行的，因此需要把作用在单元体上的非 节点荷载化成等效节点荷载。集中力作用点通常都被取做节点，因此荷载移置 的任务是将体力和面力按照静力等效的原则化成节点荷载。所谓静力等效就 是，原荷载与移置后的节点荷载在虚位移上的虚功相等。 当单元发生虚位移时，体力所做的虚功应当等于其等效节点荷载所做的虚 功，即 iaPr=ou fd 其中，P为单元的体力等效节点荷载。将=Na°代入上式，并注意到虚节 点位移的任意性，得 P=NTfan 1.19) 在平面问题中，有 Nndzdy (a 当单元发生虚位移时，面力所做的虚功应当等于其等效节点荷载所做的虚 功，即
（１）ｋｉｊ的物理意义 ｋｉｊ即单元节点位移向量中第ｊ个自由度发生单位位移而其他位移分量为零 时，在第ｉ个自由度方向引起的节点力。 （２）Ｋｅ 的对称性 根据功的互等定理，ｕ１＝１时产生的节点力 Ｖ１ 等于ｖ１＝１时产生的节点 力Ｕ１，即ｋ２１＝ｋ１２。更一般地，有ｋｉｊ＝ｋｊｉ。故 Ｋｅ 是对称矩阵。事实上，从式 （１１８）很容易看出 Ｋｅ 的对称性。 （３）Ｋｅ 的奇异性 当单元所有节点位移分量均发生单位位移时，由式（ｄ）得 Ｕ１＝ｋ１１＋ｋ１２＋ｋ１３＋ｋ１４＋ｋ１５＋ｋ１６ … Ｖ３＝ｋ６１＋ｋ６２＋ｋ６３＋ｋ６４＋ｋ６５＋ｋ 烍 烌 ６６烎 由于上述位移属于单元刚体移动，故所有节点力分量应均为零。于是，单元刚度 矩阵中任何一行元素的代数和等于零。由对称性可知，任何一列元素的代数和 也等于零。可见，单元刚度矩阵的行列式为零，即是奇异的。 １３４ 等效节点荷载 离散体系的平衡分析是在节点上进行的，因此需要把作用在单元体上的非 节点荷载化成等效节点荷载。集中力作用点通常都被取做节点，因此荷载移置 的任务是将体力和面力按照静力等效的原则化成节点荷载。所谓静力等效就 是，原荷载与移置后的节点荷载在虚位移上的虚功相等。 当单元发生虚位移时，体力所做的虚功应当等于其等效节点荷载所做的虚 功，即 δａｅＴＰｅ ｆ ＝∫Ω ｅδｕＴ ｆｄΩ 其中，Ｐｅ ｆ为单元的体力等效节点荷载。将δｕ＝Ｎδａｅ 代入上式，并注意到虚节 点位移的任意性，得 Ｐｅ ｆ ＝∫Ω ｅＮＴ ｆｄΩ （１１９） 在平面问题中，有 Ｐｅ ｆ ＝Ａ ｅＮＴ ｆｔｄｘｄｙ （ａ） 当单元发生虚位移时，面力所做的虚功应当等于其等效节点荷载所做的虚 功，即 ０１ 第１章 有限单元法基本程式
1.3单元分析 11 par 其中，P为单元的面力等效节点荷载。将M=Na代入上式得 NTdr 1.20) 在平面问题中，有 ed 6) 对于T3单元，根据图1.4所示的几何意义（参见习题1.6），不难看出 IN:drdy =A/3 (i=1.2.3) (c) Nids Ih (i=1,2,3) (d) 【例题1.1】求常体力的等效节点荷载。 图1.4N的几何意父 解当体力为常数时，式)中的体力向量可提到积分号外。注意到式。)有 片-心0。。g8 X -9001010rx) 3010101y =兽K Y X Y X YT 当体力为重力，即f=0-Pg]r(图1.5a)时，其等效节点荷载为 乃=-寸gA010101r 2 (a)常体力 图1.5外力 (6)均布面力
δａｅＴＰｅ ｐ ＝∫Γ ｅ σ δｕＴ珔ｐｄΓ 其中，Ｐｅ ｐ 为单元的面力等效节点荷载。将δｕ＝Ｎδａｅ 代入上式得 Ｐｅ ｐ ＝∫Γ ｅ σ ＮＴ珔ｐｄΓ （１２０） 在平面问题中，有 Ｐｅ ｐ ＝∫ ｌ ０ ＮＴ珔ｐｔｄｓ （ｂ） 对于Ｔ３单元，根据图１４所示的几何意义（参见习题１６），不难看出 图１４ Ｎｉ的几何意义 Ａ ｅＮｉｄｘｄｙ＝Ａ／３ （ｉ＝１，２，３） （ｃ） ∫ｌ Ｎｉｄｓ＝ｌ／２ （ｉ＝１，２，３） （ｄ） 【例题１１】 求常体力的等效节点荷载。 解 当体力为常数时，式（ａ）中的体力向量可提到积分号外。注意到式（ｃ）有 Ｐｅ ｆ ＝ｔＡ ｅ Ｎ１ ０ Ｎ２ ０ Ｎ３ ０ ［ ］ ０ Ｎ１ ０ Ｎ２ ０ Ｎ３ Ｔ ｄｘｄｙ Ｘ ｛ ｝Ｙ ＝ｔＡ ３ １ ０ １ ０ １ ０ ［ ］ ０ １ ０ １ ０ １ Ｔ Ｘ ｛ ｝Ｙ ＝ｔＡ ３［Ｘ Ｙ Ｘ Ｙ Ｘ Ｙ］Ｔ 当体力为重力，即ｆ＝［０ －ρｇ］Ｔ（图１５ａ）时，其等效节点荷载为 Ｐｅ ｆ＝－１ ３ρｇｔＡ［０ １ ０ １ ０ １］Ｔ 图１５ 外力 １３ 单 元 分 析 １１
12 第1章有限单元法基本程式 【例题1.2】设31边的长度为l,其上作用沿x方向的均布荷载p.=p(图 1.5b),求其等效节点荷载。 解在31边上，N2=0。注意到式(d,p=[p0]T的等效节点荷载为 -西&0&08'a =00010r 2 1.4整体分析 整体分析的基本任务是建立整体刚度方程，形成整体刚度矩阵和整体节点 荷载向量。 1.4.1节点平衡 根据每个节点处力（各相关单元的节点力和节点荷载）的平衡条件，可得到一 组以节点位移分量为未知量的代数方程组。从结构中取出一个节点i,环绕节点i 有若干个单元，节点i承受的节点荷载为X,Y:(图1.6)，记为P,=[X;Y;]严。 单元作用于节点上的力与节点作用于单元上的节点力大小相等方向相反。 取节点i为脱离体，那么该节点在节点荷载和各单元所施加的节点力之间保持 平衡，即 ∑Ue=X,∑Ve=Y 显然，与节点i无关的单元不进入上述求和式。如果把单元e对节点i的作用 力记为F=[UeVe]T,则上式可用矩阵表示为 F=P; 图1.6节点的平衡
【例题１２】 设３１边的长度为ｌ，其上作用沿ｘ方向的均布荷载ｐｘ＝ｐ（图 １５ｂ），求其等效节点荷载。 解 在３１边上，Ｎ２＝０。注意到式（ｄ），珔ｐ＝［ｐ ０］Ｔ 的等效节点荷载为 Ｐｅ ｐ ＝∫ ｌ ０ Ｎ１ ０ Ｎ２ ０ Ｎ３ ０ ［ ］ ０ Ｎ１ ０ Ｎ２ ０ Ｎ３ Ｔ ｐ ｛ ｝０ ｔｄｓ ＝ｐｌｔ ２ ［１ ０ ０ ０ １ ０］Ｔ １４ 整 体 分 析 整体分析的基本任务是建立整体刚度方程，形成整体刚度矩阵和整体节点 荷载向量。 １４１ 节点平衡 根据每个节点处力（各相关单元的节点力和节点荷载）的平衡条件，可得到一 组以节点位移分量为未知量的代数方程组。从结构中取出一个节点ｉ，环绕节点ｉ 有若干个单元，节点ｉ承受的节点荷载为Ｘｉ，Ｙｉ （图１６），记为Ｐｉ＝［Ｘｉ Ｙｉ］Ｔ。 单元作用于节点上的力与节点作用于单元上的节点力大小相等方向相反。 取节点ｉ为脱离体，那么该节点在节点荷载和各单元所施加的节点力之间保持 平衡，即 ∑ｅ Ｕｉｅ ＝Ｘｉ， ∑ｅ Ｖｉｅ ＝Ｙｉ 显然，与节点ｉ无关的单元不进入上述求和式。如果把单元ｅ对节点ｉ的作用 力记为Ｆｉｅ＝［ ］ Ｕｉｅ Ｖｉｅ Ｔ，则上式可用矩阵表示为 ∑ｅ Ｆｉｅ ＝Ｐｉ 图１６ 节点ｉ的平衡 ２１ 第１章 有限单元法基本程式