第十一章節(jié)平面問題有限元分析及三角形單元應(yīng)用

1、第十一講 平面問題的有限元分析及三 角形單元的應(yīng)用 第一節(jié)第一節(jié) 概述概述 第二節(jié)第二節(jié) 單元分析單元分析 第三節(jié)第三節(jié) 等效結(jié)點荷載等效結(jié)點荷載 第四節(jié)第四節(jié) 整體剛度矩陣整體剛度矩陣 第五節(jié)第五節(jié) 平面問題分析舉例平面問題分析舉例 第六節(jié)第六節(jié) 單元網(wǎng)格的劃分和計算成果的整理單元網(wǎng)格的劃分和計算成果的整理 第一節(jié) 概述 分析彈性力學(xué)平面問題時，最簡單的單元式由三個結(jié)點組成的三 角形單元。當(dāng)用以分析平面應(yīng)力問題時，可將其視為三角板；當(dāng)用以 分析平面應(yīng)變問題時，則可式為三棱柱。各單元在結(jié)點處為鉸結(jié)。圖 8-1 所示位移懸臂梁離散為三角形單元的組合體 以矩陣形式列出彈性力學(xué)平面問題的基本量和基本

2、方程。 談形體所受體力分量可表示為 T yx y x pp p p p (8-1) 所受面力分量可表示為 T yx y x pp p p p （8-2） 體內(nèi)任一點應(yīng)力分量可表示為 T xyyx （8-3） 任一點的應(yīng)變分量可表示為 T xyyx （8-4） 任一點的位移分量可表示為 T vu （8-5） 彈性力學(xué)平面問題的幾何方程的矩陣表達(dá)式為 x u y v y v x u xy y x （8-6） 平面應(yīng)力問題的物理方程的矩陣表達(dá)式為 xy y x xy y x E 2 1 00 01 01 1 2 （8-7） 或簡寫成為 D （8-8） 式中 2 1 00 01 01 1 2 E D

3、（8-9） 稱為彈性矩陣。 平面應(yīng)變問題的物理方程也可寫成式（8-8） ，但須將式（8-9）中 的 E 換成 2 1 E ，換成 2 1 ，因此得出 )1 (2 21 00 01 1 0 1 1 )21)(1 ( )1 ( 2 2 E D （8-10） 平衡微分方程及邊界條件也可以用矩陣表示， 但彈性力學(xué)有限元位移法 中，通常用虛功方程代替平衡微分方程和應(yīng)力邊界條件。虛功方程的矩陣表 達(dá)式為 tdxdytdspfptdxdyf TT （8-11） 式中： T vuf ，表示虛位移； T xyxx ，表示與虛位移相對應(yīng)的虛應(yīng)變。 為了便于計算， 作用于彈性體上的體力和面力替換為作用在結(jié)點上的集

4、中力，即等效結(jié)點荷載。設(shè)作用于各個結(jié)點上的外力分量用如下列陣來表示 T nn VUVUVUF 2211 與這些結(jié)點外力分量相對應(yīng)得結(jié)點虛位移分量列陣為 T nn vuvuvu 2211 則外力在虛位移上做的虛功為 FvVuUvVuUvVuU T nnnn 22221111 如平面彈性體的厚度為 t，該虛功除以 t，即可得出單位厚度薄 板上的外力虛功。于是，式（8-11）所示虛功方程可寫成 tdxdyF TT （8-11） 虛功方程不僅僅應(yīng)用于彈性力學(xué)，也可用于塑性力學(xué)。其應(yīng)用 條件是：只要變形體的全部外力和應(yīng)力滿足平衡方程；位移是微小 的，并滿足邊界條件，位移與應(yīng)變滿足幾何方程。所以，通常稱為

5、 變形體虛功方程。 第二節(jié) 單元分析 圖 8-2 所示為一個三角形單元。三個結(jié)點按逆時針順序編號分別為 i、 j、m，結(jié)點坐標(biāo)分別為 ),(),(),( mmjjii yxmyxjyxi、 。 圖 8-2 由 于 每 個 結(jié) 點 有 兩 個 位 移 分 量 ， 單 元 共 有 六 個 結(jié) 點 位 移 分 量 ： mmjjii vuvuvu、 ，如圖 8-2a）所示，因此三角形單元的結(jié)點位移分量 e 可表示為 T mmjjii e vuvuvu （8-13） 與這六個結(jié)點位移分量相對應(yīng)得結(jié)點力也有六個分量，如圖 8-2b)所示 T mmjjii e VUVUVUF （8-14） 在每個單元上，都

6、可以把結(jié)點力用結(jié)點位移來表示，即建立 如下關(guān)系式 eee kF （8-15） 式中 k e 稱為單元剛度矩陣。尋求 k e 的過程稱為單元分析。單 元分析按如下步驟 一 、 位 移 函 數(shù) 為 了 求 單 元 內(nèi) 任 一 點 （ x， y） 的 位 移 ， 設(shè) 該 點 的 位 移 u 、 v 為 其 坐 標(biāo) x、 y 的 某 種 函 數(shù) ， 單 元 有 六 個 結(jié) 點 位 移 分 量 ， 在 位 移 函 數(shù) 中 取 六 個 任 意 參 數(shù) i（ i=1 ， 2， ， 6） ， 并 將 位 移 函 數(shù) 取 為 線 性 函 數(shù) ， 即 yxyxv yxyxu 654 321 ),( ),( （ 8

7、-16） 一 般 情 況 下 ， 一 個 彈 性 變 形 體 在 外 界 作 用 下 ， 內(nèi) 部 點 的 位 移 變 化 比 較 復(fù) 雜 ， 不 能 用 簡 單 的 線 性 函 數(shù) 描 述 。 但 是 ， 當(dāng) 把 彈 性 體 離 散 為 許 多 微 笑 單 元 時 ， 在 每 一 個 單 元 內(nèi) 部 有 限 小 的 局 部 內(nèi) ，各 點 位 移 可 以 用 線 性 函 數(shù) 描 述 。式（ 8-16）可 寫 成 矩 陣 形 式 6 5 4 3 2 1 1000 0001 yx yx v u f （8-17） 為了求出內(nèi)部結(jié)點位移 f 與結(jié)點位移 e 之間的關(guān)系，需求出 e 與間的關(guān) 系。降格結(jié)點

8、坐標(biāo)和位移代入式（8-16） ，可得 3 2 1 1 1 1 mm jj ii m j i yx yx yx u u u 6 5 4 1 1 1 mm jj ii m j i yx yx yx v v v (a) （b） 三角形單元的面積為 mm jj ii yx yx yx A 1 1 1 2 1 （8-18） 求解方程組(a)得 m j i mji mji mji u u u ccc bbb aaa A2 1 3 2 1 （c） 求解方程組(b)得 m j i mji mji mji v v v ccc bbb aaa A2 1 6 5 4 (d) 式中， mmmiii cbacba、 由

9、下式計算 mji mji jmmji xxc yyb yxyxa (i、j、m) 上式中的(i、j、m)表示腳標(biāo)依次輪換，可寫出計算 aj、bj、cj以及 am、bm、cm的 另兩組公式。將式(c)和(d)代入(8-16)并展開，得到以結(jié)點位移表示的位移函數(shù) j i j i j i mji mji v u v u v u yxNyxNyxN yxNyxNyxN yxv yxu ),(0),(0),(0 0),(0),(0),( ),( ),( （8-20） 式中， mji NNN、 反映了單元的位移形態(tài)， 故稱為單元位移的形態(tài)函數(shù)或形函數(shù)。 矩陣 N 稱為形函數(shù)矩陣。 選取得位移函數(shù)是否合理，

10、要看隨著單元網(wǎng)格的逐步細(xì)分，有限元解是否 逼近于精確解。為了保證收斂型所選擇的單元位移函數(shù)應(yīng)滿足以下條件： （1） 包含單元的剛體位移； （2） 包含單元的常量應(yīng)變； （3） 保證相鄰單元在公共邊界處位移的連續(xù)性。 一、 單元的應(yīng)變和應(yīng)力 選擇了位移函數(shù)并以結(jié)點位移表示單元內(nèi)點的位移后， 重新寫出平面問 題的幾何方程 v u xy y x xy y x 0 0 （f） 由式（8-20）得 mmjjii mmijii vNvNvNv uNuNuNu （g） 將式（g）代入式（f） ，并利用下式 A c y N A b x N ii ii 2 2 （i、j、m） (h) 得單元應(yīng)變 m m j j

11、 i i mmjjii mji mji xy y x v u v u v u bcbcbc ccc bbb A 000 000 2 1 （8-24） 或簡寫成 e B (8-25) 式中 mmjjii mji mji bcbcbc ccc bbb A B000 000 2 1 (8-26) 式（8-25）就是由結(jié)點位移求應(yīng)變的轉(zhuǎn)換式，其轉(zhuǎn)換矩陣 B 稱為幾何矩陣。 將（8-25）代入平面問題的物理方程式（8-8）有 ee DBD (8-27) 或?qū)懗?ee S (8-28) S=DB (8-29) 稱為應(yīng)力矩陣。 二、 單元剛度矩陣 在有限單元法中，常利用虛功方程代替平衡方程。圖 8-3a)所

12、示為 三角形單元的實際力系，其結(jié)點力為 F e，應(yīng)力為；圖 8-3b)所示 為單元虛位移狀態(tài)，其結(jié)點位移為 *e,應(yīng)變?yōu)?。利用式（8-12） 可得 圖 8-3 tdxdyF TeeT （8-30） 式中： 為單元虛應(yīng)變。為單元結(jié)點虛位移； e 由式（8-25）可知 e B 因此 TeTT B 將此公式代入式（8-30） ，由于 e 中的元素是常量，公式右邊的 eT 可 以提到積分號的前面，得 tdxdyBF TeTeeT 由于虛位移 e 是任意的，則 tdxdyBF e 因為 B 和都是常量矩陣，并且積分 Adxdy ，所以 tABF Te （8-31） 利用（8-27） ，可得 tADBB

13、F eTe （8-32） 令 DBtABk T （8-33） 則式（8-32）就變成式（8-15） ，即 eee kF 單元剛度矩陣 e k 為一個 66 矩陣，它時單元結(jié)點位移與單元結(jié)點力之 間的轉(zhuǎn)換矩陣，具有以下性質(zhì)： （1） e k 示對稱矩陣，其元素jiij kk ； （2） e k 是奇異矩陣，由它的元素組成的行列式等于零，即它不存在逆 矩陣； （3） e k 具有分快性質(zhì)。 第三節(jié) 等效結(jié)點荷載 為簡化各單元得受力情況，便于分析計算，應(yīng)將單元所受各種載荷向 結(jié)點移置，化為結(jié)點荷載，荷載的移置應(yīng)安靜力等效原則進(jìn)行。靜力等效 是令原來的荷載與移置后的荷載在任意虛位移上的虛功相等。 一、

14、集中荷載 圖 8-5a）所示為單元內(nèi)的任意一點 M 受到集中荷載 P 的作用，沿 x、 y 方向的分量分別為 Px、Py，用矩陣表示為 T yx PPP 。設(shè)移置到該單元結(jié) 點上的等效荷載列陣為 T mmjjii e YXYXYXR 圖 8-5 所示單元發(fā)生虛位移，其中單元內(nèi)部任意點的虛位移為 T vuf 設(shè)單元各結(jié)點虛位移為 T mmjjii e vuvuvu 由式（8-22） ，則 e Nf （a） 圖 8-5 如果設(shè) Re為單元得等效結(jié)點荷載，則 Re在結(jié)點虛位移上所做的 虛功應(yīng)與原來集中荷載在其作用點的虛位移上做的虛功相等，即 PfR TeeT 將式(a)代如上式，得 PNR TeTe

15、eT 由于虛位移可以式任意的，因此 PNR Te （8-36） 或?qū)懗?yii xii PNY PNX （i、j、m） 二、 分布體力 設(shè)單元受分布體力的作用 T yx ppp 將微分體積 tdxdy 上的體力當(dāng)作集中荷載 P，利用式（8-36）的 積分得出 ptdxdyNR Te （8-38） 如果單元上作用的分布體力為自重，為材料的重度，即 Tp 0 那么利用式（8-38）得到等效結(jié)點荷載列陣為 T e At R101010 3 （8-39） 三、分布面力 設(shè)單元在邊界 ij 上首有分布面力的作用 T yx ppp 將微分面積 tds 上的面力當(dāng)作集中荷載 P，利用式（3-36）的積分 可

16、以得出 tdspNR Te （8-40） 下面給出作用于單元邊界上荷載的兩種簡單情況。如圖8-6所示的 邊界ij上沿x方向作用一集中力 ，其作用點距i，j兩點的距離 分別為 和 。其等效結(jié)點和在列陣為 （8-41） 0 P i l j l 0 0000 T je i l l RP ll 如圖8-7所示單元的邊界ij上受到沿x方向作用的按三角形分布的 荷載，在i點的荷載集度為q,則其等效結(jié)點荷載列陣為 （8-42） 如上述兩種荷載為任意方向，可分別按其在x方向和y方向的分量進(jìn) 行計算。 21 0000 233 T e qtl R 如單元上同時作用有集中力，分布體力和分布面力，則結(jié) 點荷載應(yīng)為三部

17、分荷載之和，即： （8-43） eTTT Al RN pN ptdxdyN ptds 第四節(jié) 整體剛度矩陣 有限元求解彈性力學(xué)問題也要用結(jié)點平衡方程求解作為基本為 質(zhì)量的整體結(jié)點位移列陣。在求出各單元的單元剛度矩陣 ke和解點荷 載列陣后，就可以用幾何的方法建立其平面彈性體的整體剛度矩陣和 整體平衡方程。如果彈性體劃分為 m 各單元，n 各結(jié)點，則有 121222 nnnn RK （8-44） 式中：整體剛度矩陣 K 為 2n2n 階方程；、R 分別為整體結(jié) 點位移列陣和整體結(jié)點荷載列陣，都是 2n1 階列陣。 用集成剛度法將e單元的單元貢獻(xiàn)矩陣表示如下, （8-45） 1 m e e KK

18、下舉例說明用單元剛度矩陣形成整體剛度矩陣的方法。圖8-8 表示一彈性體劃分為三個單元，并對結(jié)點進(jìn)行局部編號，整體 編號為1，2，3，4，55單元結(jié)點局部編號為I,j,m按逆時針方向 排序。整體剛度矩陣為10階方陣。 首先將單元e的剛度矩陣 擴大為10階方陣，成為單位貢 獻(xiàn)矩陣。單元局部編號應(yīng)與其所在位置的整體編號相對應(yīng)。以 單元3為例，其結(jié)點局部編號I,j,m與整體編號2，4，5對應(yīng)，形 成單元3的貢獻(xiàn)矩陣。 整體號 1 2 3 4 5 i j m 局部號 e k 333 3 333 333 1 2 3 4 5 iiijim jijjjm mimjmm Kkki K kkkj kkkm 然后將各個單元貢獻(xiàn)矩陣疊加起來，形成整體剛度矩陣。 在實際編程中，為節(jié)省容量，將各個單元剛度矩陣 中的子 矩陣 逐個搬到整體剛度矩陣對應(yīng)的位置上 去，建立整體矩陣如下： e k 2 2( , , ,) rs kr si j m 111 1 1 2 123 2 3 3 3 1 11 12 1 22 333 32 23 33 1 2 3 4 5 iiijim jj