有限元分析基礎(chǔ)平面問題的有限元

第三章平面問題的有限元§3.1三角形單元§3.2矩形單元§3.3平面等參元§3.4平面單元應(yīng)用比較分析第三章平面問題的有限元3.1三角形單元§3.1.1位移函數(shù)§3.1.2應(yīng)變與應(yīng)力矩陣§3.1.3形函數(shù)的性質(zhì)§3.1.4單元剛度矩陣及其特點(diǎn)§3.1.5等效節(jié)點(diǎn)載荷§3.1.6整體剛度矩陣及其特點(diǎn)§3.1.7邊界條件§3.1.8算例§3.1.9收斂準(zhǔn)則§3.1.10從能量原理推導(dǎo)剛度矩陣本節(jié)主要研究三節(jié)點(diǎn)三角形單元，首先建立以單元節(jié)點(diǎn)位移表示單元內(nèi)各點(diǎn)位移的關(guān)系式。設(shè)單元e的節(jié)點(diǎn)編號為i、j、m，如圖所示。由彈性力學(xué)平面問題可知，每個節(jié)點(diǎn)在其單元平面內(nèi)的位移有兩個分量，所以整個三角形單元將有六個節(jié)點(diǎn)位移分量，即六個自由度。用列陣可表示為式中的子向量為

（i、j、m輪換）式中：——節(jié)點(diǎn)i在x軸方向的位移分量；

——節(jié)點(diǎn)i在y軸方向的位移分量。3.1.1位移函數(shù)在有限單元法中，雖然是用離散化模型來代替原來的連續(xù)體，但每一個單元體仍是一個彈性體，所以在其內(nèi)部依然是符合彈性力學(xué)基本假設(shè)的，彈性力學(xué)的基本方程在每個單元內(nèi)部同樣適用。從彈性力學(xué)平面問題的解析解法中可知，如果彈性體內(nèi)的位移分量函數(shù)已知，則應(yīng)變分量和應(yīng)力分量也就確定了。但是，如果只知道彈性體中某幾個點(diǎn)的位移分量的值，那么就不能直接求得應(yīng)變分量和應(yīng)力分量。因此，在進(jìn)行有限元分析時，必須先假定一個位移模式。由于在彈性體內(nèi)，各點(diǎn)的位移變化情況非常復(fù)雜，很難在整個彈性體內(nèi)選取一個恰當(dāng)?shù)奈灰坪瘮?shù)來表示位移的復(fù)雜變化，但是如果將整個區(qū)域分割成許多小單元，那么在每個單元的局部范圍內(nèi)就可以采用比較簡單的函數(shù)來近似地表示單元的真實(shí)位移，將各單元的位移模式連接起來，便可近似地表示整個區(qū)域的真實(shí)位移函數(shù)。這種化繁為簡、聯(lián)合局部逼近整體的思想，正是有限單元法的絕妙之處。3.1.1位移函數(shù)基于上述思想，可以選擇一個單元位移模式，單元內(nèi)各點(diǎn)的位移可按此位移模式由單元節(jié)點(diǎn)位移通過插值而獲得。假設(shè)三節(jié)點(diǎn)三角形單元內(nèi)任意點(diǎn)的位移是x、y的線性函數(shù)，如下是待定常數(shù)，因三角形單元共有六個自由度，且位移函數(shù)在三個節(jié)點(diǎn)處的數(shù)值應(yīng)該等于這些點(diǎn)處的位移分量的數(shù)值。假設(shè)節(jié)點(diǎn)i、j、m的坐標(biāo)分別為（式中～、）、（）、（），代入位移函數(shù)，得3.1.1位移函數(shù)從上式中可以解得式中：——三角形i、j、m的面積。為保證求得的面積為正值，節(jié)點(diǎn)i、j、m的編排次序必須是逆時針方向，如圖所示。的計算公式，如下3.1.1位移函數(shù)將式所求得得系數(shù)代入位移函數(shù)，得式中：（i、j、m輪換）3.1.1位移函數(shù)同理可以求出單元內(nèi)任意點(diǎn)的y方向位移，如下令

（i、j、m輪換）則單元內(nèi)任意點(diǎn)的位移可以表示成如下形式式中：——單元內(nèi)位移矢量3.1.1位移函數(shù)位移函數(shù)可以進(jìn)一步寫成如下形式式中：

（i、j、m輪換）是坐標(biāo)的函數(shù)，它們反映了單元的位移狀態(tài)，所以一般稱之為形狀函數(shù)，簡稱形函數(shù)。顯然，形函數(shù)決定了單元內(nèi)的“位移模式”，反映了i節(jié)點(diǎn)位移對單元內(nèi)任意點(diǎn)位移的貢獻(xiàn)。3.1.1位移函數(shù)有了單元的位移模式，就可以利用平面問題的幾何方程求出單元內(nèi)任意點(diǎn)的應(yīng)變，式中：——單元應(yīng)變列陣；

——單元應(yīng)變轉(zhuǎn)換矩陣，或稱單元幾何矩陣。3.1.2應(yīng)變與應(yīng)力矩陣單元應(yīng)變轉(zhuǎn)換矩陣，為3×6矩陣，其子矩陣為（i、j、m輪換）可以看出，表示i節(jié)點(diǎn)位移對單元應(yīng)變的貢獻(xiàn)率。當(dāng)單元確定后，也就確定了，此時單元內(nèi)的應(yīng)變僅依賴于節(jié)點(diǎn)的位移。對于三節(jié)點(diǎn)三角形單元，由于面積和

等都是常量，所以單元應(yīng)變轉(zhuǎn)換矩陣中的諸元素都是常量，因而單元中各點(diǎn)的應(yīng)變分量也都是常量，即單元內(nèi)各點(diǎn)的應(yīng)變相同，通常稱這種單元為常應(yīng)變單元。、3.1.2應(yīng)變與應(yīng)力矩陣求得應(yīng)變之后，再將應(yīng)變代入物理方程，即公式（2.19），便可推導(dǎo)出以節(jié)點(diǎn)位移表示的應(yīng)力，如下式中：——單元應(yīng)力矩陣。令則式中：——單元應(yīng)力轉(zhuǎn)換矩陣。3.1.2應(yīng)變與應(yīng)力矩陣的子矩陣為單元應(yīng)力轉(zhuǎn)換矩陣可以看出，表示i節(jié)點(diǎn)位移對單元應(yīng)力的貢獻(xiàn)率。當(dāng)單元確定后，也就確定了，此時單元內(nèi)的應(yīng)力僅依賴于節(jié)點(diǎn)的位移。對于三節(jié)點(diǎn)三角形單元，由于面積和

等都是常量，所以單元應(yīng)力轉(zhuǎn)換矩陣中的諸元素都是常量，因而單元中各點(diǎn)的應(yīng)力分量也都是常量，即單元內(nèi)各點(diǎn)的應(yīng)力相同，通常稱這種單元為常應(yīng)力單元。、3.1.2應(yīng)變與應(yīng)力矩陣可見，對于三節(jié)點(diǎn)三角形單元，由于所選取的位移模式是線性的，因而導(dǎo)致單元內(nèi)各點(diǎn)的應(yīng)變、應(yīng)力相同，但相鄰單元將具有不同的應(yīng)力和應(yīng)變，即在單元的公共邊界上應(yīng)力和應(yīng)變的值將會有突變，但位移卻是連續(xù)的。3.1.2應(yīng)變與應(yīng)力矩陣由前文可知，形函數(shù)決定了單元內(nèi)的“位移模式”，反映了i節(jié)點(diǎn)位移對單元內(nèi)任意點(diǎn)位移的貢獻(xiàn)，形函數(shù)具有如下性質(zhì)：1.形函數(shù)在各單元節(jié)點(diǎn)上的值，具有“本點(diǎn)是1、它點(diǎn)為零”的性質(zhì)即

（i、j、m輪換）證明如下：3.1.3形函數(shù)的性質(zhì)在節(jié)點(diǎn)j、m上類似有在節(jié)點(diǎn)i上3.1.3形函數(shù)的性質(zhì)2.在單元內(nèi)任意點(diǎn)上，三個形函數(shù)之和等于1即證明：將單元內(nèi)任意點(diǎn)的坐標(biāo)代入上式，得，因此容易得到此結(jié)論。由此可見，三個形函數(shù)中只有二個是獨(dú)立的?？梢宰C明，，3.1.3形函數(shù)的性質(zhì)例如：假設(shè)單元各節(jié)點(diǎn)的位移均相同并等于，則根據(jù)此性質(zhì)可得單元內(nèi)任意點(diǎn)的位移為，顯然，該性質(zhì)反映了單元剛體位移，即當(dāng)單元做剛體運(yùn)動時，單元內(nèi)任意點(diǎn)的位移均等于剛體位移。3.1.3形函數(shù)的性質(zhì)3.三角形單元任意一條邊上的形函數(shù)，僅與該邊的兩端節(jié)點(diǎn)坐標(biāo)有關(guān)。例如在ij邊上，有事實(shí)上，因ij

邊的直線方程方程為3.1.3形函數(shù)的性質(zhì)代入及，得則3.1.3形函數(shù)的性質(zhì)另外，由可以求得利用形函數(shù)的這一性質(zhì)可以證明，相鄰單元的位移分別進(jìn)行線性插值之后，在其公共邊上將是連續(xù)的。例如，對圖3.2所示的單元具有ij公共邊，可知，在ij邊上有這樣，不論按哪個單元來計算，公共邊ij上的位移均由下式表示3.1.3形函數(shù)的性質(zhì)由此可見，在公共邊上的位移u、v將完全由公共邊上的兩個節(jié)點(diǎn)i、j的位移所確定，因而相鄰單元的位移是連續(xù)的。例3.1求如圖3.3所示等邊三角形的形函數(shù)。三角形的面積將以上數(shù)據(jù)代入式（3.10），得3.1.3形函數(shù)的性質(zhì)①將各節(jié)點(diǎn)坐標(biāo)代入以上三式，有即滿足“本點(diǎn)是1、它點(diǎn)為零”的性質(zhì)。3.1.3形函數(shù)的性質(zhì)即滿足“在單元內(nèi)任意點(diǎn)上，三個形函數(shù)之和等于1”的性質(zhì)。②將三個形函數(shù)相加，得，即在ij邊上任意位置的位移，將完全由節(jié)點(diǎn)i、j的位移所確定。取其它邊同樣可以得出其結(jié)論，即滿足“三角形單元任意一條邊上的形函數(shù)，僅與該邊的兩端節(jié)點(diǎn)坐標(biāo)有關(guān)”的性質(zhì)。③以ij邊為例，在ij邊上，則3.1.3形函數(shù)的性質(zhì)單元在節(jié)點(diǎn)處受力，單元會發(fā)生變形，也就是說單元在節(jié)點(diǎn)處受到的力與單元節(jié)點(diǎn)位移之間有必然的聯(lián)系。單元間正是通過與節(jié)點(diǎn)間的相互作用力連接起來成為整體，而一個單元仍是一個彈性體，如果將每一個單元的受力－位移關(guān)系找到，則整體的受力－位移關(guān)系也容易清楚。為了推導(dǎo)單元的節(jié)點(diǎn)力和節(jié)點(diǎn)位移之間的關(guān)系，可應(yīng)用虛位移原理對單元進(jìn)行分析。單元在節(jié)點(diǎn)處受到的力稱為單元節(jié)點(diǎn)力，單元在節(jié)點(diǎn)力的作用下處于平衡，節(jié)點(diǎn)力可以表示為

式中：——節(jié)點(diǎn)i沿x方向的節(jié)點(diǎn)力分量；

——節(jié)點(diǎn)i沿y方向的節(jié)點(diǎn)力分量。3.1.4單元剛度矩陣及其特點(diǎn)單元節(jié)點(diǎn)的虛位移為相應(yīng)的單元內(nèi)的虛位移場為

單元內(nèi)的虛應(yīng)變?yōu)?/p>

設(shè)單元僅在節(jié)點(diǎn)上受力，則單元節(jié)點(diǎn)力在節(jié)點(diǎn)虛位移上的虛功為

3.1.4單元剛度矩陣及其特點(diǎn)單元內(nèi)的應(yīng)力在虛應(yīng)變上所做的功為式中：——在體積上積分。令單元剛度矩陣為由虛功原理可知，式（3.27）和（3.28）相等，再由虛位移的任意性，可得3.1.4單元剛度矩陣及其特點(diǎn)上式為單元剛度矩陣的定義式，其推導(dǎo)過程與形函數(shù)的具體表達(dá)形式、節(jié)點(diǎn)個數(shù)均無關(guān)，這說明該表達(dá)式具有普遍意義，則上式是單元節(jié)點(diǎn)力與單元節(jié)點(diǎn)位移之間的關(guān)系式，即單元節(jié)點(diǎn)平衡方程，或稱單元剛度方程，具有普遍意義。對于三節(jié)點(diǎn)三角形單元，由于、中的元素都是常量，這里假設(shè)單元厚度均勻，則式中：——單元面積；

——單元厚度。3.1.4單元剛度矩陣及其特點(diǎn)上式進(jìn)一步可以表示為將幾何矩陣式（3.16）和彈性矩陣式（2.19）代入，可得每個子塊的矩陣表達(dá)式，如下3.1.4單元剛度矩陣及其特點(diǎn)單元剛度矩陣具有如下特點(diǎn)：1.中的每一個元素都是一個剛度系數(shù)假設(shè)單元的節(jié)點(diǎn)位移為，由，得到節(jié)點(diǎn)力為

3.1.4單元剛度矩陣及其特點(diǎn)由此可得：表示i節(jié)點(diǎn)在水平方向產(chǎn)生單位位移時，在節(jié)點(diǎn)i的水平方向上需要施加的節(jié)點(diǎn)力；表示i節(jié)點(diǎn)在水平方向產(chǎn)生單位位移時，在節(jié)點(diǎn)i的垂直方向上需要施加的節(jié)點(diǎn)力。選擇不同的單元節(jié)點(diǎn)位移，可以得到單元剛度矩陣中每個元素的物理含義，如下：表示s節(jié)點(diǎn)在水平方向產(chǎn)生單位位移時，在節(jié)點(diǎn)r的水平方向上需要施加的節(jié)點(diǎn)。表示s節(jié)點(diǎn)在水平方向產(chǎn)生單位位移時，在節(jié)點(diǎn)r的垂直方向上需要施加的節(jié)點(diǎn)力；表示s節(jié)點(diǎn)在垂直方向產(chǎn)生單位位移時，在節(jié)點(diǎn)r的水平方向上需要施加的節(jié)點(diǎn)力；表示s節(jié)點(diǎn)在水平方向產(chǎn)生單位位移時，在節(jié)點(diǎn)r的垂直方向上需要施加的節(jié)點(diǎn)力；因此，中的每一個元素都是一個剛度系數(shù)，表示單位節(jié)點(diǎn)位移分量所引起的節(jié)點(diǎn)力分量，中的每一個子塊表示單元某個節(jié)點(diǎn)位移矢量對單元某個節(jié)點(diǎn)力矢量的貢獻(xiàn)率。3.1.4單元剛度矩陣及其特點(diǎn)2.是對稱矩陣由式（3.32）可知，是對稱矩陣。根據(jù)對稱性，可以減少計算量和存儲量。

3.是奇異矩陣，即根據(jù)這一性質(zhì)，當(dāng)已知單元節(jié)點(diǎn)位移時，可以從式（3.31）中求出單元節(jié)點(diǎn)力。反之，由于單元剛陣奇異，不存在逆矩陣，因此，當(dāng)已知單元節(jié)點(diǎn)力時不能求出單元上的節(jié)點(diǎn)位移。從物體變形的實(shí)際情況來說，單元剛陣奇異是必須的，因?yàn)閱卧水a(chǎn)生變形外，還會產(chǎn)生剛體位移，即僅僅依靠節(jié)點(diǎn)力是無法唯一地確定剛體位移的。事實(shí)上，當(dāng)單元的節(jié)點(diǎn)力為零時，單元仍可做剛體運(yùn)動。3.1.4單元剛度矩陣及其特點(diǎn)例如，假定單元產(chǎn)生了x方向的單位剛體位即并假設(shè)此時對應(yīng)的單元節(jié)點(diǎn)力為零，則由

得可以得到，在單元剛度矩陣中1，3，5列中對應(yīng)行的系數(shù)相加為零，由行列式的性質(zhì)可知，。3.1.4單元剛度矩陣及其特點(diǎn)

4.單元的剛度不隨單元或坐標(biāo)軸的平行移動而改變單元的剛度取決于單元的形狀、大小、方向和彈性系數(shù)，而與單元的位置無關(guān)，即不隨單元或坐標(biāo)軸的平行移動而改變。因此，只要單元的形狀、大小、方向和彈性系數(shù)相同，無論單元出現(xiàn)在任何位置均有相同的單元剛度矩陣，根據(jù)對稱性，可以減少計算量。，厚度為t。例3.2求圖3.4所示等腰直角三角形單元的剛度矩陣，設(shè)解根據(jù)1（0，a）、2（0，0）、3（a，0）的坐標(biāo)和式（3.8）可得3.1.4單元剛度矩陣及其特點(diǎn)單元面積，將以上各數(shù)據(jù)代入幾何矩陣式（3.16）為由于，則平面應(yīng)力與平面應(yīng)變問題的彈性矩陣相同，為3.1.4單元剛度矩陣及其特點(diǎn)由式（3.18）可得應(yīng)力轉(zhuǎn)換矩陣，如下則單元剛度矩陣為3.1.4單元剛度矩陣及其特點(diǎn)從以上求得的結(jié)果可以看出，單元剛度矩陣是對稱的。1、3、5列相加等于零，2、4、6列相加也等于零，單元剛度矩陣是奇異的。另外在幾何矩陣中所涉及到的元素均為“坐標(biāo)差量”，而與坐標(biāo)的具體值無關(guān)，從而也說明了單元剛度矩陣不隨單元或坐標(biāo)軸的平行移動而改變。1）獲取單元節(jié)點(diǎn)信息；2）計算單元面積及參數(shù)、等；3）計算幾何矩陣中各元素；4）根據(jù)問題類型（平面應(yīng)力、平面應(yīng)變）計算彈性矩陣；5）計算應(yīng)力轉(zhuǎn)換矩陣；6）計算單元剛度矩陣。從單元剛度矩陣的計算公式可以看出，在求解單元剛度矩陣時，可以不必求得應(yīng)力轉(zhuǎn)換矩陣，但在后期求解單元應(yīng)力時要用到該矩陣，所以在此處一并求出，也可以根據(jù)實(shí)際情況具體考慮。3.1.4單元剛度矩陣及其特點(diǎn)實(shí)際上作用在物體上的外力可以直接作用在節(jié)點(diǎn)上，也可以不作用在節(jié)點(diǎn)上，然而在有限元法中要求作用在物體上的各種外力必須用作用在節(jié)點(diǎn)上的力表示。這一點(diǎn)體現(xiàn)了有限法中“離散化”這一概念，即將連續(xù)體離散成只有在節(jié)點(diǎn)處相連的單元體，單元與單元之間的聯(lián)系只能通過節(jié)點(diǎn)。如前文所述，單元內(nèi)任意點(diǎn)的位移、應(yīng)力、應(yīng)變等變量最終都用單元節(jié)點(diǎn)位移來表示。同樣，作用在物體上的各種外力也必須用作用在節(jié)點(diǎn)上的力表示，這一過程稱為外力的靜力等效移置，所得到的節(jié)點(diǎn)力稱為等效節(jié)點(diǎn)力。

3.1.5等效節(jié)點(diǎn)載荷靜力等效原則：

對于剛體來說，所謂靜力等效原則就是單元上原有的外力系和將外力系向各節(jié)點(diǎn)移置所得的等效節(jié)點(diǎn)力，二者向同一點(diǎn)簡化應(yīng)具有相同的主失和主矩；對于彈性體來說，所謂靜力等效原則就是指單元上的外力系和將該力系向各節(jié)點(diǎn)移置后的等效節(jié)點(diǎn)力，二者在虛位移上的虛功相等，也即外力作用在單元上所引起的變形能和移置后的等效節(jié)點(diǎn)力在單元上引起的變形能相等，在一定的位移模式下這種移置是唯一的。3.1.5等效節(jié)點(diǎn)載荷1.彈性體靜力等效原則——虛功原理虛功等效的思路是：就一個單元來說，把作用在單元上的外力移置到節(jié)點(diǎn)上后，應(yīng)當(dāng)與原來的實(shí)際外力所作虛功等效。其計算方法是：對任意允許的微小虛位移，原外力所作虛功等于移置后的等效節(jié)點(diǎn)力所作虛功。即式中：——節(jié)點(diǎn)虛位移；——單元虛位移場；——集中力作用處產(chǎn)生的虛位移。3.1.5等效節(jié)點(diǎn)載荷等式左邊表示等效節(jié)點(diǎn)力在節(jié)點(diǎn)虛位移上所作虛功；右邊第一項(xiàng)表示作用在單元任意點(diǎn)c處的集中力所作虛功、第二項(xiàng)表示表面力所作的虛功、第三項(xiàng)表示體積力所作的虛功。將代入上式，并考慮的任意性，得式中：——由集中力等效的單元等效節(jié)點(diǎn)力；——由表面力等效的單元等效節(jié)點(diǎn)力；——由體積力等效的單元等效節(jié)點(diǎn)力；——形函數(shù)在集中力作用處的值。3.1.5等效節(jié)點(diǎn)載荷單元等效節(jié)點(diǎn)力、和的表達(dá)式分別為

以上三式不但適用三角形單元，同樣適用其它單元，具有普遍意義。下面根據(jù)以上三式，分別進(jìn)行討論。3.1.5等效節(jié)點(diǎn)載荷可以看出，當(dāng)外力作用在節(jié)點(diǎn)上時，如作用在i節(jié)點(diǎn)上，此時，，則計算所得到的等效節(jié)點(diǎn)力與原外力在數(shù)值上是相同的。3.1.5等效節(jié)點(diǎn)載荷（1）作用在單元任意點(diǎn)的集中力3.1.5等效節(jié)點(diǎn)載荷（2）作用在單元上的表面力現(xiàn)舉例說明，如圖所示的單元e，在ij邊上作用有表面力，假設(shè)ij邊的長度為l，單元厚度t為常數(shù)，其上任一點(diǎn)P距節(jié)點(diǎn)i的距離為s。有，，，則例3.3如圖所示，設(shè)31邊的長度為l，單元厚度為t，其上作用沿x方向的均布荷載，求其等效節(jié)點(diǎn)荷載。解：在31邊上，，

，則等效節(jié)點(diǎn)荷載為從計算結(jié)果來看，等效后的節(jié)點(diǎn)力是將原來作用在31邊上的均布面力平分到1、3節(jié)點(diǎn)上。3.1.5等效節(jié)點(diǎn)載荷（3）作用在單元上的體積力設(shè)單元所受體積力為，單元厚度t為常數(shù)，如所示，則3.1.5等效節(jié)點(diǎn)載荷當(dāng)單元體是均質(zhì)、等厚、比重為時，則、，有即當(dāng)單元有均勻自重時，在三節(jié)點(diǎn)三角形單元中每個節(jié)點(diǎn)沿y方向的等效節(jié)點(diǎn)力等于單元總重量的1/3。3.1.5等效節(jié)點(diǎn)載荷2.剛體靜力等效原則——力系簡化（1）作用在單元任意點(diǎn)的集中力若集中力作用在節(jié)點(diǎn)上，則可直接將集中力認(rèn)為是等效節(jié)點(diǎn)力即可。若集中力作用于單元的某一邊界上，如圖所示，則先將分解為、，然后按線段的比例把分別移置到i，j兩點(diǎn)上。即3.1.5等效節(jié)點(diǎn)載荷（2）作用在單元上的表面力如圖所示，已知在ij邊受有均布面力，單元厚度t為常數(shù)，則移置到i、j節(jié)點(diǎn)上的等效節(jié)點(diǎn)力為如圖所示，當(dāng)某一邊上有三角形分布的面力時，則等效節(jié)點(diǎn)力為

3.1.5等效節(jié)點(diǎn)載荷（3）作用在單元上的體積力如圖所示，如果三角形單元ijm的重心c上受有體積力，單元厚度t為常數(shù)，則按剛體靜力等效原理可把體積力直接移置到i，j，m三個節(jié)點(diǎn)上，即根據(jù)圣維南原理可知，按剛體靜力等效原則移置后的等效節(jié)點(diǎn)力所引起的應(yīng)力誤差只是局部的不會影響到整體，即物體內(nèi)部距實(shí)際外力作用處相當(dāng)遠(yuǎn)的各點(diǎn)處應(yīng)力，與其外力的具體分布情況關(guān)系很小。在線性位移模式下，彈性體和剛體按各自的靜力等效原則移置外力的結(jié)果是一致的。3.1.5等效節(jié)點(diǎn)載荷

整體節(jié)點(diǎn)載荷列陣：由各節(jié)點(diǎn)所受等效節(jié)點(diǎn)力按節(jié)點(diǎn)號碼以從小到大的順序排列組成的列陣。等效節(jié)點(diǎn)力是由集中力、表面力和體積力共同移置構(gòu)成的，其中集中力包括直接作用在彈性體上的外力和邊界約束力，如支座反力。前文對單元體進(jìn)行了分折，得到了單元剛度方程，但要解決問題，還必須進(jìn)一步建立整個計算模型的整體剛度方程。完成這一步的關(guān)鍵，在于怎樣將單元的剛度矩陣和節(jié)點(diǎn)荷載列陣，分別“組裝”成整體剛度矩陣和整體節(jié)點(diǎn)荷載列陣。這里通過研究任意節(jié)點(diǎn)的平衡來建立整體剛度矩陣，該方法不但比較直觀、易懂，而且對怎樣編寫計算機(jī)程序是很有幫助的。為了研究整體剛度矩陣的組裝過程，先引入兩個概念。

整體節(jié)點(diǎn)位移列陣：由各節(jié)點(diǎn)位移按節(jié)點(diǎn)號碼以從小到大的順序排列組成的列陣。3.1.6整體剛度矩陣及其特點(diǎn)式中：，

3.1.6整體剛度矩陣及其特點(diǎn)不失一般性，僅考慮計算模型中有4個單元，如圖所示。四個單元的整體節(jié)點(diǎn)位移列陣為對每個單元都可以寫出相應(yīng)的單元剛度方程，即單元節(jié)點(diǎn)平衡方程。例如，對①號單元，有式中：——①號單元中第i（i=1,2,3）節(jié)點(diǎn)所受力。為了便于組裝整體剛度矩陣，將上式以整體節(jié)點(diǎn)位移表示，即

3.1.6整體剛度矩陣及其特點(diǎn)——①號單元的擴(kuò)大剛度矩陣或稱為單元貢獻(xiàn)矩陣。同理，對于②單元，有式中：——②號單元中第i（i=1,3,4）節(jié)點(diǎn)所受力；——②號單元的擴(kuò)大剛度矩陣。3.1.6整體剛度矩陣及其特點(diǎn)對于④單元，有式中：——④號單元中第i（i=3,4,5）節(jié)點(diǎn)所受力；——④號單元的擴(kuò)大剛度矩陣。3.1.6整體剛度矩陣及其特點(diǎn)對于任意一個節(jié)點(diǎn)，可能承受兩種力的作用，一種是其它單元給予該節(jié)點(diǎn)的反作用力；另一種是作用在節(jié)點(diǎn)上的等效節(jié)點(diǎn)力。對整體而言，前者屬于內(nèi)力，后者屬于外力，每個節(jié)點(diǎn)在兩種力的作用下處于平衡。將各單元剛度方程左邊相加，即將各節(jié)點(diǎn)所受力相加，由于對于整體而言，單元給予節(jié)點(diǎn)的反作用力屬于內(nèi)力，在相加過程中相互抵消，所以各節(jié)點(diǎn)所受力相加的結(jié)果只有外力，即等效節(jié)點(diǎn)力，從而得到整體節(jié)點(diǎn)荷載列陣，如下3.1.6整體剛度矩陣及其特點(diǎn)將各單元剛度方程右邊相加，從而得到整體剛度矩陣，如下3.1.6整體剛度矩陣及其特點(diǎn)通過以上分析得，整體節(jié)點(diǎn)載荷與整體節(jié)點(diǎn)位移之間的關(guān)系式，即結(jié)構(gòu)整體有限元方程，如下式中：——整體剛度矩陣。3.1.6整體剛度矩陣及其特點(diǎn)整體剛度矩陣組裝的基本步驟：1）將單元剛度矩陣中的每個子塊放到在整體剛度矩陣中的對應(yīng)位置上，得到單元的擴(kuò)大剛度矩陣。注意對于單元剛度矩陣是按照局部編碼排列的，即對應(yīng)單元剛度矩陣中的i、j、m；對于整體剛度矩陣是按照整體編碼排列的，即按節(jié)點(diǎn)號碼以從小到大的順序排列。在組裝過程中，必須知道單元節(jié)點(diǎn)的局部編碼與該節(jié)點(diǎn)在整體結(jié)構(gòu)中的整體編碼之間的關(guān)系，才能得到單元剛度矩陣中的每個子塊在整體剛度矩陣中的位置。將單元剛度矩陣中的每個子塊按總體編碼順序重新排列后，可以得到單元的擴(kuò)大矩陣。例如在圖中，單元②的局部編碼為i、j、m，對應(yīng)整體編碼為1、3、4，然后將單元②剛度矩陣中的每個子塊按總體編碼順序重新排列后，可以得到單元的擴(kuò)大矩陣。注意有些書籍中將局部編碼表示為1、2、3或1，2，3等；2）將全部單元的擴(kuò)大矩陣相加得到整體剛度矩陣。3.1.6整體剛度矩陣及其特點(diǎn)通過以上組裝過程可以得到組裝整體剛度矩陣的一般規(guī)則：1）結(jié)構(gòu)中的等效節(jié)點(diǎn)力是相關(guān)單元結(jié)點(diǎn)力的疊加，整體剛度矩陣的子矩陣是相關(guān)單元的單元剛度矩陣子矩陣的集成；2）當(dāng)整體剛度矩陣中的子矩陣中r=s時，該節(jié)點(diǎn)（節(jié)點(diǎn)r或s）被哪幾個單元所共有，則就是這幾個單元的剛度矩陣中的子矩陣的相加。如應(yīng)該是單元①－④中對應(yīng)子矩陣的集成，即3.1.6整體剛度矩陣及其特點(diǎn)3.1.6整體剛度矩陣及其特點(diǎn)3）當(dāng)中時，若rs邊是組合體的內(nèi)邊，則就是共用該邊的兩相鄰單元剛度矩陣中的子矩陣的相加。如13邊為單元①和②的共用邊，則4）當(dāng)中r和s不同屬于任何單元時，則=0。如節(jié)點(diǎn)r＝1和s＝5不同屬于任何單元，此時=0。上述組裝基本步驟和規(guī)則具有普遍意義，對于不同類型、不同形式的單元，只是相應(yīng)節(jié)點(diǎn)的子矩陣的階數(shù)（節(jié)點(diǎn)自由度×節(jié)點(diǎn)自由度）可能不同，至于組裝整體剛度矩陣的規(guī)律仍是相同的。正是應(yīng)為有了這種組裝規(guī)律，使得有限元法能夠很方便地應(yīng)用電子計算機(jī)進(jìn)行計算。例3.4如圖所示有限元模型，彈性模量為，厚度為，為簡化計算取，求整體剛度矩陣。3.1.6整體剛度矩陣及其特點(diǎn)單元編號①②③④整體編碼1、2、32、4、55、3、23、5、6局部編碼i、j、mi、j、mi、j、mi、j、m以整體編碼表示的單元剛度矩陣子塊解：該模型中共有6個節(jié)點(diǎn)，4個單元，各單元的信息如表3.1所示。表3.1各單元信息3.1.6整體剛度矩陣及其特點(diǎn)同例3.2類似的分析，得根據(jù)單元剛度矩陣的性質(zhì)，可知，由于幾何矩陣，所以有，整體剛度矩陣中的各子塊是對所有單元相應(yīng)的子塊求和得到的（實(shí)際只是對相關(guān)單元求和），其中各子塊矩陣均為2行×2列，整體剛度矩陣用子塊矩陣可以表示為3.1.6整體剛度矩陣及其特點(diǎn)上式中任意一子塊矩陣均為2行×2列，在計算過程中，無需將每個單元剛度矩陣進(jìn)行擴(kuò)大，只需判斷整體剛度矩陣子塊的下標(biāo)，然后利用組裝整體剛度矩陣的一般規(guī)則進(jìn)行計算，如，由圖形可知，節(jié)點(diǎn)2由單元①、②和③所共有，則3.1.6整體剛度矩陣及其特點(diǎn)，由圖形可知，25邊為單元②和③的共用邊，則，由圖形可知，節(jié)點(diǎn)1、5不同屬于任何單元，則采用同樣的方法進(jìn)行計算，得到整體剛度矩陣為3.1.6整體剛度矩陣及其特點(diǎn)3.1.6整體剛度矩陣及其特點(diǎn)4.是奇異矩陣，在排除剛體位移后，它是正定陣3.是稀疏矩陣，非零元素呈帶狀分布用有限元方法分析復(fù)雜工程問題時，節(jié)點(diǎn)的數(shù)目比較多，整體剛度矩陣的階數(shù)通常也是很高的。那么在進(jìn)行計算時，如果存儲整體剛度矩陣的全部元素，將會浪費(fèi)較大的資源、降低計算效率。如果根據(jù)整體剛度矩陣的特點(diǎn)進(jìn)行編寫程序，可以大大節(jié)省資源、并提高計算效率。因此有必要了解和掌握整體剛度矩陣的特點(diǎn)，整體剛度矩陣具有以下幾個顯著的特點(diǎn)：1.是對稱矩陣2.中主對角元素總是正的3.1.6整體剛度矩陣及其特點(diǎn)3.1.6整體剛度矩陣及其特點(diǎn)1.是對稱矩陣由單元剛度矩陣的對稱性和整體剛度矩陣的組裝過程，可知整體剛度矩陣必為對稱矩陣。利用對稱性，在計算機(jī)編寫程序時，只存儲整體剛度矩陣上三角或下三角部分即可。2.中主對角元素總是正的例如，剛度矩陣中的元素表示節(jié)點(diǎn)2在x方向產(chǎn)生單位位移，而其它位移均為零時，在節(jié)點(diǎn)2的x方向上必須施加的力；表示節(jié)點(diǎn)2在y方向產(chǎn)生單位位移，而其它位移均為零時，在節(jié)點(diǎn)2的y方向上必須施加的力。很顯然在此情況下力的方向應(yīng)該與位移方向一致，故應(yīng)為正號。3.1.6整體剛度矩陣及其特點(diǎn)3.是稀疏矩陣，非零元素呈帶狀分布如果遵守一定的節(jié)點(diǎn)編號規(guī)則，就可使矩陣的非零元素都集中在主對角線附近呈帶狀。整體剛度矩陣中的子矩陣只有當(dāng)下標(biāo)s等于r或者s與r同屬于一個單元時才不為零，這就說明，在第r雙行中非零子矩陣的塊數(shù)，應(yīng)該等于節(jié)點(diǎn)r周圍直接相鄰的節(jié)點(diǎn)數(shù)目加1。可見，中元素一般都不是填滿的，是稀疏矩陣，且非零元素呈帶狀分布。以下圖所示的單元網(wǎng)格為例，其整體剛度矩陣中的非零子塊（每個子塊為2行2列）的分布情況如下圖所示，圖中陰影部分表示該子塊不為零，其它子塊部位均為零。3.1.6整體剛度矩陣及其特點(diǎn)3.1.6整體剛度矩陣及其特點(diǎn)顯然，帶狀剛度矩陣的帶寬取決于單元網(wǎng)格中相鄰節(jié)點(diǎn)號碼的最大差值D。把半個斜帶形區(qū)域中各行所具有的非零元素的最大個數(shù)叫做剛度矩陣的半帶寬（包括主對角元），用B表示，如下

B=2(D+1)通常的有限元程序，一般都利用剛度矩陣的對稱和稀疏帶狀的特點(diǎn)，在計算求解中，只存儲上半帶的元素，即所謂的半帶存儲。因此，在劃分完有限元網(wǎng)格進(jìn)行節(jié)點(diǎn)編號時，要采用合理的編碼方式，使同一單元中相鄰兩節(jié)點(diǎn)的號碼差盡可能小，以便節(jié)省存儲空間、提高計算效率。3.1.6整體剛度矩陣及其特點(diǎn)對于同樣的有限元單元網(wǎng)格，按照圖（a）的結(jié)點(diǎn)編碼，最大的半帶寬為B=2×(10－4+1)=14；按照圖（b）的結(jié)點(diǎn)編碼，最大的半帶寬為B=2×(10－2+1)＝18；按照圖（c）的結(jié)點(diǎn)編碼，最大的半帶寬為B=2×(10－6+1)＝10。(a)(b)(c)3.1.6整體剛度矩陣及其特點(diǎn)4.是奇異矩陣，在排除剛體位移后，它是正定陣無約束的彈性體（或結(jié)構(gòu)物）的整體剛度矩陣是奇異的，不存在逆矩陣，即關(guān)于位移的解不唯一。這是因?yàn)閺椥泽w在外力的作用下處于平衡，外力的分量應(yīng)該滿足三個靜力平衡方程。這反映在整體剛度矩陣中就意味著存在三個線性相關(guān)的行或列，所以是個奇異陣，不存在逆矩陣。例如：設(shè)彈性體在外力的作用下處于平衡，這時相應(yīng)的解為，然后在給予彈性體以剛體位移而相應(yīng)的節(jié)點(diǎn)位移，這時，仍是問題的解，因?yàn)閯傮w位移不會破壞平衡。注：當(dāng)排除剛體位移后，整體剛度矩陣是正定矩陣。前文已經(jīng)提到在排除剛體位移后，整體剛度矩陣是正定的，方程才可求得唯一解。排除剛體位移可以通過引入邊界約束條件來實(shí)現(xiàn)，這里介紹兩種比較簡單的引入已知位移的方法。1．代入法2．乘大數(shù)法3.1.7邊界條件3.1.7邊界條件1．代入法該方法保持方程組仍為2n×2n系統(tǒng)，僅對整體剛度矩陣和整體載荷列陣進(jìn)行修正。下面以一個只有四個方程的簡單例子加以說明，方程如下假定系統(tǒng)中節(jié)點(diǎn)位移、，則當(dāng)引入這些節(jié)點(diǎn)的已知位移之后，方程就變成3.1.7邊界條件若，則然后，用這組維數(shù)不變的方程來求解所有的節(jié)點(diǎn)位移。顯然，其解答仍為原方程的解答。在手算時，可直接將零位移約束所對應(yīng)的整體剛度矩陣中的行和列直接劃去，使得整體剛陣的維數(shù)變小，更便于手算。3.1.7邊界條件2．乘大數(shù)法將中與指定的節(jié)點(diǎn)位移相對應(yīng)的主對角元素乘上一個大數(shù)，同時將中的對應(yīng)元素?fù)Q成指定的節(jié)點(diǎn)位移值、該大數(shù)與節(jié)點(diǎn)位移相對應(yīng)的主對角元素三者的乘積。若把此方法用于上面的例子，則方程就變成該方程組的第一個方程為解得，這種方法就是使中相應(yīng)行的修正項(xiàng)遠(yuǎn)大于非修正項(xiàng)。3.1.7邊界條件1．代入法2．乘大數(shù)法在以上的兩種方法中，代入法接近人工解法，雖然該方法比較直觀，但該方法對剛度矩陣改變較多，程序效率不高。乘大數(shù)法對剛度矩陣改變較少，工作量較小，但相乘的“大數(shù)”若取得過大，求解時會發(fā)生“溢出”、若取得太小則會引起較大的誤差。根據(jù)前面的討論，現(xiàn)以三角形常應(yīng)變單元為例來說明應(yīng)用有限元法求解彈性力學(xué)平面問題的具體步驟。1）力學(xué)模型的確定。根據(jù)工程實(shí)際情況確定問題的力學(xué)模型，并按一定比例繪制結(jié)構(gòu)圖、注明尺寸、載荷和約束情況等；2）結(jié)構(gòu)進(jìn)行離散化。將彈性體劃分為許多單元，并對節(jié)點(diǎn)進(jìn)行編號，確定全部節(jié)點(diǎn)的坐標(biāo)值；對單元進(jìn)行編號，并列出各單元三個節(jié)點(diǎn)的節(jié)點(diǎn)號；3）計算等效節(jié)點(diǎn)載荷，形成整體載荷列陣；4）計算單元剛度矩陣；5）組裝整體剛度矩陣；6）處理約束，引入邊界條件；7）求解線性方程組，得到節(jié)點(diǎn)位移和節(jié)點(diǎn)力；8）計算單元應(yīng)力；9）整理計算結(jié)果（后處理部分）。3.1.8算例例3.5如圖所示為一厚度t=1cm的均質(zhì)正方形薄板，邊長為2m，上下受均勻拉力=106N/m2，材料彈性模量為E，泊松比，不記自重，試用有限元法求其應(yīng)力分量。

3.1.8算例解：1）力學(xué)模型的確定。由于此結(jié)構(gòu)長、寬遠(yuǎn)大于厚度，而載荷作用于板平面內(nèi)，且沿板厚均勻分布，故可按平面應(yīng)力問題處理，考慮到結(jié)構(gòu)和載荷的對稱性，可取結(jié)構(gòu)的1/4來研究。3.1.8算例表3.2節(jié)點(diǎn)坐標(biāo)值節(jié)點(diǎn)坐標(biāo)1234x0110y00112）結(jié)構(gòu)進(jìn)行離散化。為了計算簡便并能說明整個計算過程，將該1/4結(jié)構(gòu)離散為兩個三角形單元。節(jié)點(diǎn)編號、單元劃分如圖所示，各節(jié)點(diǎn)的坐標(biāo)值如表3.2所示，各單元編碼如表3.3所示。3.1.8算例表3.3單元編碼單元號①②整體編碼1、2、33、4、1局部編碼i、j、mi、j、m3）計算等效節(jié)點(diǎn)載荷，形成整體載荷列陣。按照彈性體靜力等效原則，應(yīng)用公式（3.41），并參考例3.3可得單元②由均布載荷等效的節(jié)點(diǎn)力為3.1.8算例若按照剛體靜力等效原則，參考公式（3.45），同樣可得到以上結(jié)果?？紤]到支座反力，整體載荷列陣為4）計算單元剛度矩陣。對于單元①有，局部編碼i、j、m對應(yīng)整體編碼1、2、3，利用公式（3.34）并參考例3.2，可得單元①的剛度矩陣為3.1.8算例對于單元②有，局部編碼i、j、m對應(yīng)整體編碼3、4、1，由單元剛度矩陣的性質(zhì)可知，單元②的剛度矩陣與單元①的剛度矩陣相同。

5）組裝整體剛度矩陣。在組裝整體剛度矩陣時，要注意單元局部編碼與整體編碼之間的關(guān)系，按照組裝整體剛度矩陣的一般規(guī)則和基本步驟，可得整體剛度矩陣為3.1.8算例

6）處理約束，引入邊界條件。由以上分析可知，整體有限元方程為3.1.8算例根據(jù)邊界約束條件：，采用代入法引入邊界條件，劃去整體剛度矩陣中1、2、4、7的行和列，得3.1.8算例

7）求解線性方程組，得到節(jié)點(diǎn)位移。求解上面方程組可得出相應(yīng)得節(jié)點(diǎn)位移為3.1.8算例所以，整體節(jié)點(diǎn)為列陣為代入整體有限元方程，便可得整體載荷列陣為若要求得每個單元所受的節(jié)點(diǎn)力，則可以從整體節(jié)點(diǎn)位移列陣中提取單元所對應(yīng)的節(jié)點(diǎn)位移，然后在利用，便可求得單元節(jié)點(diǎn)所受的力。如單元①所對應(yīng)的節(jié)點(diǎn)位移為3.1.8算例求得單元①所對應(yīng)的節(jié)點(diǎn)力為單元②所對應(yīng)的節(jié)點(diǎn)位移為求得單元②所對應(yīng)的節(jié)點(diǎn)力為注意，這里求得的單元節(jié)點(diǎn)位移、和單元所對應(yīng)的節(jié)點(diǎn)力、均對應(yīng)局部編碼。很顯然，每個單元均滿足力的平衡。8）計算單元應(yīng)力。先求出各單元的應(yīng)力矩陣、，然后再求得各單元的應(yīng)力分量，如下3.1.8算例9）整理計算結(jié)果（后處理部分）。對于三節(jié)點(diǎn)三角形單元，單元內(nèi)各點(diǎn)的應(yīng)力值相等，算出的應(yīng)力一般作為單元形心處的應(yīng)力。由于單元應(yīng)力為常數(shù)，整個結(jié)構(gòu)的應(yīng)力場呈階梯狀，在單元之間不連續(xù)。而工程上往往更加關(guān)心邊界和節(jié)點(diǎn)上的受力情況，因此，必須對所得到的應(yīng)力再次進(jìn)行處理，得到更加合理的應(yīng)力場，并得到所需點(diǎn)上的應(yīng)力值。這里介紹兩種簡單的方法，一種方法稱為節(jié)點(diǎn)平均法，即把環(huán)繞某一節(jié)點(diǎn)的各單元的應(yīng)力加以平均作為該節(jié)點(diǎn)的應(yīng)力值。例如圖中節(jié)點(diǎn)3的應(yīng)力為3.1.8算例為了使通過這樣平均得來的應(yīng)力比較接近實(shí)際情況，要求環(huán)繞節(jié)點(diǎn)的各單元尺寸不應(yīng)相差太大。這種做法，對內(nèi)節(jié)點(diǎn)比較好，對邊界點(diǎn)則可能很差。因此，邊界節(jié)點(diǎn)處的應(yīng)力不宜直接由單元應(yīng)力平均來獲得，而要根據(jù)內(nèi)節(jié)點(diǎn)的應(yīng)力構(gòu)造插值函數(shù)推算出來。例如圖中邊界點(diǎn)1的應(yīng)力，可以先用節(jié)點(diǎn)平均法求得節(jié)點(diǎn)2、3、4處的應(yīng)力，在構(gòu)造相應(yīng)的插值函數(shù)推算邊界點(diǎn)1的應(yīng)力，如常用的拋物線插值公式如下3.1.8算例另一種方法稱為單元平均法，即把兩相鄰單元的應(yīng)力加以平均，用以表示公共邊界中點(diǎn)的應(yīng)力。為了由這樣平均所得到的應(yīng)力具有較好的精度，兩相鄰單元的面積不應(yīng)相差太大。如圖中單元②和③邊界的中點(diǎn)處的應(yīng)力為3.1.8算例在不同的有限元軟件中均具有各自的后處理方法，但無論后期怎樣處理，應(yīng)力場來源于應(yīng)力的計算，應(yīng)力的精度主要依賴于單元的尺寸、單元的類型（位移模式）。對于有限元這種數(shù)值計算方法，一般總是希望隨著網(wǎng)格的逐步細(xì)分所得到的解能夠收斂于問題的精確解。根據(jù)前面的分析，可知在有限元分析中，一旦確定了單元的形狀之后，位移模式的選擇將是非常關(guān)鍵的。由于載荷的移置、應(yīng)力矩陣和剛度矩陣的建立等等，都依賴于單元的位移模式，所以，如果所選擇的位移模式與真實(shí)的位移分布有很大的差別，那么就很難獲得良好的數(shù)值解。為了保證解答的收斂性，要求位移模式必須滿足以下三個條件，即1）位移模式必須包含單元的剛體位移。2）位移模式必須包含單元的常應(yīng)變。3）位移模式在單元內(nèi)要連續(xù)、且在相鄰單元之間的位移必須協(xié)調(diào)。3.1.9收斂準(zhǔn)則3.1.9收斂準(zhǔn)則

1）位移模式必須包含單元的剛體位移。也就是說，當(dāng)節(jié)點(diǎn)位移是由某個剛體位移所引起時，彈性體內(nèi)將不會產(chǎn)生應(yīng)變。所以，位移模式不但要具有描述單元本身形變的能力，而且還要具有描述由于其它單元形變而通過節(jié)點(diǎn)位移引起單元剛體位移的能力。2）位移模式必須包含單元的常應(yīng)變。每個單元的應(yīng)變一般都是包含著兩個部分：一部分是與該單元中各點(diǎn)的坐標(biāo)位置有關(guān)的應(yīng)變（即所謂各點(diǎn)的變應(yīng)變）；另一部分是與位置坐標(biāo)無關(guān)的應(yīng)變（即所謂的常應(yīng)變）。從物理意義上看，當(dāng)單元尺寸無限縮小時，每個單元中的應(yīng)變應(yīng)該趨于常量。因此，在位移模式中必須包含有這些常應(yīng)變，否則就不可能使數(shù)值解收斂于正確解。3.1.9收斂準(zhǔn)則3）位移模式在單元內(nèi)要連續(xù)、且在相鄰單元之間的位移必須協(xié)調(diào)。當(dāng)選擇多項(xiàng)式來構(gòu)成位移模式時，單元內(nèi)的連續(xù)性要求總是得到滿足的，單元間的位移協(xié)調(diào)性，就是要求單元之間既不會出現(xiàn)開裂也不會出現(xiàn)重疊的現(xiàn)象。通常，當(dāng)單元交界面上的位移取決于該交界面上節(jié)點(diǎn)的位移時，就可以保證位移的協(xié)調(diào)性。在有限單元法中，把能夠滿足條件1和2的單元，稱為完備單元，完備單元是收斂的必要條件。滿足條件3的單元，叫做協(xié)調(diào)單元或保續(xù)單元。前面討論過的三節(jié)點(diǎn)三角形單元，均能同時滿足上述三個條件，因此屬于完備的協(xié)調(diào)單元，完備的協(xié)調(diào)單元是收斂的充分條件。在某些梁、板及殼體分析中，要使單元滿足條件3比較困難，所以實(shí)踐中有時也出現(xiàn)一些只滿足條件1和2的單元，其收斂性往往也能夠令人滿意。特別是放松條件3的單元，即完備而不協(xié)調(diào)的單元，這時要保證其收斂必須通過分片試驗(yàn)的測試，關(guān)于分片試驗(yàn)的內(nèi)容請參見相關(guān)書籍。目前，完備而不協(xié)調(diào)的單元，已獲得了很多成功應(yīng)用。3.1.9收斂準(zhǔn)則在有限元分析中，將實(shí)際連續(xù)體分成許多單元體的組合后，根據(jù)線性或非線性位移的假定，人為地選擇一個位移場，通過這些措施所得到的模型比實(shí)際連續(xù)體的剛性要高。因而，近似模型的剛度是實(shí)際連續(xù)體剛度的上界。若選擇不協(xié)調(diào)單元，那么這種模型可能由于單元分離、疊加或單元之間形成鉸而降低剛性，這種影響就有可能使得不協(xié)調(diào)元比應(yīng)用協(xié)調(diào)元所得的結(jié)果要好。不過，應(yīng)用不協(xié)調(diào)單元事先不能肯定所得的剛度是真實(shí)剛度的上界。換句話說，不協(xié)調(diào)元不一定象協(xié)調(diào)單元那樣剛硬，可能比較柔軟，因此有可能會比協(xié)調(diào)單元收斂得快。3.1.9收斂準(zhǔn)則經(jīng)驗(yàn)證明，根據(jù)巴斯卡（Pascal）三角形來選擇二維多項(xiàng)式的各項(xiàng)。在二維多項(xiàng)式中，如果包含有對稱軸一邊的某一項(xiàng)，那么就必須同時包含有另一邊的對稱項(xiàng)。選擇多項(xiàng)式位移模式時，還應(yīng)該要考慮的一個因素是，多項(xiàng)式中的項(xiàng)數(shù)必須等于或稍大于單元邊界上的外節(jié)點(diǎn)的自由度數(shù)。通常是取項(xiàng)數(shù)與單元的外節(jié)點(diǎn)的自由度數(shù)相等，取過多的項(xiàng)數(shù)是不恰當(dāng)?shù)?。從?jié)點(diǎn)平衡得到結(jié)構(gòu)整體有限元方程式的方法物理概念清晰，但這種方法在數(shù)學(xué)上不夠嚴(yán)謹(jǐn)。本節(jié)從勢能原理導(dǎo)出有限元方程，該方法是一種更加嚴(yán)謹(jǐn)、適用性更加廣泛的方法。由第2章可知，系統(tǒng)的總勢能為如果有集中力作用等式右邊還應(yīng)加一項(xiàng)集中力所作的功，如下3.1.10從能量原理推導(dǎo)剛度矩陣當(dāng)把連續(xù)體離散為有限個單元，并在各單元內(nèi)假設(shè)了位移模式，在單元內(nèi)應(yīng)力與應(yīng)變用節(jié)點(diǎn)位移（或單元中）來描述，、，在連續(xù)體上的整體積分等于在每個單元上分別積分再求和，如下由勢能原理可知，對于處于平衡狀態(tài)的彈性體，位移真實(shí)解應(yīng)總使總勢能取極小值，反之依然，則由，得3.1.10從能量原理推導(dǎo)剛度矩陣上式進(jìn)一步可寫為如果有集中力作用，有上式進(jìn)一步可寫為上式同式（3.55），即為整體結(jié)構(gòu)有限元方程。

3.1.10從能量原理推導(dǎo)剛度矩陣3.2矩形單元雖然三角形單元具有很好的“適應(yīng)性”，幾乎任何復(fù)雜邊界的彈性體總可以劃分為三角形，并且三角形單元計算公式簡單，但精度較低。三角形單元間雖然能夠保證位移連續(xù)，但應(yīng)力的精度較差，不能很好的反映彈性體內(nèi)應(yīng)力的準(zhǔn)確分布規(guī)律。為了提高計算精度，準(zhǔn)確反映彈性體內(nèi)的應(yīng)力狀態(tài)，可以采用一些較精密的單元類型。本節(jié)將介紹常用的矩形單元，它采用了比常應(yīng)變?nèi)切螁卧螖?shù)更高的位移模式，因而可以更好地反映彈性體中的位移狀態(tài)和應(yīng)力狀態(tài)。另外，對一些邊界比較規(guī)則且呈直線的平面結(jié)構(gòu)的分析，采用矩形單元較合適。這時單元總數(shù)可以減少，相應(yīng)的原始數(shù)據(jù)準(zhǔn)備工作和單元特征計算工作均可節(jié)省。3.2.1位移函數(shù)如圖所示的矩形單元，不失一般性，令矩形單元的長、寬分別為2a、2b。矩形單元有4個節(jié)點(diǎn)，共8個自由度，即共有8個節(jié)點(diǎn)位移，采用類似三角形單元的分析方法，同樣可以完成對矩形單元的力學(xué)特性分析。這里引入一個局部坐標(biāo)系、，這樣可以推出比較簡潔的結(jié)果。如圖所示，取矩形單元的形心o為局部坐標(biāo)系的原點(diǎn)，和軸分別與整體坐標(biāo)軸x和y平行，兩坐標(biāo)系存在有以下的坐標(biāo)變換關(guān)系式中：、——矩形形心處坐標(biāo)。矩形形心處坐標(biāo)以及矩形長、寬可由下式計算3.2.1位移函數(shù)在局部坐標(biāo)系中，節(jié)點(diǎn)i的坐標(biāo)是，其值分別為±1。如節(jié)點(diǎn)1在局部坐標(biāo)系下的坐標(biāo)為（－1，－1）。由于矩形有4個節(jié)點(diǎn)，共8個自由度，可以選擇有8個待定參數(shù)的位移模式，如下采用類似三角形的分析方法，將各節(jié)點(diǎn)的局部坐標(biāo)代入上式，便可以求得相應(yīng)的待定參數(shù)，然后將帶回式位移模式，可得到用節(jié)點(diǎn)位移表示的位移模式，如下3.2.1位移函數(shù)式中：——矩形單元的形函數(shù)，i＝1，2，3，4；——形函數(shù)矩陣；——單元節(jié)點(diǎn)位移列陣，，i＝1，2，3，4。3.2.1位移函數(shù)（i＝1，2，3，4）形函數(shù)的表達(dá)式為3.2.1位移函數(shù)引入符號，，i＝1，2，3，4，則上式可以統(tǒng)一寫為可以看出，矩形單元的形函數(shù)具有和三角形單元形函數(shù)同樣的性質(zhì)，即：形函數(shù)在各單元節(jié)點(diǎn)上的值，具有“本點(diǎn)是1、它點(diǎn)為零”的性質(zhì)；在單元內(nèi)任意點(diǎn)上，四個形函數(shù)之和等于1；單元任意一條邊上的形函數(shù)，僅與該邊的兩端節(jié)點(diǎn)坐標(biāo)有關(guān)。3.2.1位移函數(shù)有了單元的位移模式，就可以利用平面問題的幾何方程求出單元內(nèi)任意點(diǎn)的應(yīng)變，將位移代入幾何方程，得式中的應(yīng)變轉(zhuǎn)換矩陣的子塊（i＝1，2，3，4）為3.2.2應(yīng)變與應(yīng)力矩陣求得應(yīng)變之后，再將應(yīng)變代入物理方程，便可推導(dǎo)出以節(jié)點(diǎn)位移表示的應(yīng)力，如下式中，應(yīng)力矩陣為其子塊（i＝1，2，3，4）為3.2.2應(yīng)變與應(yīng)力矩陣由前面的討論可以發(fā)現(xiàn)，四邊形單元的位移模式，式（3.77）比常應(yīng)變?nèi)切螁卧捎玫木€性位移模式增添了項(xiàng)（即相當(dāng)于xy項(xiàng)），把這種位移模式稱為雙線性模式。在這種模式下，單元內(nèi)的應(yīng)變分量將不再是常量，這一點(diǎn)可以從的表達(dá)式中看出。另外四邊形單元的位移模式中的與三角形單元相同，它反映了剛體位移和常應(yīng)變，而且在單元的邊界上(=±1或

=±1)，位移是按線性變化的，顯然在兩個相鄰單元的公共邊界上，其位移是連續(xù)的。3.2.2應(yīng)變與應(yīng)力矩陣由單元的應(yīng)力矩陣表達(dá)式還可以看出，矩形單元中的應(yīng)力分量也都不是常量。正應(yīng)力、和剪應(yīng)力均沿、兩個方向線性變化，即沿x、y兩個方向線性變化。正因?yàn)槿绱耍粼趶椥泽w中采用相同數(shù)目的節(jié)點(diǎn)時，矩形單元的精度要比常應(yīng)變?nèi)切螁卧木雀摺５?，矩形單元也有一些明顯的缺點(diǎn)，矩形單元不能適應(yīng)斜交的邊界和曲線邊界，不便于對不同部位采用不同大小的單元，以便提高有限元分析計算的效率和精度。3.2.2應(yīng)變與應(yīng)力矩陣矩形單元剛度矩陣的推導(dǎo)過程與三節(jié)點(diǎn)三角形單元類似，同樣可以表示為式（3.30）的形式，即，由前文可知的推導(dǎo)過程與形函數(shù)的具體表達(dá)形式、節(jié)點(diǎn)個數(shù)均無關(guān)，該表達(dá)式具有普遍意義。若單元厚度t

為常量，則可以進(jìn)一步表示為將單元剛度矩陣寫成子塊的形式，如下3.2.3單元剛度矩陣上式中每一個子塊矩陣均為2行×2列，單元剛度矩陣中的子塊矩陣的表達(dá)式為

（r、s＝1，2，3，4）將應(yīng)變轉(zhuǎn)換矩陣子塊和彈性矩陣，代入上式，得3.2.3單元剛度矩陣式中：（r、s＝1，2，3，4）3.2.3單元剛度矩陣?yán)?.6如圖3.22所示，該模型中有兩個四邊形單元，彈性模量為＝210GPa，厚度為＝0.025m，泊松比＝0.3，＝1kN，求單元所受應(yīng)力。3.2.4算例解：單元①所對應(yīng)的節(jié)點(diǎn)為2、4、3、1，單元②所對應(yīng)的節(jié)點(diǎn)為4、6、5、3。在有限元分析過程中，首先求解單元①和②的剛度矩陣；然后組裝整體剛度矩陣，組裝的過程同三角形單元，此處不再給出；最后引入邊界條件（，＝0）并結(jié)合受力情況，求得整體節(jié)點(diǎn)位移列陣＝10-5×{0，0，0，0，0.1162，-0.1674，-0.1149，-0.1628，0.1514，-0.4707，-0.1568，-0.4978}T。為了求解單元應(yīng)力，需要先求得每個單元所對應(yīng)的節(jié)點(diǎn)位移列陣，結(jié)合單元的節(jié)點(diǎn)編號，可從整體位移列陣中提取單元①和單元②的位移列陣，如下＝10-5×{0，0，-0.1149，-0.1628，0.1162，-0.1674，0，0}T＝10-5×{-0.1149，-0.1628，-0.1568，-0.4978，0.1514，-0.4707，0.1162，-0.1674}T3.2.4算例注意：整體節(jié)點(diǎn)位移列陣是按照節(jié)點(diǎn)編號由小到大排列的，而單元位移列陣是按照單元節(jié)點(diǎn)編號排列的，如單元①所對應(yīng)的節(jié)點(diǎn)為2、4、3、1，則單元①的位移列陣中的前兩個數(shù)則表示節(jié)點(diǎn)2的x和y方向位移。將單元①和單元②的位移列陣代入應(yīng)力矩陣，可求得單元①和單元②的應(yīng)力，如下3.2.4算例通過以上兩式便可以求得單元①和單元②內(nèi)任意點(diǎn)的應(yīng)力，以單元①為例，若取，，則表示單元①的13邊中點(diǎn)處的應(yīng)力；若取，，則表示24邊中點(diǎn)處的應(yīng)力。若計算單元形心處的應(yīng)力，則取，為通過分析結(jié)果可知，單元內(nèi)任意點(diǎn)的應(yīng)力是坐標(biāo)的函數(shù)，若在彈性體中采用相同數(shù)目的節(jié)點(diǎn)時，矩形單元的精度顯然要高于常應(yīng)變?nèi)切螁卧木取?.2.4算例3.3平面等參元前一節(jié)討論的矩形單元具有精度較高、形狀規(guī)則、便于計算自動化等優(yōu)點(diǎn)，對于結(jié)構(gòu)復(fù)雜的曲邊外形，只能通過減小單元尺寸，增加單元數(shù)量進(jìn)行逐漸逼近。這樣，自由度的數(shù)目隨之增加，計算時間長，工作量大。另外，這些單元的位移模式是線性模式，是實(shí)際位移模式的最低級逼近形式，問題的求解精度受到限制。為了克服以上缺點(diǎn)，人們試圖找出這樣一種單元：一方面，單元能很好地適應(yīng)曲線邊界和曲面邊界，準(zhǔn)確地模擬結(jié)構(gòu)形狀；另一方面，這種單元要具有較高次的位移模式，能更好地反映結(jié)構(gòu)的復(fù)雜應(yīng)力分布情況，即使單元網(wǎng)格劃分比較稀疏，也可以得到比較好的計算精度。等參數(shù)單元（等參元）就具備了以上兩條優(yōu)點(diǎn)，等參元為任意四邊形單元，其網(wǎng)格劃分不受邊界形狀的限制，單元大小可以不相等，是一種精度高而且應(yīng)用廣泛的單元。3.3平面等參元然而直接對任意四邊形進(jìn)行單元分析是困難的，這是由于它的幾何形狀不規(guī)則，沒有統(tǒng)一的形狀，對各個單元逐個按不同的公式計算，因其工作量大而難以進(jìn)行。為了解決這一問題，人們采用坐標(biāo)變換的方法，如圖所示，把坐標(biāo)系ξη（局部坐標(biāo)系）內(nèi)形狀規(guī)則單元（母單元）變換為另一坐標(biāo)系xy（整體坐標(biāo)系）中的任意四邊形單元（子單元）。然后導(dǎo)出關(guān)于局部坐標(biāo)系形狀規(guī)整的單元（母單元）的高階位移模式的形函數(shù)，然后利用形函數(shù)進(jìn)行坐標(biāo)變換，得到關(guān)于整體坐標(biāo)系的復(fù)雜形狀的單元（子單元）。a）母單元b）子單元3.3平面等參元從圖形變換的角度看，坐標(biāo)系ξη和坐標(biāo)系xy可以分別看成是母單元和子單元這兩個不同單元的坐標(biāo)系，它們都是直角坐標(biāo)系。而從另一角度看，坐標(biāo)系ξη和坐標(biāo)系xy又可以看成是同一單元（子單元）的兩種不同的坐標(biāo)系。坐標(biāo)系xy是子單元的直角坐標(biāo)系，而坐標(biāo)系ξη可看成是子單元的曲線坐標(biāo)系?？梢钥闯鲎鴺?biāo)系xy始終扮演同一角色，即子單元的直角坐標(biāo)；而坐標(biāo)系ξη則扮演兩種角色，它既是母單元的直角坐標(biāo)，又是子單元的曲線坐標(biāo)。在有限元分析中，兩者的作用是不同的。直角坐標(biāo)系xy在整個結(jié)構(gòu)的所有子單元中共同采用，所以稱為整體坐標(biāo)系。而曲線坐標(biāo)系ξη則只適用于單個獨(dú)立的子單元，所以稱為局部坐標(biāo)。整體坐標(biāo)在整體分析中采用，局部坐標(biāo)則在單元分析中采用。平面問題中常用的等參元有四節(jié)點(diǎn)四邊形單元、八節(jié)點(diǎn)曲邊四邊形單元和6~8可變節(jié)點(diǎn)曲邊四邊形單元等，本節(jié)以八節(jié)點(diǎn)曲邊四邊形等參元為例介紹平面等參元的分析過程。八節(jié)點(diǎn)曲邊四邊形等參元的母單元是二維二次單元。八個節(jié)點(diǎn)分別為正方形的四個角點(diǎn)和四個邊中點(diǎn)，母單元（邊長為2的正方形）采用直角坐標(biāo)系ξη。由于母單元為16個自由度，因此單元的位移模式為3.3.1等參元剛度矩陣式中：式中的16個常數(shù)用8個節(jié)點(diǎn)的位移（，）表示后，則上式為3.3.1等參元剛度矩陣采用坐標(biāo)變換可使母單元的八個節(jié)點(diǎn)坐標(biāo)與等參元的八個節(jié)點(diǎn)坐標(biāo)建立一一對應(yīng)的關(guān)系，整體坐標(biāo)和局部坐標(biāo)的變換式為式中：——母單元的形函數(shù)?？芍竼卧先我庵边?，例如2－3邊，有ξ＝1，則有關(guān)的形函數(shù)N2、N3和N6均是η的二次函數(shù)，其余的形函數(shù)均為零。這樣變換就成為二次非線性變換，它可將母單元的直邊2－3映射成子單元的曲邊2－3。根據(jù)等參元的思想，利用上述的形函數(shù)，可得等參元的位移模式為3.3.1等參元剛度矩陣有了單元的位移模式，就可以利用平面問題的幾何方程求出單元內(nèi)任意點(diǎn)的應(yīng)變，將上式代入式（2.12），得等參元的應(yīng)變矩陣為式中：——i節(jié)點(diǎn)位移，（i=1，2，…，8）；3.3.1等參元剛度矩陣——單元幾何矩陣，其中為（i=1，2，…，8）在上式中，是ξ、η的函數(shù)，因此必須用坐標(biāo)變換式來轉(zhuǎn)換導(dǎo)數(shù)關(guān)系，根據(jù)復(fù)合函數(shù)的求導(dǎo)法則，有（i=1，2，…，8）3.3.1等參元剛度矩陣式中的J稱為雅可比（Jacobi）矩陣，為則3.3.1等參元剛度矩陣雅可比（Jacobi）矩陣的逆矩陣為將單元應(yīng)變代入平面問題的物理方程式，就得到平面八節(jié)點(diǎn)等參元的應(yīng)力列陣，為（i=1，2，…，8）3.3.1等參元剛度矩陣等參元剛度矩陣仍可用式（3.30）計算，即在上式中，應(yīng)變矩陣為局部坐標(biāo)系ξη下的顯式，因此在計算上式時應(yīng)將微面積dxdy變換成dξ和dη所包圍的面積dA，由于所以，有3.3.1等參元剛度矩陣則等參元剛度矩陣為應(yīng)該指出，上式是對ξ和η的重積分，盡管其積分區(qū)域十分簡單，但其被積函數(shù)卻比較復(fù)雜，一般該積分沒有顯式，需要采用數(shù)值積分法求解。至于等參元的等效節(jié)點(diǎn)力的計算一般也要采用數(shù)值積分的方法求解，這里不再討論，可參考相關(guān)書籍。由以上分析可知，在劃分單元時，只需確定單元節(jié)點(diǎn)的整體坐標(biāo)值，而不必畫出等參元的具體形狀，因?yàn)樵谟嬎阒袑?shí)際使用的只有單元八個節(jié)點(diǎn)在整體坐標(biāo)系下的位移值。等參元變換的條件為，因此在有限元網(wǎng)格劃分時，要特別注意這一點(diǎn)。3.3.1等參元剛度矩陣如圖所示，在圖b中3、4點(diǎn)退化為一個點(diǎn)，在該點(diǎn)，同理在圖c中2、3點(diǎn)退化為一個點(diǎn)，在該點(diǎn)，在圖d在點(diǎn)1、2、3處>0，在點(diǎn)4處<0，又因?yàn)樵趩卧獌?nèi)連續(xù)變換，所以單元內(nèi)肯定存在＝0，這是由于單元過分歪曲造成的。a）正常b）不正常c）不正常