數(shù)值模擬：第五講平面問題(二)——三角形單元分析

1、5.1 離散化要解決的問題離散化要解決的問題如何進行離散結(jié)構(gòu)的求解？如何進行離散結(jié)構(gòu)的求解？如何得到小單元的剛度特性？如何得到小單元的剛度特性？為什么離散化結(jié)構(gòu)的解能作為原連續(xù)問題的近似解？為什么離散化結(jié)構(gòu)的解能作為原連續(xù)問題的近似解？ 研究上述問題就涉及到有限元法的基本原理和基本理論。研究上述問題就涉及到有限元法的基本原理和基本理論。5 5 三節(jié)點三角形單元解平面問題三節(jié)點三角形單元解平面問題 目 標: 掌握平面問題 簡單三角形單元位移模式。5.2 5.2 三節(jié)點三角形單元的特性分析三節(jié)點三角形單元的特性分析5.2.1 單元作為分析對象單元作為分析對象按前面結(jié)構(gòu)矩陣位移法分析思想，要求解按前

2、面結(jié)構(gòu)矩陣位移法分析思想，要求解平面問題的有限元離散結(jié)構(gòu)，需要知道單平面問題的有限元離散結(jié)構(gòu)，需要知道單元（三角形薄片）在節(jié)點自由度上受力時元（三角形薄片）在節(jié)點自由度上受力時的彈性特性或剛度特性。這是一個新問題，的彈性特性或剛度特性。這是一個新問題，一個特殊的彈性力學問題。一個特殊的彈性力學問題。 三角形頂點設為節(jié)點，其三角形頂點設為節(jié)點，其局部編號局部編號為為l，m，n（逆時針）。逆時針）。每節(jié)點有總體坐標每節(jié)點有總體坐標x，y方向兩個待求位移分量：方向兩個待求位移分量：u，v。單元共單元共有有6個位移分量個位移分量6個自由度。個自由度。下面研究有限元法中特有的求解該特下面研究有限元法中特

3、有的求解該特殊彈性力學問題的方法。殊彈性力學問題的方法。 有限元離散結(jié)構(gòu)受力平衡后，取出一個典型三節(jié)點三角形單元有限元離散結(jié)構(gòu)受力平衡后，取出一個典型三節(jié)點三角形單元e。 單元平衡時要在節(jié)點處受到節(jié)點力（節(jié)點對單元的作用單元平衡時要在節(jié)點處受到節(jié)點力（節(jié)點對單元的作用力），每節(jié)點有力），每節(jié)點有2個節(jié)點力分量，單元有個節(jié)點力分量，單元有6個節(jié)點力分量。個節(jié)點力分量。 （2） 單元節(jié)點力列陣單元節(jié)點力列陣 Tnnmmllnmlevuvuvu Tynxnymxmylxlnmlepppppppppp下面要研究的問題是下面要研究的問題是該三角形薄片彈性體在保持平衡時所受節(jié)該三角形薄片彈性體在保持平衡時

4、所受節(jié)點力和節(jié)點位移的關系點力和節(jié)點位移的關系。（1）單元節(jié)點位移列陣）單元節(jié)點位移列陣5.2.2 單元位移模式單元位移模式按彈性力學位移法求近似解的思路，位移作為基本未知量時，需按彈性力學位移法求近似解的思路，位移作為基本未知量時，需要要對單元上位移的分布作出假設對單元上位移的分布作出假設，即構(gòu)造含待定參量的簡單位移即構(gòu)造含待定參量的簡單位移函數(shù)函數(shù)位移模式位移模式。yaxaayxvyaxaayxu654321),(),(為待定系數(shù)，稱為廣義坐標為待定系數(shù)，稱為廣義坐標。61 aa通常用多項式函數(shù)作位移模式，通常用多項式函數(shù)作位移模式，對三節(jié)點三角形單元對三節(jié)點三角形單元，有，有6個待定個待

5、定節(jié)點位移分量，所以單元上的位移函數(shù)只能是節(jié)點位移分量，所以單元上的位移函數(shù)只能是含含6個待定系數(shù)的完個待定系數(shù)的完全一次多項式全一次多項式：3211aaayxu321111aaayxyxyxuuunnmmllnml654111aaayxyxyxvvvnnmmllnml位移多項式寫成矩陣形式位移多項式寫成矩陣形式：代入各節(jié)點取值條件后代入各節(jié)點取值條件后：坐標取節(jié)點值坐標取節(jié)點值6541aaayxv坐標取節(jié)點值坐標取節(jié)點值 由于有限元法中未知量是節(jié)點位移，所以上面單由于有限元法中未知量是節(jié)點位移，所以上面單元位移模式需要轉(zhuǎn)換為以節(jié)點位移分量為待定參量的元位移模式需要轉(zhuǎn)換為以節(jié)點位移分量為待定參

6、量的形式。過程如下形式。過程如下：321111aaayxyxyxuuunnmmllnmlnmlnmlnmlnmlnmlnnmmlluuucccbbbaaauuuyxyxyxaaa211111321其中：其中：nnmmllyxyxyx1112為三角形面積由第一組方程求解由第一組方程求解:31aa節(jié)點坐標行列式節(jié)點坐標行列式分別解出分別解出6個待定系數(shù)個待定系數(shù)：654111aaayxyxyxvvvnnmmllnmlnmlnmlnmlnmlnmlnnmmllvvvcccbbbaaavvvyxyxyxaaa211111654為常數(shù)。元素的代數(shù)余子式，均個，行第，的第分別是節(jié)點坐標行列式321)(,n

7、mlkkcbakkk由第二組方程求解由第二組方程求解:64aa上面求出的待定系數(shù)上面求出的待定系數(shù) 代回位移多項式，得到代回位移多項式，得到：61 aanmlnmlnmlnmluuucccbbbaaayxaaayxu2111321nmliiinnmmllvNvNvNvNv,nmlnnnmmmllluuuycxbaycxbaycxba21nmlnmluuuNNNnmliiinnmmlluNuNuNuN,至此，單元位移模式已轉(zhuǎn)換為節(jié)點位移的插值形式。至此，單元位移模式已轉(zhuǎn)換為節(jié)點位移的插值形式。上式中上式中：)(21ycxbaNiiii），（nmli 位移插值基函數(shù)，稱為位移插值基函數(shù)，稱為形狀函

8、數(shù)（形函數(shù)）形狀函數(shù)（形函數(shù)） enmlnmlNNNNNNNvuvu000000合并后用矩陣表示為：和 N 稱為形函數(shù)矩陣，是對單元節(jié)點位移進行插值得到單元位移稱為形函數(shù)矩陣，是對單元節(jié)點位移進行插值得到單元位移分布函數(shù)的轉(zhuǎn)換矩陣分布函數(shù)的轉(zhuǎn)換矩陣。 用節(jié)點位移插值表示單元位移模式是有限元法中除了離用節(jié)點位移插值表示單元位移模式是有限元法中除了離散化之外最具代表性，最重要的步驟散化之外最具代表性，最重要的步驟!5.2.3 形函數(shù)及其性質(zhì)形函數(shù)及其性質(zhì)對于單元位移模式對于單元位移模式：nnmmllnnmmllvNvNvNvuNuNuNu)(21ycxbaNiiii），（nmli 假設假設：01n

9、mluuu得到得到：lNyxu),(形函數(shù)。對應一個形函數(shù)。單元每個節(jié)點節(jié)點的形狀函數(shù)，簡稱因此稱為的位移分布，節(jié)點固定不動時，單元節(jié)點發(fā)生單位位移，是lnmlNl,顯然，顯然，形函數(shù)決定了單元上位移分布的形態(tài)形函數(shù)決定了單元上位移分布的形態(tài)。事實上，單元位。事實上，單元位移模式就是所有形函數(shù)的線性組合。移模式就是所有形函數(shù)的線性組合。一個單元的位移模式?jīng)Q定了該單元描述局部位移場的能力，決定一個單元的位移模式?jīng)Q定了該單元描述局部位移場的能力，決定求解的精度、收斂性等，而形函數(shù)是最重要的因素。求解的精度、收斂性等，而形函數(shù)是最重要的因素。 針對三節(jié)點三角形單元，可以導出單元形函數(shù)針對三節(jié)點三角形

10、單元，可以導出單元形函數(shù)的下列性質(zhì)。的下列性質(zhì)。性質(zhì)性質(zhì)1：單元上某節(jié)點的形函數(shù)在該節(jié)點的值為單元上某節(jié)點的形函數(shù)在該節(jié)點的值為1，在其它節(jié)點，在其它節(jié)點 的值為零的值為零。性質(zhì)性質(zhì)2：單元上所有形函數(shù)之和等于單元上所有形函數(shù)之和等于1。0),(0),(1),(nnlmmllllyxNyxNyxN)(nml，1nmlNNN（簡單三角形單元的形函數(shù)只有（簡單三角形單元的形函數(shù)只有2個獨立）個獨立）性質(zhì)性質(zhì)3（推論）（推論）：簡單三角形單元的形函數(shù)在邊界上的性質(zhì)。簡單三角形單元的形函數(shù)在邊界上的性質(zhì)。 某節(jié)點的形函數(shù)在該點鄰邊上呈線性分布，取值在某節(jié)點的形函數(shù)在該點鄰邊上呈線性分布，取值在01之間

11、，之間，在該點對邊上值為零在該點對邊上值為零。簡單三角形單元形函數(shù)的幾何意義簡單三角形單元形函數(shù)的幾何意義v 由形函數(shù)表達式和性質(zhì)由形函數(shù)表達式和性質(zhì)1可畫出下列形函數(shù)幾何圖形?？僧嫵鱿铝行魏瘮?shù)幾何圖形。v根據(jù)位移模式表達式及其形函數(shù)的性質(zhì)，可以推斷出兩個相鄰根據(jù)位移模式表達式及其形函數(shù)的性質(zhì)，可以推斷出兩個相鄰三角形單元上位移分布形狀和公共邊界上位移的情況：三角形單元上位移分布形狀和公共邊界上位移的情況：nnmmllnnmmllvNvNvNvuNuNuNu兩個單元上位移線性連續(xù)分布，各單兩個單元上位移線性連續(xù)分布，各單元在公共邊界上位移線性分布，數(shù)值元在公共邊界上位移線性分布，數(shù)值相同相同邊

12、界位移協(xié)調(diào)！邊界位移協(xié)調(diào)！v 由圖形幾何性質(zhì)可以推斷簡單三角形單元形函數(shù)的下列結(jié)論：由圖形幾何性質(zhì)可以推斷簡單三角形單元形函數(shù)的下列結(jié)論：1）三角形形心上三角形形心上：31nmlNNNlmllmlldsNdsNlmlm212）形函數(shù)在邊界上的積分形函數(shù)在邊界上的積分：3）形函數(shù)在單元上的積分形函數(shù)在單元上的積分：3),(AdxdyyxNi),(nmli 5.2.4 單元應變和應力單元應變和應力已知節(jié)點位移插值形式的單元位移模式已知節(jié)點位移插值形式的單元位移模式： eNvu代入平面問題幾何方程（應變代入平面問題幾何方程（應變位移關系）得到單元應變：位移關系）得到單元應變： eenmlennmml

13、lnmlnmlBBBBxNyNxNyNxNyNyNyNyNxNxNxN000000 enmlnmlxyyxNNNNNNxyyxvuxyyx0000000000上式簡寫為上式簡寫為： eB nnmmllnmlnmlbcbcbccccbbbB00000021常應變?nèi)切螁卧?。常?shù)，故該單元又稱為應變在單元上是是常數(shù)矩陣，該單元的所以該單元的應變矩陣關的常數(shù)，均為與單元節(jié)點坐標有、由于上式中),(nmlicbii 為應變矩為應變矩陣，其一個陣，其一個子塊的計算子塊的計算式為式為：B),(00nmlixNyNyNxNBiiiii該式建立了用單元節(jié)點位移表該式建立了用單元節(jié)點位移表達單元上應變分布的關系

14、。達單元上應變分布的關系。對簡單三角形單元，應變矩陣為對簡單三角形單元，應變矩陣為：把單元應變代入平面問題物理方程，即得到單元應力把單元應變代入平面問題物理方程，即得到單元應力： Dxyyx eS 平面問題彈性矩陣單元應力矩陣DBDS eBv應力矩陣的分塊形式應力矩陣的分塊形式： nmlnmlSSSBBBDS iiBDS ),(nmli eeSBDv對于平面應力問題，應力矩陣的子塊為： iiiiiiibccbcbES2121)1 (22),(nmli v由于均質(zhì)材料的彈性系數(shù)為常數(shù)，應力矩陣也是常數(shù)矩陣，即該單元上應力分布也是常數(shù)。對平面應變問題的應力矩陣，只要把上式中彈性系數(shù)作相應變換：21

15、EEu1u 平面問題簡單三角形單元應力的討論平面問題簡單三角形單元應力的討論單元上應力和應變均為常數(shù)單元上應力和應變均為常數(shù)，其大小與單元幾何和節(jié)點位移，其大小與單元幾何和節(jié)點位移有關，一般各單元應力值不相同，應力在單元之間不連續(xù)，有關，一般各單元應力值不相同，應力在單元之間不連續(xù)，這是這是有限元解的近似性的反映有限元解的近似性的反映。一般把單元上這個常應力值。一般把單元上這個常應力值作為單元中心的應力值較合理，有較高的精度。作為單元中心的應力值較合理，有較高的精度。理論和數(shù)值試驗均可證明，用該單元求解問題時，誤差隨單元理論和數(shù)值試驗均可證明，用該單元求解問題時，誤差隨單元尺寸的減小和單元數(shù)目

16、的增加而減小，這就是尺寸的減小和單元數(shù)目的增加而減小，這就是有限元解的收斂有限元解的收斂性性。事實上，當單元變得越來越小時，結(jié)構(gòu)在一個小單元區(qū)域。事實上，當單元變得越來越小時，結(jié)構(gòu)在一個小單元區(qū)域上的應力趨于均勻，而大量小單元拼接在一起就可能很精確地上的應力趨于均勻，而大量小單元拼接在一起就可能很精確地模擬結(jié)構(gòu)中變化的應力。模擬結(jié)構(gòu)中變化的應力。5.2.5 單元剛度方程與剛度矩陣單元剛度方程與剛度矩陣前面已經(jīng)得到了用單元節(jié)點前面已經(jīng)得到了用單元節(jié)點位移表示的單元內(nèi)部位移、位移表示的單元內(nèi)部位移、應變、應力的分布：應變、應力的分布： Tnnmmllevuvuvu Tynxnymxmylxlepp

17、ppppp eNvu eBD eB結(jié)構(gòu)受力平衡時，每個單元在節(jié)點力結(jié)構(gòu)受力平衡時，每個單元在節(jié)點力 作作用下也保持平衡，節(jié)點產(chǎn)生相應位移用下也保持平衡，節(jié)點產(chǎn)生相應位移 。 e ep下面下面用彈性體虛位移原理（虛功原理）建立單元的平衡關系用彈性體虛位移原理（虛功原理）建立單元的平衡關系。v設單元節(jié)點虛位移：設單元節(jié)點虛位移： Tnnmmllevuvuvu*則相應有單元內(nèi)部虛位移則相應有單元內(nèi)部虛位移： eNvu*相應單元虛應變?yōu)椋合鄳獑卧搼優(yōu)椋?eB*v外力虛功計算外力虛功計算： eTeynnxnnymmxmmyllxllppvpupvpupvpuW)(*對一個單元而言，外力即單元節(jié)點力對

18、一個單元而言，外力即單元節(jié)點力。v按虛功原理，虛位移下單元的外力虛功等于虛應變能按虛功原理，虛位移下單元的外力虛功等于虛應變能。v虛應變能計算虛應變能計算單元虛應變能即單元內(nèi)應力在虛應變上所作虛功。單元虛應變能即單元內(nèi)應力在虛應變上所作虛功。 dVdVUTVxyxyyyxVxee*)( eB* dVBDBdVUeTVeVTee)(* eBD eVeTTedVBDB)(* eVTTeedVBDB)(*v由虛功原理得到由虛功原理得到：根據(jù)上面推導，有：根據(jù)上面推導，有： eTepW)(* eVTTeedVBDBU)(*外力虛功外力虛功：虛應變能虛應變能：UW eVTTeeTeedVBDBp)()(

19、*考慮到虛位移考慮到虛位移 的任意性，由上式立即得到下列方程：的任意性，由上式立即得到下列方程：e* eVTeedVBDBp上式簡寫成： eeekp其中單元剛度矩陣： eVTedVBDBk該方程建立了單元節(jié)點力和節(jié)點位移之間的關系，為單元剛度方程。至此，我們已經(jīng)得到了關于彈性實體小單元的力學特性方程。由于彈性矩陣 是對稱的，由上式可以看出，單元剛度矩陣 是對稱方陣。對3節(jié)點三角形單元，單元剛度矩陣是66 方陣。 D ek上面單元剛度方程和單元剛度矩陣的推導雖然在平面問題的簡單三角形單元下得到，但推導的過程和公式的形式具有一般性，推導原理和結(jié)果適用一般的連續(xù)體力學問題的所有單元。提提 示：示：單

20、元剛度方程和單元剛度矩陣的建立是單元分析的核心內(nèi)容。一般情況下，單元應變矩陣是坐標的函數(shù)矩陣，所以單元剛度矩陣的計算需要進行積分運算。所建立的單元剛度矩陣反映了一般彈性體小單元近似的彈性性質(zhì)，是單元特性的核心。 eeekp單元剛度矩陣的計算單元剛度矩陣的計算彈性力學平面問題的單元剛度矩陣通式彈性力學平面問題的單元剛度矩陣通式： eVTedVBDBk（單元剛度矩陣通式） ehdxdyBDBkTe（平面問題單剛通式）單元厚度（常數(shù)）h由于應變矩陣是常數(shù)矩陣，所以由平面問題剛度矩陣通式得到：平面問題簡單三角形單元剛度矩陣計算 nmlTnTmTlTeBBBDBBBhhBDBk nnnmnlmnmmml

21、lnlmllekkkkkkkkkk33分塊形式對平面應力問題，上述剛度矩陣的一個22子矩陣為： srsrsrsrsrsrsrsrsTrersbbcccbbcbccbccbbEhhBDBk21212121)1 (42）（nmlsr, 單元剛度方程的分塊形式 eeekp nnnmnlmnmmmllnlmllekkkkkkkkkk nmle nmlepppp簡單三角形單元剛度方程按節(jié)點分塊形式：nmlnnnmnlmnmmmllnlmllnmlkkkkkkkkkppp展開的單元剛度方程： ),(nmlikkkpninmimlili每行各節(jié)點剛度子塊反映了節(jié)點位移對節(jié)點力的貢獻。5.2.6 簡單三角形單

22、元分析計算實例簡單三角形單元分析計算實例(1)(1) 形函數(shù)和形函數(shù)矩陣前面針對一維桿單元和平面問題簡單三角形單元引入的單元形函數(shù)，應變、應力矩陣，單元剛度矩陣等是結(jié)構(gòu)有限元分析中最基本、最重要的概念。dd厚度h下面通過分析計算一個具體的單元加以說明。按公式計算各形函數(shù)按公式計算各形函數(shù))(21ycxbaNiiii），（nmli 200101011112dddyxyxyxnnmmlldddddcccbbbaaanmlnmlnml00002代數(shù)余子式矩陣：按形函數(shù)公式，得到三個形函數(shù)：計算節(jié)點坐標行列式：dd厚度h容易驗證其滿足形函數(shù)性質(zhì)dydxNdyNdxNnml1按形函數(shù)性質(zhì)和幾何意義，直接

23、寫出形函數(shù)表達式按形函數(shù)性質(zhì)和幾何意義，直接寫出形函數(shù)表達式dd厚度h形函數(shù) 的圖形是一個平面，其方程是：lNdxNl上述兩個函數(shù)定義在三角形區(qū)域的部分就是形函數(shù)。同理，直接寫出形函數(shù) 圖形的方程為：mNdyNm由形函數(shù)的性質(zhì)，得到第三個形函數(shù)：dydxNNNmln11顯然，與公式計算得到的結(jié)果相同！3）形函數(shù)矩陣形函數(shù)矩陣根據(jù)計算所得單元形函數(shù)，寫出單元的形函數(shù)矩陣： dydxdydxdydxdydxN10000100(2) (2) 應變矩陣應變矩陣由公式計算單元應變矩陣： 1101101010000100011dB),(00nmlixNyNyNxNBiiiii nnmmllnmlnmlbcbcbccccbbbB00000021dddddcccbbbaaanmlnmlnml00002(3) (3) 應力矩陣應力矩陣 212102121011001001)1 (2dESdddddcccbbbaaanmlnmlnml00002 iiiiiiibccbcbES2121)1 (22),(nmli 對平面應力問題：(4) (4) 單元剛度矩陣計算單元剛度矩陣計算dddddcccbbbaaanmlnmlnml00002 srsrsrsrsrsrsrsrsTrersbbcccbbcbccbccbbEhhBDBk21212121)1 ( 42）（nmlsr, 平