1、第二章第二章 結構分析的有限元法結構分析的有限元法 2.1 2.1 有限元法發展簡況有限元法發展簡況 利用定義在三角形區域上的分片連續函數和最小位能原理St.Venant扭轉問題的近似解有限元法的研究現代有限元法1943Courant應用數學家、物理學家、工程師1960Tumer、Clough第一次用三角形單元平面應力問題解答提出了有限單元法的名稱各種非線性問題多物理場耦合問題多尺度問題商品化有限元軟件20世紀70年代國外幾何非線性：因幾何變形引起結構剛度改變材料非線性：彈性（超彈和多線性彈性）、粘彈性、非彈性狀態非線性：接觸問題2.1 有限元法發展簡況 固體力學流體力學傳熱學電磁學學科應用力

2、學計算結構優化計算功能計算技術純粹數值技術前、后處理技術的高度智能化和與CAD的集成化2.2 2.2 有限元法的基本思路及其求解步驟有限元法的基本思路及其求解步驟 經典的解析法 從連續體的微分方程入手，尋求滿足微分方程和定解條件的適合全域的解析解，一旦得到解析解，就可知道域內任意點的解大多數問題，特別是實際問題 很難甚至無法用解析法得到問題的解析解在整個求解域上滿足控制方程在邊界上滿足邊界條件的場函數尋找很困難很困難有限元法有限元法單元節點有限元模型2.2 有限元法的基本思路及其求解步驟 有限元法基本思路拋棄尋找一個滿足整個求解域的場函數的思路把求解域劃分成有限個四邊形單元對每一個單元通過插值

3、的方法，用其節點上的位移建立該單元的位移函數123每個單元都有與其對應的位移函數表達式用全部單元域之和代替整個求解域，用全部單元的位移函數之和代替滿足整個求解域的位移函數4對單元進行力學特性分析，建立單元節點力與單元節點位移的關系，并將結構的外載荷等效移植到節點上，再在節點上建立力的平衡方程，求解后得到節點上的位移，繼而得到各個單元的應力5以以節點位移節點位移為未知量，通過為未知量，通過求解求解力的平衡方程力的平衡方程獲得節點獲得節點位移，然后按位移，然后按單元單元計算應力計算應力2.2 有限元法的基本思路及其求解步驟 有限元法求解步驟離散化離散化將結構（求解域）劃分為有限個單元，讓全部單元的

4、集合與原結構近似等價劃分單元時，二者在幾何形體上越逼近越好，特別是在位移和應力急劇變化的地方選擇單元位移函數選擇單元位移函數在有限元法中，需要用單元節點位移通過插值方法建立單元位移函數（單元位移模式），即用單元節點位移來描述單元位移。單元位移函數的合理與否，直接關系到有限元分析的計算精度、效率和收斂性。通常取為多項式形式2.2 有限元法的基本思路及其求解步驟 單元特性分析單元特性分析（1）依照應變與位移之間的幾何關系，根據所選擇的單元位移函數，建立單元應變與單元節點位移之間的關系式。據此式，在求出節點位移后，可以求得單元應變。（2）依照物理關系（胡克定律），建立單元應力與單元節點位移之間的關系

5、式。據此式，在求出節點位移后，可以求得單元應力。（3）根據虛位移原理或最小勢能原理，建立單元剛度方程，即單元節點力與單元節點位移之間的關系式。此步驟核心是計算單元剛度矩陣。外載荷處理外載荷處理將外載荷（體力、面力等）等效移植到節點上。2.2 有限元法的基本思路及其求解步驟 建立節點上的力平衡方程建立節點上的力平衡方程按照有限元法的統一格式，形成如下形式的以節點位移為未知量的代數方程組 KF K由各個單元的剛度矩陣組裝成的總體剛度矩陣 待求的節點位移列陣 F按節點編號順序形成的節點載荷列陣處理邊界條件、解算節點位移處理邊界條件、解算節點位移(2.1)按照實際位移邊界條件，對式(2.1)進行整理，

6、解之，可得單元節點位移。有了節點位移，即可根據單元特性分析中建立的關系式，求應力、應變、內力等。有了節點位移，即可根據單元特性分析中建立的關系式，求應力、應變、內力等。后處理后處理：對所選應力、應變等，以：對所選應力、應變等，以彩色云圖彩色云圖或或圖表圖表的形式顯示的形式顯示計算結果。計算結果。2.3 2.3 有限元程序的結構簡介有限元程序的結構簡介 對一個題目或一個實際工程問題進行有限元分析，大體上對一個題目或一個實際工程問題進行有限元分析，大體上分分3個個主要步驟主要步驟有限元建模有限元求解計算結果分析與整理前處理求解器后處理程程序序結結構構2.3 有限元程序的結構簡介 前處理幾何模型的建

7、立定義約束條件網格剖分確定材料參數和載荷有限元建模形成有限元分析所需用的有限元計算數據形成有限元分析所需用的有限元計算數據可視化可視化 有限元模型有限元模型前處理中，可以用圖形顯示所建立的幾何模型、單元網格、約束條件等求解器2.3 有限元程序的結構簡介 有限元程序的核心部分主要完成有限元模型的力學計算，即根據前處理形成的有限元計算數據，完成以下工作：計算單元剛度矩陣計算節點載荷組裝總體剛度矩陣將載荷等效簡化到節點上形成總體有限元平衡方程求解節點位移計算應力、應變、內力等2.3 有限元程序的結構簡介 后處理根據計算者的要求對計算結果進行檢查、分析、整理、打印輸出等進行數據檢索響應量合成繪制變形圖

8、、應力圖、應變圖、曲線圖等可視化可視化的方式分析、觀察計算結果的方式分析、觀察計算結果計算者計算者進行有限元分析的工作量主要體進行有限元分析的工作量主要體現在現在前處理前處理和和后處理后處理方面方面2.4 2.4 算例算例 2.4.1 平面三角形常應變單元任意區域三角形單元網 格剖分示意圖mvmujvjuiviuijm典型三角形單元yxo單元內任意點(x,y)的位移uv、坐標x和y的函數建立單元位移函數通過插值方法建立，即用單元的節點位移來表示單元內任意點的位移1. 單元位移函數2.4 算例 mvmujvjuiviuijm典型三角形單元yxo單元位移函數選用坐標x和y的一次多項式123456u

9、xyvxy(1)123456、待定系數123123123iiijjjmmmuxyuxyuxy(2)未知量求解(2)(3)得到123456、123456、456456456iiijjjmmmvxyvxyvxy(3)2.4 算例 123123123iiijjjmmmuxyuxyuxy(2)1121iijjmmxyDxyAxy123112111121111121iiijjjiijjmmmmmiijjiijjmmmmiijjiijjmmmmuxyuxyaua ua uDAuxyuyuybub ub uDAuyxuxucuc uc uDAxu456456456iiijjjmmmvxyvxyvxy(3)1

10、121iijjmmxyDxyAxy456112111121111121iiijjjiijjmmmmmiijjiijjmmmmiijjiijjmmmmvxyvxyava va vDAvxyvyvybvb vb vDAvyxvxvcvc vc vDAxv2.4 算例 ijmmjjmiimmijjiax yx yax yx yax yx yijmjmimijbyybyybyyijmjmimijcxxcxxcxx 123456uxyvxy(1)代入123456、將求得的得到用單元的節點位移表示的單元位移函數iijjmmiijjmmuN uN uN uvN vN vN v(4)式中是單元形狀函數，簡稱形

11、函數121212iiiijjjjmmmmNab xc yANab xc yANab xc yA（5）,ijmN NN是常數，取決于單元的三個節點坐標, , ,iiima b cc，返回P242.4 算例 11211111.()( 1).() 1.()()222iijjmmjmmjimmiijjiijmxyDxyAxyADx yx yx yx yx yx yaaa 1.()( 1).()1.()()()()()()()0.jmmjimmiijjijmimijijmjjmmiijjjmjmmmiiimimjjmiimjmimjiijmijjix yx yx yx yx yx yxyyxyyx yy

12、b xb xb xb xb xb xb xb xb xb xbxxb xxbbb xb cbcxbcb c 1()2ijjiAbcb c三角形單元的面積A單元位移函數表達式iijjmmiijjmmuN uN uN uvN vN vN v(4)2.4 算例 000000iijeijmijmjmmuvuNNNuuNNNvvuv eeuN寫成矩陣形式簡寫為其中， eu表示單元內任意點處位移的單元位移函數列陣 000000ijmijmNNNNNNN為形函數矩陣返回P282.4 算例 形函數的性質在節點上形函數的值是式（6）表示形函數Ni在其自身節點上的值等于1，在其他節點上的值等于0，即1(,)( ,

13、 ,)0ijjijjiN xyi j mji（6）( ,)1(,)0(,)0iiiijjimmN x yN xyN xy，( ,)0(,)1(,)0jiijjjjmmNx yNxyNxy，( ,)0(,)0(,)1miimjjmmmNx yNxyNxy，單元中任意一點上的各個形函數之和等于1，即2.4 算例 ( , )( , )( , )1ijmN x yNx yNx y由121212iiiijjjjmmmmNab xc yANab xc yANab xc yA（5） ( , )( , )( , )11122212ijmiiijjjmmmijmijmijmN x yNx yNx yab xc

14、yab xc yab xc yAAAaaabbbxcccyAijmjmimijbyybyybyyijmjmimijcxxcxxcxx 1()2ijmAaaa120.0.21AxyA2.4 算例 小結小結（1 1）本節的三角形單元，形函數是）本節的三角形單元，形函數是線性線性的，為的，為x x、y y的一的一次函數；次函數；（2 2）在單元內部和各條單元邊上，位移也是）在單元內部和各條單元邊上，位移也是線性線性的，可的，可由兩個節點的位移由兩個節點的位移唯一唯一確定；確定；（3 3）相鄰單元的公共節點的節點位移是相等的，因此，）相鄰單元的公共節點的節點位移是相等的，因此，能保證相鄰單元在公共邊界

15、上以及單元內部的能保證相鄰單元在公共邊界上以及單元內部的位移連續位移連續性性。2.4 算例 單元位移函數確定后，根據幾何方程xyxyuxvyuvyx求得單元內任意點處的應變，即單元應變 000000ijimixjejimyjxyjjiimmmmuNNNuvxxxxuNNNvvyyyyuvNNNNNNuyxyxyxyxv（7）iijjmmiijjmmuN uN uN uvN vN vN v(4)2. 單元應變和單元應力2.4 算例 由121212iiiijjjjmmmmNab xc yANab xc yANab xc yA（5）形函數對坐標變量求偏導111222111222jimijmjimij

16、mNNNbbbxAxAxANNNcccyAyAyA，（8），式（8）代入式（7）中，得到 0001000(9)2iixijmjeeyijmjiijjmmxymmuvbbbucccBvAcbcbcbuv =返回P28由式（9）得到 0001000(10)2ijmijmijmiijjmmbbbBcccBBBAcbcbcb2.4 算例 單元應變矩陣（幾何矩陣）000111000222ijmiijjmmiijjmmbbbBcBcBcAAAcbcbcb，（11）分塊矩陣參數,; ,ijmijmb b bc c c由單元的節點坐標確定，因此，它們取決于單元形狀，當單元的節點坐標確定后，它們都是常量，所以，

17、3節點三角形單元的應變矩陣B是常數矩陣2.4 算例 根據物理方程 xeeeeyxyDDBS（12）其中 000200101011002ED為平面應力（平面應變）問題的彈性矩陣平面應力問題00,EE平面應變問題002,11EE 0000000200000002(1)111111222222iijjmmiijjmmiijjmmbcbcbcESDBbcbcbcAcbcbcb（13）應力矩陣S也是常數矩陣單元應力矩陣返回P28（1 1）3 3節點平面三角形單元，節點平面三角形單元，應變矩陣應變矩陣和和應力矩陣應力矩陣都為都為常數矩陣常數矩陣；（2 2） 3節點平面三角形單元，各點的應變和應力都是相節點

18、平面三角形單元，各點的應變和應力都是相同的，且是常數，所以同的，且是常數，所以3節點三角形單元是節點三角形單元是常應變單元常應變單元，也是也是常應力單元常應力單元；（3 3）采用）采用3 3節點三角形單元時，在應力變化劇烈或應力節點三角形單元時，在應力變化劇烈或應力梯度較大的部位，單元劃分應適當加密。梯度較大的部位，單元劃分應適當加密。2.4 算例 小結小結3. 單元剛度矩陣2.4 算例 分析結構在載荷作用下產生變形和應力，于是在各單元之間就產生相互作用。實際上，各單元之間的相互作用是通過相鄰邊界上（即，單元的邊，實際是面）的分布力而產生的。按照有限元方法，結構離散化為一個個單元后，單元之間的

19、相互作用就由單元的節點力來實現，即用單元節點力等效代替相鄰邊界上的相互作用力，這樣，節點力就與單元應力相關，而單元應力與節點位移相關，因此，單元節點力與單元節點位移相關。建立單元節點力與單元節點位移之間的關系 eR表示單元節點力mvmujvjuiviuijm (14)Teiijjmmuvuvuu (15)TeixiyjxjymxmyRRRRRRR單元節點位移由（9）式2.4 算例 (9)eeB= eeDB（12）由（12）式把一個單元作為分析對象時，可以把節點力看作外力。單元節點力和單元節點位移之間的關系可由虛位移原理導出在外力作用下，處于平衡狀態的變形體，當發生約束允許的任意微小的虛位移時，

20、外力在虛位移上所做的虛功等于整個體積內的應力在虛應變上所做的虛功。推導令單元的節點虛位移為 *eTiijjmmuvuvuv eeuN由P18*(16)eeuN eeB由(9)式 *(17)eeB2.4 算例 (15)TeixiyjxjymxmyRRRRRRR *eTiijjmmuvuvuv節點力在虛位移上所做的虛功為 *(18)Teeiixiiyjjxjjymmxmmyu Rv Ru Rv Ru Rv RR單元應力在虛應變上所做的虛功為 *(19)eTeeVdV單元體積 *(17)eeB將 eeDB（12）代入(19)式 *eeeTeeVTTeeTeeVVdVBDBdVBDBdV由于節點虛位移

21、是任意的 *eTeTeVBDB dV建立單元的虛功方程為 *(20)eTTeeTeeVRBDB dV由單元的虛功方程 *(20)eTTeeTeeVRBDB dV2.4 算例 任意性相等 (21)eTeeVRBD B dV節點力和節點位移之間的關系式就建立起來了 (22)eeTVKBDB dV令 (23)eeeRK則(21)式變為單元剛度方程注意：這里，節點力不是結構上的外載荷，而是按虛位移原理把單元邊界上的分布力近似等效到單元節點上的一種節點力。節點力在實際結構中是不存在的總結：式（22）、（23）是由三角形常應變單元推導得到的，但是，這兩式及其推導過程所基于的原理和方法具有普遍性。原則上說（

22、22）式是位移有限元分析中普遍適用的單元剛度矩陣表達式，對于不同單元，只是其中的具體計算細節不同。三角形常應變單元的剛度矩陣分析2.4 算例 0001000(10)2ijmijmiijjmmbbbBcccAcbcbcb一般情況，單元應變矩陣B是坐標的函數矩陣。 000200101011002ED這里，三角形常應變單元，B是常數矩陣。如果材料是線性的、勻質的，矩陣D也是常數矩陣。單元厚度t是常量,則dV=tdxdy,因此，三角形常應變單元的剛度矩陣可以寫成： (24)eeTTTTVAAKBDB dVBDB tdxdyBDB tdxdyBDB tA=將單元剛度矩陣寫成分塊形式： 222221111

23、112222221111112222221124(1)eTiiiiiiijijijijimimimimiiiiiiijijijijimimimimjijijiKBDB tAbcbccbbbccbccbbbccbccbcbbccbcbbcccbbcbbcccbbb bc cb cEtA222221111222221111112222221111122222jijjjjjjjmjmjmjmjijijijijjjjjjjmjmjmjmm im im im imjmjmjmjmmc bbcb cc bb bc cb cc bc bb cc cb bb cc bcbc bb cc cb bb bc cb

24、 cc bb bc cb cc bbc22212111111222222immmmm im im im imjmjmjmjmmmmmmiiijimjijjjmmimjmmcb cc bc bb cc cb bc bb cc cb bc bb ccbkkkkkkkkk（25）2.4 算例 21122, ,114(1)22rsrsrsrsrsrsrsrsrsb bc cb cc bEtkr si j mAc bb cc cb b2.4 算例 （1 1）對稱性對稱性單元剛度矩陣是單元剛度矩陣是對稱矩陣對稱矩陣；（2 2）奇異性奇異性單元剛度矩陣是單元剛度矩陣是奇異矩陣奇異矩陣，它不存在逆矩，它不存在

25、逆矩陣陣；（3 3）主元恒正主元恒正單元剛度矩陣對角元素的數值恒大于單元剛度矩陣對角元素的數值恒大于0 0，可由（可由（2525）式看出；）式看出；（4 4）單元剛度矩陣的元素具有明確的物理意義，如圖所示）單元剛度矩陣的元素具有明確的物理意義，如圖所示單元，其剛度矩陣第一列的單元，其剛度矩陣第一列的6 6個元素個元素k ki1i1( (i i=1,2,3,4,5,6)=1,2,3,4,5,6)的物理意義是，當單元的的物理意義是，當單元的第第1 1個節點位移個節點位移( (節點節點i i的的u ui i) )為為1 1，而而其他節點位移其他節點位移全為全為0 0時，需要在時，需要在6 6個節點位

26、移方向上施加個節點位移方向上施加的的節點力節點力的的大小大小。單元剛度矩陣的特性單元剛度矩陣的特性mvmujvjuiviuijm2.4 算例 4. 單元等效節點載荷應用虛位移原理進行載荷等效移植：移植后的節點載荷和移植前的載荷在約束允許的任意虛位移上所做的功相等。單元的節點虛位移為 *eTiijjmmuvuvuv*(16)eeuN單元的虛位移為劃分單元時，一般都將作用有集中力的地方劃分為節點，集中力即可直接施加到節點上。以下說明體積力和表面力向節點移植：F2.4 算例 (1)體積力等效移植體積力等效移植 令單位體積的力為： (26)Txyggg單位體積力在x軸和y軸方向的分量令單元體積力等效的

27、移植到單元節點上的等效節點載荷為： (27)eTgigxigyjgxjgymgxmgyFFFFFFF由虛位移原理， egF g和在虛位移上所做的虛功相等，即 *(28)TTeeegAFug tdxdy單元面積單元厚度(16)式代入(28)式： *TTTTeeeeeTgAAAFug tdxdyNg tdxdyNg tdxdy2.4 算例 000000ixiyTeTjxTijmgxyAAAijmjymxmyN gN gN gNNNFNg tdxdytggdxdytdxdyNNNN gN gN g(29)特殊地，體積力是重力，且重力方向為負y方向，單元的單位體積力是 0TTxyggg= 000000

28、00100113001TeTTijmgxyAAijmixiyijxAAjjymxmmyNNNFNg tdxdytggdxdyNNNN gN gNN gtdxdytdxdytANN gN gNN g (30)2.4 算例 (2)表面力等效移植表面力等效移植 2.4 算例 工程問題中，表面力一般都垂直于其作用面，所以，在有限元法中，要求定義的表面力垂直于其作用面，這樣，可以將表面力分解到沿坐標軸方向。令q為表面力矢量，則它可以表示為： (31)Txyqqq表面力在x軸和y軸方向的分量令表面力等效的移植到單元節點上的等效節點載荷為： (32)eTqiqxiqyjqxjqymqxmqyFFFFFFF由

29、虛位移原理， eqF q和在虛位移上所做的虛功相等，即 *(33)eTTeeeqSFuq tds表面力作用的邊界*(16)eeuN單元的虛位移為代入(33)式得到： *eeeTTTTeeeeeTqSSSFuq tdsNq tdsNq tds 000000eeeixiyTeTjxTijmqxySSSijmjymxmyN qN qN qNNNFNq tdsqqtdstdsNNNN qN qN q(34)2.4 算例 2.4 算例 xy0ijmq沿單元邊界均勻分布的表面力mi邊上有垂直于邊界的均勻分布的表面力q,將分解為x和y方向的均布力qx和qy，這樣，mi邊上的表面力可以表示為： 11sinco

30、s11mijxyimjyybqqllqqqqqxxcll邊的邊長為mil 002eTTxqxyxylyqtlFNtdsqqqqq(35)2.4 算例 5. 總體平衡方程的建立厚度為t的正方形板，左邊固定，在其右上角分別作用有x、y方向的集中力F1x，F1y，建立其總體有限元平衡方程。(1)劃分單元，所有節點總體編號1xF1yF121234(2)對號入座12324123( ) i( ) i( ) j( )m( )m單元號局部編碼整體編碼i14j32m2312( ) j單元剛度矩陣，其分塊形式為2.4 算例 11111111112222222222iiijimjijjjmmimjmmiiijimj

31、ijjjmmimjmmkkkkkkkkkkkkkkkkkkkk整體剛度矩陣，其分塊形式為 1111112131411212221222324112122313233342224142434400iiimijmimmjjmjjmjijijmmjjjmmmiijimiikkkKKKKkkkkkkKKKKKKKKKkkkkkkKKKKkkk單元號局部編碼整體編碼i14j32m23122.4 算例 (3)總體平衡方程 1111111212222112122332224400iiimijmimmjjmjjmjijijmmjjjmmmiijimiikkkUFkkkkkkUFUFkkkkkkUFkkk 1x

32、F1yF121234單元號局部編碼整體編碼i14j32m2312 KF總體平衡方程 K整體剛度矩陣 整體節點位移 F整體節點力2.4 算例 1111111212222112122332224400iiimijmimmjjmjjmjijijmmjjjmmmiijimiikkkUFkkkkkkUFUFkkkkkkUFkkk 11111121121222211212233222344400iiimijmimmjjmjjmjijijmmjjjmmmiijimiiuvkkkUukkkkkkUvKUuUkkkkkkvUkkkuv， 111222333444xyxyxyxyFFFFFFFFFFFFF，1xF1yF1212346. 位移邊界條件的處理及總體平衡方程的求解2.4 算例 (1)零位移邊界條件的處理1112131415161718212223242526272831323334353637384142434445464748515253545556575861626364656667687172737