第05章 有限元法-1

第5章有限元法(1)

ⅤFiniteElementMethod

(FEM)第5章有限元法內(nèi)容簡(jiǎn)介

有限元法是結(jié)構(gòu)分析的一種數(shù)值計(jì)算方法。它在20世紀(jì)50年代初期隨著計(jì)算機(jī)的發(fā)展應(yīng)運(yùn)而生。

這一方法的理論基礎(chǔ)牢靠，物理概念清晰，解題效率高，適應(yīng)性強(qiáng)，目前已成為機(jī)械產(chǎn)品動(dòng)、靜、熱特性分析的重要手段，它的程序包是機(jī)械產(chǎn)品計(jì)算機(jī)輔助設(shè)計(jì)方法庫(kù)中不可缺少的內(nèi)容之一。

本章介紹了如下內(nèi)容:

■

有限元法的基本思想及應(yīng)用

■

平面問(wèn)題有限元分析原理及步驟

■

有限元分析的前后處理

■

有限元法的設(shè)計(jì)應(yīng)用及計(jì)算實(shí)例5.1概述

在工程分析和科學(xué)研究中，常常會(huì)遇到大量的由常微分方程、偏微分方程及相應(yīng)的邊界條件描述的場(chǎng)問(wèn)題，如位移場(chǎng)、應(yīng)力場(chǎng)和溫度場(chǎng)等問(wèn)題。目前求解這類場(chǎng)問(wèn)題的方法主要有兩種：

●

用解析法求得精確解；

●

用數(shù)值解法求其近似解。其中，

能用解析法求出精確解的只能是方程性質(zhì)比較簡(jiǎn)單且?guī)缀芜吔缦喈?dāng)規(guī)則的少數(shù)問(wèn)題。

而對(duì)于絕大多數(shù)問(wèn)題，則很少能得出解析解。這就需要研究它的數(shù)值解法，以求出近似解。

目前，工程中實(shí)用的數(shù)值解法主要有三種：

有限差分法

有限元法

邊界元法其中，以有限元法通用性最好，解題效率高，工程應(yīng)用最廣。目前它已成為機(jī)械產(chǎn)品動(dòng)、靜、熱特性分析的重要手段，它的程序包是機(jī)械產(chǎn)品計(jì)算機(jī)輔助設(shè)計(jì)方法庫(kù)中不可缺少的內(nèi)容之一。“有限元法

”的基本思想早在20世紀(jì)40年代初期就有人提出，但真正用于工程中則是電子計(jì)算機(jī)出現(xiàn)以后?！坝邢拊?/p>

”這一名稱是1960年美國(guó)的克拉夫（Clough，R.W.）在一篇題為“平面應(yīng)力分析的有限元法”論文中首先使用。此后，有限元法的應(yīng)用得到蓬勃發(fā)展。到20世紀(jì)80年代初期國(guó)際上較大型的結(jié)構(gòu)分析有限元通用程序多達(dá)幾百種，從而為工程應(yīng)用提供了方便條件。由于有限元通用程序使用方便，計(jì)算精度高，其計(jì)算結(jié)果已成為各類工業(yè)產(chǎn)品設(shè)計(jì)和性能分析的可靠依據(jù)。

有限元法的分析過(guò)程可概括如下：●

連續(xù)體離散化●

單元分析●

整體分析●

確定約束條件●

有限元方程求解●

結(jié)果分析與討論

1.連續(xù)體離散化

連續(xù)體：是指所求解的對(duì)象（如物體或結(jié)構(gòu)）。

離散化（劃分網(wǎng)格或網(wǎng)絡(luò)化）：是將所求解的對(duì)象劃分為有限個(gè)具有規(guī)則形狀的微小塊體，把每個(gè)微小塊體稱為單元，相鄰兩個(gè)單元之間只通過(guò)若干點(diǎn)互相連接，每個(gè)連接點(diǎn)稱為節(jié)點(diǎn)。相鄰單元只在節(jié)點(diǎn)處連接，載荷也只通過(guò)節(jié)點(diǎn)在各單元之間傳遞，這些有限個(gè)單元的集合體，即原來(lái)的連續(xù)體。

單元?jiǎng)澐趾螅o每個(gè)單元及節(jié)點(diǎn)進(jìn)行編號(hào)；

選定坐標(biāo)系，計(jì)算各個(gè)節(jié)點(diǎn)坐標(biāo)；

確定各個(gè)單元的形態(tài)和性態(tài)參數(shù)以及邊界條件等。

圖5-1所示為將一懸臂梁建立有限元分析模型的例子。圖中將該懸臂梁劃分為許多三角形單元；

三角形單元的三個(gè)頂點(diǎn)都是節(jié)點(diǎn)。圖5-1懸臂梁及其有限元模型

2.單元分析

連續(xù)體離散化后，即可對(duì)單元體進(jìn)行特性分析，簡(jiǎn)稱為單元分析。

單元分析工作主要有兩項(xiàng)：

(1)選擇單元位移模式(位移函數(shù))

用節(jié)點(diǎn)位移來(lái)表示單元體內(nèi)任一點(diǎn)的位移、應(yīng)變和應(yīng)力，就需搞清各單元中的位移分布。

一般是假定單元位移是坐標(biāo)的某種簡(jiǎn)單函數(shù)，用其模擬內(nèi)位移的分布規(guī)律，這種函數(shù)就稱為位移模式或位移函數(shù)。通常采用的函數(shù)形式多為多項(xiàng)式。

根據(jù)所選定的位移模式，就可以導(dǎo)出用節(jié)點(diǎn)位移來(lái)表示單元體內(nèi)任一點(diǎn)位移的關(guān)系式。

進(jìn)行單元力學(xué)特性分析，將作用在單元上的所有力（表面力、體積力、集中力）等效地移置為節(jié)點(diǎn)載荷；

采用有關(guān)的力學(xué)原理建立單元的平衡方程，求得單元內(nèi)節(jié)點(diǎn)位移與節(jié)點(diǎn)力之間的關(guān)系矩陣單元?jiǎng)偠染仃嚒?/p>

(2)

分析單元的特性，建立單元?jiǎng)偠染仃?/p>

3.

整體分析

把各個(gè)單元的剛度矩陣集成為總體剛度矩陣，以及將各單元的節(jié)點(diǎn)力向量集成總的力向量，求得整體平衡方程。

集成過(guò)程所依據(jù)的原理是節(jié)點(diǎn)變形協(xié)調(diào)條件和平衡條件。

4.

確定約束條件

由上述所形成的整體平衡方程是一組線性代數(shù)方程，在求解之前，必修根據(jù)具體情況分析，確定求解對(duì)象問(wèn)題的邊界約束條件，并對(duì)這些方程進(jìn)行適當(dāng)修正。5.

有限元方程求解通過(guò)求解整體平衡方程，即可求得各節(jié)點(diǎn)的位移，進(jìn)而根據(jù)位移可計(jì)算單元的應(yīng)力及應(yīng)變。6.

結(jié)果分析與討論

應(yīng)用有限元法求解機(jī)械結(jié)構(gòu)應(yīng)力類問(wèn)題時(shí)，根據(jù)未知量和分析方法的不同，有三種基本解法：

位移法

力法

混合法

(1)位移法

此法是以節(jié)點(diǎn)位移作為基本未知量，通過(guò)選擇適當(dāng)?shù)奈灰坪瘮?shù)，進(jìn)行單元的力學(xué)特性分析。在節(jié)點(diǎn)處建立單元?jiǎng)偠确匠?，再組合成整體剛度矩陣，求解出節(jié)點(diǎn)位移后，進(jìn)而由節(jié)點(diǎn)位移求解出應(yīng)力。

位移法優(yōu)點(diǎn)是比較簡(jiǎn)單，規(guī)律性強(qiáng)，易于編寫(xiě)計(jì)算機(jī)程序。所以得到廣泛應(yīng)用，其缺點(diǎn)是精度稍低。

(2)力法

該法是以節(jié)點(diǎn)力作為基本未知量，在節(jié)點(diǎn)處建立位移連續(xù)方程，求解出節(jié)點(diǎn)力后，再求解節(jié)點(diǎn)位移和單元應(yīng)力。

力法的特點(diǎn)是計(jì)算精度高。

(3)混合法

此法是取一部分節(jié)點(diǎn)位移和一部分節(jié)點(diǎn)力作為基本未知量，建立平衡方程進(jìn)行求解。

有限元法最初用于飛機(jī)結(jié)構(gòu)的強(qiáng)度設(shè)計(jì)，由于它在理論上的通用性，因而它可用于解決工程中的許多問(wèn)題。目前，它可以解決幾乎所有的連續(xù)介質(zhì)和場(chǎng)的問(wèn)題，包括熱傳導(dǎo)、電磁場(chǎng)、流體動(dòng)力學(xué)、地質(zhì)力學(xué)、原子工程和生物醫(yī)學(xué)等方面的問(wèn)題。

機(jī)械設(shè)計(jì)中，從齒輪、軸、軸承等通用零部件到機(jī)床、汽車、飛機(jī)等復(fù)雜結(jié)構(gòu)的應(yīng)力和變形分析（包括熱應(yīng)力和熱變形分析）。

有限元法不僅可以解決工程中的線性問(wèn)題、非線性問(wèn)題，而且對(duì)于各種不同性質(zhì)的固體材料，如各向同性和各向異性材料，粘彈性和粘塑性材料以及流體均能求解；對(duì)于工程中最有普遍意義的非穩(wěn)態(tài)問(wèn)題也能求解。

到20世紀(jì)80年代初期，國(guó)際上已開(kāi)發(fā)出了多種用于結(jié)構(gòu)分析的有限元通用程序，其中著名的有NASTRAN、ANSYS、ASKA、ADINA、SAP等。

表5-1列出了幾種國(guó)際上流行的商用有限元程序的應(yīng)用范圍。

表5-1

幾種有限元程序及其應(yīng)用范圍

程序名應(yīng)用范圍ADINAANSYSASKAMARCNASTRANSAP非線性分析√√√√√塑性分析√√√√斷裂力學(xué)√熱應(yīng)力與蠕變√√√√√厚板厚殼√√√√√√管道系統(tǒng)√√√√船舶結(jié)構(gòu)√√√√√焊接接頭√粘彈性材料√√表5-1(續(xù))幾種有限元程序及其應(yīng)用范圍

程序名應(yīng)用范圍ADINAANSYSASKAMARCNASTRANSAP熱傳導(dǎo)√√√√√√薄板薄殼√√√√√√復(fù)合材料√結(jié)構(gòu)穩(wěn)定性√√流體力學(xué)√√氣動(dòng)彈性力學(xué)√電場(chǎng)√5.2

單元特性的推導(dǎo)方法5.2.1單元?jiǎng)澐址椒霸瓌t

連續(xù)體離散化時(shí)，要根據(jù)設(shè)計(jì)對(duì)象的具體情況（結(jié)構(gòu)物的形狀、載荷特性、邊界條件等），確定單元（網(wǎng)格）的大小和形狀、單元的數(shù)目以及劃分方案。

圖5-2桿狀單元圖5-2所示為桿狀單元?？捎幸痪S、二維和三維梁?jiǎn)卧D5-3

平面單元和多面體單元常用的平面單元和多面體單元，如圖5-3所示。

單元的劃分基本上是任意的，一個(gè)結(jié)構(gòu)體可以有多種劃分結(jié)果。但應(yīng)遵循以下劃分原則：

(1)分析清楚所討論對(duì)象的性質(zhì)，例如，是桁架結(jié)構(gòu)還是結(jié)構(gòu)物，是平面問(wèn)題還是空間問(wèn)題等等。

(2)單元的幾何形狀取決于結(jié)構(gòu)特點(diǎn)和受力情況，單元的幾何尺寸(大小)要按照要求確定。一般來(lái)說(shuō)，單元幾何形體各邊的長(zhǎng)度比不能相差太大。

(3)有限元模型的網(wǎng)格劃分越密，其計(jì)算結(jié)果越精確，但計(jì)算工作量就越大。因此，在保證計(jì)算精度的前提下，單元網(wǎng)格數(shù)量應(yīng)盡量少。

(4)在進(jìn)行網(wǎng)格疏密布局時(shí)，應(yīng)力集中或變形較大的部位，單元網(wǎng)格應(yīng)取小一些，網(wǎng)格應(yīng)劃分得密一些，而其他部分則可疏一些。

(5)在設(shè)計(jì)對(duì)象的厚度或者彈性系數(shù)有突變的情況下，應(yīng)該取相應(yīng)的突變線作為網(wǎng)格的邊界線；

(6)相鄰單元的邊界必須相容，不能從一單元的邊或者面的內(nèi)部產(chǎn)生另一個(gè)單元的頂點(diǎn)。

(7)網(wǎng)格劃分后，要將全部單元和節(jié)點(diǎn)按順序編號(hào)，不允許有錯(cuò)漏或者重復(fù)；

(8)劃分的單元集合成整體后，應(yīng)精確逼近原設(shè)計(jì)對(duì)象。原設(shè)計(jì)對(duì)象的各個(gè)頂點(diǎn)都應(yīng)該取成單元的頂點(diǎn)。所有網(wǎng)格的表面頂點(diǎn)都應(yīng)該在原設(shè)計(jì)對(duì)象的表面上。所有原設(shè)計(jì)對(duì)象的邊和面都應(yīng)被單元的邊和面所逼近。

5.2.2

單元特性的推導(dǎo)方法

單元?jiǎng)偠染仃嚨耐茖?dǎo)是有限元分析的基本步驟之一。目前，建立單元?jiǎng)偠染仃嚨姆椒ㄖ饕幸韵滤姆N：直接剛度法虛功原理法能量變分法加權(quán)殘數(shù)法

1.直接剛度法

直接剛度法是直接應(yīng)用物理概念來(lái)建立單元的有限元方程和分析單元特性的一種方法。這一方法僅能適用于簡(jiǎn)單形狀的單元，如梁?jiǎn)卧５梢詭椭斫庥邢拊ǖ奈锢砀拍睢?/p>

圖5-4所示是xoy平面中的一簡(jiǎn)支梁簡(jiǎn)圖，現(xiàn)以它為例，來(lái)說(shuō)明用直接剛度法建立單元?jiǎng)偠染仃嚨乃枷牒瓦^(guò)程。圖5-4

平面簡(jiǎn)支梁元及其計(jì)算模型

梁在橫向外載荷（可以是集中力或分布力或力矩等）作用下產(chǎn)生彎曲變形，在水平載荷作用下產(chǎn)生線位移。

對(duì)于該平面簡(jiǎn)支梁?jiǎn)栴}：

梁上任一點(diǎn)受有三個(gè)力的作用：水平力Fx，

剪切力Fy，

和彎矩Mz。

相應(yīng)的位移為：水平線位移u，

撓度v，

和轉(zhuǎn)角z

。由上圖可見(jiàn)：

水平線位移和水平力向右為正，

撓度和剪切力向上為正，

轉(zhuǎn)角和彎矩逆時(shí)針?lè)较驗(yàn)檎?/p>

通常規(guī)定：

為使問(wèn)題簡(jiǎn)化，可把圖示的梁看作是一個(gè)梁?jiǎn)卧?。如圖5-4所示，當(dāng)令左支承點(diǎn)為節(jié)點(diǎn)i

，右支承點(diǎn)為節(jié)點(diǎn)

j

時(shí)，則該單元的節(jié)點(diǎn)位移和節(jié)點(diǎn)力可以分別表示為：稱為單元的節(jié)點(diǎn)位移列陣。稱為單元的節(jié)點(diǎn)力列陣；若{F}

為外載荷，則稱為載荷列陣。

（5-1）（5-2）寫(xiě)成矩陣形式為顯然，梁的節(jié)點(diǎn)力和節(jié)點(diǎn)位移是有聯(lián)系的。在彈性小變形范圍內(nèi)，這種關(guān)系是線性的，可用下式表示

或（5-3b）（5-3a）

上式(5-3b)稱為單元有限元方程，或稱為單元?jiǎng)偠确匠蹋砹藛卧妮d荷與位移之間（或力與變形之間）的聯(lián)系；式中，[K](e)稱為單元?jiǎng)偠染仃嚕菃卧奶匦跃仃嚒?/p>

對(duì)于圖5-4所示的平面梁?jiǎn)卧獑?wèn)題，利用材料力學(xué)中的桿件受力與變形間的關(guān)系及疊加原理，可以直接計(jì)算出單元?jiǎng)偠染仃嘯K](e)中的各系數(shù)kst(s,t=i,j)的數(shù)值，具體方法如下：

(1)假設(shè)ui=１，其余位移分量均為零，即vi=iz=uj=vj=zj=0，此時(shí)梁?jiǎn)卧鐖D5-5(a)所示，由梁的變形公式得圖5-5(a)伸縮：撓度：轉(zhuǎn)角：由上述三式可以解得根據(jù)靜力平衡條件

由式(5-3a)解得

(2)同理，設(shè)vi=１，其余位移分量均為零，即ui=iz=uj=vj=zj=0,

此時(shí)梁?jiǎn)卧鐖D5-5(b)所示，由梁的變形公式得伸縮：撓度：轉(zhuǎn)角：圖5-5(b)由上述三式可以解得利用靜力平衡條件由式(5-3a)解得

(3)同理，設(shè)，其余位移分量均為零，即此時(shí)梁?jiǎn)卧鐖D5-5(c)所示，由梁的變形公式得伸縮：撓度：轉(zhuǎn)角：由上述三式可以解得圖5-5(c)利用靜力平衡條件，可得

由式(5-3a)解得

剩余三種情況，仿此可推出。最后可以得到平面彎曲梁元的單元?jiǎng)偠染仃嚍檫@樣，便求得式（5-3a）單元?jiǎng)偠染仃囍星叭袆偠认禂?shù)?？梢钥闯觯琜K](e)為對(duì)稱矩陣。（5-4）2.虛功原理法

下面以平面問(wèn)題中的三角形單元為例，說(shuō)明利用虛功原理法來(lái)建立單元?jiǎng)偠染仃嚨牟襟E。

如前所述，將一個(gè)連續(xù)的彈性體分割為一定形狀和數(shù)量的單元，從而使連續(xù)體轉(zhuǎn)換為有限個(gè)單元組成的組合體。單元與單元之間僅通過(guò)節(jié)點(diǎn)連結(jié)，除此之外再無(wú)其他連結(jié)。也就是說(shuō)，一個(gè)單元上的只能通過(guò)節(jié)點(diǎn)傳遞到相鄰單元。

從分析對(duì)象的組合體中任取一個(gè)三角形單元：設(shè)其編號(hào)為e，

三個(gè)節(jié)點(diǎn)的編號(hào)為i、j、m，在定義的坐標(biāo)系

xoy中，節(jié)點(diǎn)坐標(biāo)分別為(xi,yi)、(xj,yj)、(xm,ym)，如圖5-6所示。圖5-6三節(jié)點(diǎn)三角形單元

由彈性力學(xué)平面問(wèn)題的特點(diǎn)可知，單元每個(gè)節(jié)點(diǎn)有兩個(gè)位移分量，即每個(gè)單元有6個(gè)自由度，相應(yīng)有6個(gè)節(jié)點(diǎn)載荷，寫(xiě)成矩陣形式，即單元節(jié)點(diǎn)載荷矩陣：?jiǎn)卧?jié)點(diǎn)位移矩陣：圖5-6三節(jié)點(diǎn)三角形單元(1)設(shè)定位移函數(shù)

按照有限元法的基本思想：首先需設(shè)定一種函數(shù)來(lái)近似表達(dá)單元內(nèi)部的實(shí)際位移分布，稱為位移函數(shù)，或位移模式。

三節(jié)點(diǎn)三角形單元有6個(gè)自由度，可以確定6個(gè)待定系數(shù)，故三角形單元的位移函數(shù)為（5-5）式(5-5)為線性多項(xiàng)式，稱為線性位移函數(shù)，相應(yīng)的單元稱為線性單元。上式(5-5)也可用矩陣形式表示，即

式中，jxjnxrp為單元內(nèi)任意點(diǎn)的位移列陣。

（5-6）

由于節(jié)點(diǎn)i、j、m

在單元上，它們的位移自然也就滿足位移函數(shù)式(5-5)。設(shè)三個(gè)節(jié)點(diǎn)的位移值分別為(ui,vi)、(uj,vj)、(um,vm)，將節(jié)點(diǎn)位移和節(jié)點(diǎn)坐標(biāo)代入式(5-5)，得（5-7）式中（5-8）

由上可知，共有6個(gè)方程，可以求出6個(gè)待定系數(shù)。解方程，求得各待定系數(shù)和節(jié)點(diǎn)位移之間的表達(dá)式為

為三角形單元的面積。其中：

（5-9）將式(5-7)及式(5-8)、式(5-9)代入式(5-6)中，得到式中，矩陣[N]

稱為單元的形函數(shù)矩陣；

為單元節(jié)點(diǎn)位移列陣。其中，為單元的形函數(shù)，它們反映單元內(nèi)位移的分布形態(tài)，是x,y坐標(biāo)的連續(xù)函數(shù)，且有

（5-12）式(5-10)又可以寫(xiě)成（5-13）上式清楚地表示了單元內(nèi)任意點(diǎn)位移可由節(jié)點(diǎn)位移插值求出。

(2)

利用幾何方程由位移函數(shù)求應(yīng)變

根據(jù)彈性力學(xué)的幾何方程

，

線應(yīng)變

剪切應(yīng)變

則應(yīng)變列陣可以寫(xiě)成式中，[B]稱為單元應(yīng)變矩陣，它是僅與單元幾何尺寸有關(guān)的常量矩陣，即（5-14）（5-15）

上述方程(5-14)稱為單元應(yīng)變方程，它的意義在于：

單元內(nèi)任意點(diǎn)的應(yīng)變分量亦可用基本未知量即節(jié)點(diǎn)位移分量來(lái)表示。

(3)利用廣義虎克定律求出單元應(yīng)力方程根據(jù)廣義虎克定律，對(duì)于平面應(yīng)力問(wèn)題上式（5-16）也可寫(xiě)成（5-16）（5-17）式中，為應(yīng)力列陣；

[D]

稱為彈性力學(xué)平面問(wèn)題的彈性矩陣，并有則有如下單元應(yīng)力方程由式(5-19)可求單元內(nèi)任意點(diǎn)的應(yīng)力分量，它也可用基本未知量即節(jié)點(diǎn)位移分量來(lái)表示。（5-18）（5-19）(4)由虛功原理求單元?jiǎng)偠染仃?/p>

根據(jù)虛功原理，當(dāng)彈性結(jié)構(gòu)受到外載荷作用處于平衡狀態(tài)時(shí)，在任意給出的微小的虛位移上，外力在虛位移上所做的虛功

AF等于結(jié)構(gòu)內(nèi)應(yīng)力在虛應(yīng)變上所存儲(chǔ)的虛變形勢(shì)能

A，即設(shè)處于平衡狀態(tài)的彈性結(jié)構(gòu)內(nèi)任一單元發(fā)生一個(gè)微小的虛位移，則單元各節(jié)點(diǎn)的虛位移為（5-21）（5-22）（5-20）則單元內(nèi)部必定產(chǎn)生相應(yīng)的虛應(yīng)變，故單元內(nèi)任一點(diǎn)的虛應(yīng)變?yōu)轱@然，虛應(yīng)變和虛位移之間關(guān)系為設(shè)節(jié)點(diǎn)力為

則外力虛功為

（5-25）（5-23）（5-24）單元內(nèi)的虛變形勢(shì)能為

根據(jù)虛功原理

因?yàn)椋?-27）（5-26）代人式(5-26)，則有式中，，均與坐標(biāo)

x,y無(wú)關(guān)