有限元與有限差分法基礎

1、1第二講第二講 有限元與有限差分法基礎有限元與有限差分法基礎 CAE的工具：的工具： 有限元法（有限元法（FEM）、有限差分法（有限差分法（FDM）、邊界元法邊界元法（BEM）、有限體積法（、有限體積法（FVM）、無網(wǎng)格法等等）、無網(wǎng)格法等等 在材料成形的在材料成形的CAE中主要使用的是中主要使用的是有限元法和有限差分法有限元法和有限差分法2“ ” 的的早在早在20世紀世紀40年代初期就年代初期就有人提出，但真正用于工程中則是有人提出，但真正用于工程中則是電子計算機電子計算機出現(xiàn)以后。出現(xiàn)以后。 “ ” 這一名稱是這一名稱是1960年美國的年美國的克拉夫克拉夫（Clough，R.W.）在一篇題

2、為）在一篇題為 “平面應力分析的有限元平面應力分析的有限元法法” 論文中首先使用。此后，論文中首先使用。此后，有限元法有限元法的應用得到蓬勃的應用得到蓬勃發(fā)展。發(fā)展。 到到20世紀世紀80年代初期國際上較大型的年代初期國際上較大型的結構分析結構分析有限元有限元通用程序通用程序多達多達幾百種幾百種，從而為，從而為工程應用工程應用提供了方便條件。提供了方便條件。由于有限元通用程序使用方便，計算精度高，其計算結果由于有限元通用程序使用方便，計算精度高，其計算結果已成為已成為各類工業(yè)產(chǎn)品設計各類工業(yè)產(chǎn)品設計和和性能分析性能分析的可靠依據(jù)。的可靠依據(jù)。 3 最初用于飛機結構的最初用于飛機結構的，由于它在

3、，由于它在理論上的通用性，因而它可用于解決工程中的許多問題。理論上的通用性，因而它可用于解決工程中的許多問題。目前，它可以解決幾乎所有的目前，它可以解決幾乎所有的連續(xù)介質連續(xù)介質和和場場的問題，的問題，包括熱傳導、電磁場、流體動力學、地質力學、原子工程包括熱傳導、電磁場、流體動力學、地質力學、原子工程和生物醫(yī)學等方面的問題。和生物醫(yī)學等方面的問題。 中，從齒輪、軸、軸承等通用零部件到機中，從齒輪、軸、軸承等通用零部件到機床、汽車、飛機等床、汽車、飛機等復雜結構復雜結構的應力和變形分析（包括熱應的應力和變形分析（包括熱應力和熱變形分析）。力和熱變形分析）。不僅可以解決工程中的線性問題、非線性問不

4、僅可以解決工程中的線性問題、非線性問題，而且對于各種不同性質的固體材料，如各向同性和各題，而且對于各種不同性質的固體材料，如各向同性和各向異性材料，粘彈性和粘塑性材料以及流體均能求解；向異性材料，粘彈性和粘塑性材料以及流體均能求解；對于工程中最有普遍意義的對于工程中最有普遍意義的非穩(wěn)態(tài)問題非穩(wěn)態(tài)問題也能求解。也能求解。42.1 有限元法基礎有限元法基礎基本思想：基本思想： 將一個連續(xù)求解域（對象）將一個連續(xù)求解域（對象）離散離散（剖分）成（剖分）成有限有限個形狀個形狀簡單的簡單的子域（單元）子域（單元） 利用有限個利用有限個節(jié)點節(jié)點將各將各子域子域連接起來連接起來 在給定的在給定的初始條件初始

5、條件和和邊界條件邊界條件下進行綜合計算求解，從下進行綜合計算求解，從而獲得對復雜工程問題的而獲得對復雜工程問題的近似數(shù)值解近似數(shù)值解 5物理系統(tǒng)舉例 幾何體幾何體 載荷載荷 物理系統(tǒng)物理系統(tǒng)結構結構熱熱電磁電磁6有限元模型真實系統(tǒng)真實系統(tǒng)有限元模型有限元模型 有限元模型有限元模型 是真實系統(tǒng)理想化的數(shù)學抽象是真實系統(tǒng)理想化的數(shù)學抽象。定義定義7自由度（DOFs）自由度自由度(DOFs) 用于描述一個物理場的響應特性用于描述一個物理場的響應特性。結構結構 DOFs 結構結構 位移位移 熱熱 溫度溫度 電電 電位電位 流體流體 壓力壓力 磁磁 磁位磁位 方向方向 自由度自由度ROTZUYROTYU

6、XROTXUZ8節(jié)點(node)和單元(element) 網(wǎng)格（grid）節(jié)點節(jié)點: 空間中的坐標位置，具有一定自由度空間中的坐標位置，具有一定自由度和和存在相互存在相互物理作用物理作用。單元單元: 一組節(jié)點自由度間相互作用的數(shù)值、一組節(jié)點自由度間相互作用的數(shù)值、矩陣描述（稱為剛度或系數(shù)矩陣矩陣描述（稱為剛度或系數(shù)矩陣)。單元有線。單元有線、面或實體以及二維或三維的單元等種類。、面或實體以及二維或三維的單元等種類。有限元模型由一些簡單形狀的單元組成，單元之間通有限元模型由一些簡單形狀的單元組成，單元之間通過節(jié)點連接，并承受一定載荷。過節(jié)點連接，并承受一定載荷。載荷載荷載荷載荷9節(jié)點和單元節(jié)點和

7、單元信息是通過單元之間的公共節(jié)點傳遞的。信息是通過單元之間的公共節(jié)點傳遞的。分離但節(jié)點重疊的單元分離但節(jié)點重疊的單元A和和B之間沒有信息傳遞之間沒有信息傳遞（需進行節(jié)點合并處理）（需進行節(jié)點合并處理）具有公共節(jié)點的單元具有公共節(jié)點的單元之間存在信息傳遞之間存在信息傳遞 .AB.AB.1 node2 nodes每個單元的特性是通過一些線性方程式來描述的。每個單元的特性是通過一些線性方程式來描述的。作為一個整體，單元形成了整體結構的數(shù)學模型。作為一個整體，單元形成了整體結構的數(shù)學模型。10節(jié)點和單元節(jié)點和單元節(jié)點自由度是隨連接該節(jié)點節(jié)點自由度是隨連接該節(jié)點 單元類型單元類型 變化的。變化的。JII

8、JJKLILKIPOMNKJIL三維桿單元三維桿單元 (鉸接鉸接)UX, UY, UZ三維梁單元三維梁單元二維或軸對稱實體單元二維或軸對稱實體單元UX, UY三維四邊形殼單元三維四邊形殼單元UX, UY, UZ,三維實體熱單元三維實體熱單元TEMPJPOMNKJIL三維實體結構單元三維實體結構單元ROTX, ROTY, ROTZROTX, ROTY, ROTZUX, UY, UZ,UX, UY, UZ11為什么要離散？為什么要離散？ 1.無法得到復雜實際問題的解析解無法得到復雜實際問題的解析解 2.將域劃分成一些微小而形狀規(guī)則的單元后，便于在一個將域劃分成一些微小而形狀規(guī)則的單元后，便于在一個

9、單元內(nèi)得到近似解單元內(nèi)得到近似解 3.域中所有單元的解可視為該復雜問題的近似解域中所有單元的解可視為該復雜問題的近似解12有限元分析的過程有限元分析的過程 1.連續(xù)體離散化連續(xù)體離散化 2.單元分析單元分析 3.整體分析整體分析 4.確定約束條件確定約束條件 5.方程求解方程求解 6.結果分析與討論結果分析與討論131.連續(xù)體離散化連續(xù)體離散化 連續(xù)體連續(xù)體：是指所：是指所求解的對象求解的對象（如物體或結構）。（如物體或結構）。 離散化離散化（劃分網(wǎng)格或網(wǎng)絡化）：是將所求解的對象（劃分網(wǎng)格或網(wǎng)絡化）：是將所求解的對象劃劃分分為有限為有限 個具有規(guī)則形狀的微小塊體，把每個微小塊體稱為個具有規(guī)則形

10、狀的微小塊體，把每個微小塊體稱為單元單元，相鄰兩個相鄰兩個 單元單元之間只通過之間只通過若干點若干點互相連接，每個連接點稱為互相連接，每個連接點稱為節(jié)點節(jié)點。 相鄰相鄰單元單元只在只在節(jié)點處節(jié)點處連接，連接，載荷載荷也只通過節(jié)點在各單也只通過節(jié)點在各單元之間傳元之間傳 遞，這些遞，這些有限個單元有限個單元的的集合體，集合體，即原來的即原來的連續(xù)體連續(xù)體。 *單元單元劃分后，給每個劃分后，給每個單元單元及及節(jié)點節(jié)點進行編號；進行編號； *選定坐標系，計算各個選定坐標系，計算各個節(jié)點坐標節(jié)點坐標； *確定各個確定各個單元單元的的形態(tài)和性態(tài)參數(shù)形態(tài)和性態(tài)參數(shù)以及以及邊界條件邊界條件等。等。14基本上

11、是任意的，一個結構體可以有基本上是任意的，一個結構體可以有。但應遵循以下。但應遵循以下劃分原則劃分原則：(1) 分析清楚所討論對象的性質，例如，是分析清楚所討論對象的性質，例如，是桁架結構桁架結構還是還是結構物結構物，是，是平面問題平面問題還是還是空間問題空間問題等等。等等。(2) 單元單元的幾何形狀取決于結構特點和受力情況，的幾何形狀取決于結構特點和受力情況，單單元元的的幾何尺寸幾何尺寸(大小大小)要按照要求確定。一般來說，要按照要求確定。一般來說，單元幾單元幾何形體何形體各邊的長度比各邊的長度比不能相差太大。不能相差太大。 (3) 有限元模型的有限元模型的網(wǎng)格網(wǎng)格劃分越密，其計算結果越精確

12、，劃分越密，其計算結果越精確，但計算工作量就越大。因此，在保證但計算工作量就越大。因此，在保證計算精度計算精度的前提下，的前提下，單元網(wǎng)格數(shù)量應盡量少。單元網(wǎng)格數(shù)量應盡量少。(4) 在進行在進行網(wǎng)格疏密布局網(wǎng)格疏密布局時，時，應力集中應力集中或或變形較大變形較大的的部位，部位，單元網(wǎng)格單元網(wǎng)格應取小一些，應取小一些，網(wǎng)格網(wǎng)格應劃分得密一些，而其應劃分得密一些，而其他部分則可疏一些。他部分則可疏一些。15(5) 在在設計對象設計對象的的厚度厚度或者或者彈性系數(shù)彈性系數(shù)有突變有突變的情況下，的情況下，應該取相應的突變線作為網(wǎng)格的應該取相應的突變線作為網(wǎng)格的邊界線邊界線；(6) 相鄰單元的相鄰單元的

13、邊界邊界必須相容，不能從一必須相容，不能從一單元單元的邊或者的邊或者面的內(nèi)部產(chǎn)生另一個單元的頂點。面的內(nèi)部產(chǎn)生另一個單元的頂點。(7) 網(wǎng)格劃分后，要將網(wǎng)格劃分后，要將全部單元全部單元和和節(jié)點節(jié)點按順序編號，不按順序編號，不允許有錯漏或者重復。允許有錯漏或者重復。(8) 劃分的劃分的單元單元集合成集合成整體整體后，應精確逼近原設計對象。后，應精確逼近原設計對象。原設計對象的原設計對象的各個頂點各個頂點都應該取成都應該取成單元的頂點單元的頂點。 所有網(wǎng)格的所有網(wǎng)格的表面頂點表面頂點都應該在原設計對象的都應該在原設計對象的表面上表面上。所。所有原設計對象的邊和面都應被有原設計對象的邊和面都應被單元

14、單元的的邊和面邊和面所逼近。所逼近。 16 圖例圖例 將將懸臂梁懸臂梁劃分為許多劃分為許多三角形單元三角形單元 三角形單元的三角形單元的三個頂點三個頂點都是都是節(jié)點節(jié)點 載荷載荷直接施加在直接施加在節(jié)點節(jié)點上上懸臂梁及其有限元模型懸臂梁及其有限元模型 172.單元分析單元分析 連續(xù)體離散化連續(xù)體離散化后，即可對后，即可對單元體單元體進行特性分析，簡稱為進行特性分析，簡稱為單元分析單元分析。單元分析工作單元分析工作主要有主要有兩項兩項：(1)選擇單元位移模式選擇單元位移模式(位移函數(shù)位移函數(shù)) 用用節(jié)點位移節(jié)點位移來表示來表示單元體內(nèi)任一點單元體內(nèi)任一點的的位移位移、應變應變和和應力應力，就需，

15、就需搞清各單元中的搞清各單元中的位移分布位移分布。 一般是假定一般是假定單元位移單元位移是坐標的某種簡單函數(shù)，用其模是坐標的某種簡單函數(shù)，用其模擬內(nèi)位移的分布規(guī)律，擬內(nèi)位移的分布規(guī)律，這種函數(shù)這種函數(shù)就稱為就稱為位移模式位移模式或或位移函位移函數(shù)數(shù)。通常采用的函數(shù)形式多為多項式。通常采用的函數(shù)形式多為多項式。 根據(jù)所選定的根據(jù)所選定的位移模式位移模式，就可以導出用，就可以導出用節(jié)點位移節(jié)點位移來來表示單元體內(nèi)表示單元體內(nèi)任一點位移的關系式任一點位移的關系式。182.單元分析單元分析(2) (2) 分析單元的特性，建立單元剛度矩陣分析單元的特性，建立單元剛度矩陣 進行進行單元單元力學特性分析力學

16、特性分析，將作用在，將作用在單元上單元上的的所有力所有力（表面（表面 力、體積力、集中力）等效地移置為力、體積力、集中力）等效地移置為節(jié)點載荷節(jié)點載荷； 采用有關的力學原理建立采用有關的力學原理建立單元的平衡方程單元的平衡方程，求得單元，求得單元內(nèi)內(nèi)節(jié)節(jié) 點位移點位移與與節(jié)點力節(jié)點力之間的關系矩陣之間的關系矩陣單元剛度矩陣單元剛度矩陣。 19 3. 整體分析整體分析 把把各個單元各個單元的的剛度矩陣剛度矩陣集成為集成為總體剛度矩陣總體剛度矩陣，以及將，以及將各單元的各單元的節(jié)點力向量節(jié)點力向量集成集成總的力向量總的力向量，求得，求得整體平衡整體平衡方程方程。 集成過程所依據(jù)的原理是集成過程所依

17、據(jù)的原理是節(jié)點變形協(xié)調條件節(jié)點變形協(xié)調條件和和平衡條平衡條件件。 204. 確定約束條件確定約束條件 由上述所形成的由上述所形成的是是一組線性代數(shù)方程一組線性代數(shù)方程，在求解之前，必修根據(jù)具體情況分析，確定在求解之前，必修根據(jù)具體情況分析，確定的的邊界約束條件邊界約束條件，并對，并對這些方程這些方程進行適當修正。進行適當修正。215. 有限元方程求解有限元方程求解 應用應用求解機械結構應力類問題時，根據(jù)求解機械結構應力類問題時，根據(jù)未知未知量量和和分析分析 有有三種基本解法三種基本解法： 位移法位移法 力法力法 混合法混合法22 (1)位移法位移法以以節(jié)點位移節(jié)點位移作為作為基本未知量基本未知

18、量，通過選擇適當?shù)?，通過選擇適當?shù)奈灰坪灰坪瘮?shù)數(shù)，進行，進行單元單元的力學特性分析。在的力學特性分析。在節(jié)點處節(jié)點處建立建立單元剛度方單元剛度方程程，再組合成，再組合成整體剛度矩陣整體剛度矩陣，求解出，求解出節(jié)點位移節(jié)點位移后，進而由后，進而由節(jié)點位移節(jié)點位移求解出求解出應力應力。 位移法優(yōu)點位移法優(yōu)點是比較簡單，規(guī)律性強，易于編寫計算機程是比較簡單，規(guī)律性強，易于編寫計算機程序。所以得到廣泛應用，其缺點是精度稍低。序。所以得到廣泛應用，其缺點是精度稍低。 (2)力法力法以以節(jié)點力節(jié)點力作為作為基本未知量基本未知量，在，在節(jié)點處節(jié)點處建立位移連續(xù)方建立位移連續(xù)方程，求解出程，求解出節(jié)點力節(jié)點

19、力后，再求解節(jié)點位移和單元應力。后，再求解節(jié)點位移和單元應力。力法的特點力法的特點是計算精度高。是計算精度高。 (3)混合法混合法取一部分取一部分節(jié)點位移節(jié)點位移和一部分和一部分節(jié)點力節(jié)點力作為作為基本未知量基本未知量，建立建立平衡方程平衡方程進行求解。進行求解。23單元特性的推導方法單元特性的推導方法 的推導是的推導是的的基本步驟之基本步驟之一一。目前，建立。目前，建立單元剛度矩陣單元剛度矩陣的方法主要有以下的方法主要有以下四種四種： 直接剛度法直接剛度法 虛功原理法虛功原理法 能量變分法能量變分法 加權殘數(shù)法加權殘數(shù)法241. 直接剛度法直接剛度法 是直接應用物理概念來建立是直接應用物理概

20、念來建立單元的有限元單元的有限元方程方程和分析單元特性的一種方法。和分析單元特性的一種方法。這一方法這一方法僅能適用于僅能適用于簡單形狀的單元簡單形狀的單元，如，如梁單元梁單元。但它可以幫助理解。但它可以幫助理解有限元有限元法法的物理概念。的物理概念。 圖圖1所示是所示是xoy平面中的一平面中的一，現(xiàn)以它為例，現(xiàn)以它為例，來說明用來說明用建立建立單元剛度矩陣單元剛度矩陣的思想和過程。的思想和過程。圖圖1平面簡支梁元及其計算模型平面簡支梁元及其計算模型25 梁在梁在橫向外載荷橫向外載荷（可以是集中力或分布力或力矩等）作用下產(chǎn)（可以是集中力或分布力或力矩等）作用下產(chǎn)生生彎曲變形彎曲變形，在，在水平

21、載荷水平載荷作用下產(chǎn)生作用下產(chǎn)生線位移線位移。 對于對于該平面簡支梁該平面簡支梁問題：問題：梁上任一點梁上任一點受有受有三個力三個力的作用：的作用： 水平力水平力Fx， 剪切力剪切力Fy ， 彎矩彎矩Mz。相應的位移相應的位移為：為： 水平線位移水平線位移u， 撓度撓度v ， 轉角轉角 z 。 由由上圖上圖可見：可見： 水平線位移水平線位移和和水平力水平力向右為正，向右為正， 撓度撓度和和剪剪切力切力向上為正，向上為正， 轉角轉角和和彎矩彎矩逆時針方向為正。逆時針方向為正。 通常規(guī)定：通常規(guī)定：26為使為使問題簡化問題簡化，可把，可把圖示的梁圖示的梁看作是一個看作是一個梁單元梁單元。如如圖圖1

22、所示，當令所示，當令左支承點左支承點為為節(jié)點節(jié)點 i ，右支承點右支承點為為節(jié)點節(jié)點 j 時，時，則則該單元該單元的的節(jié)點位移節(jié)點位移和和節(jié)點力節(jié)點力可以分別表示為：可以分別表示為：稱為稱為單元的節(jié)點位移列陣單元的節(jié)點位移列陣。稱為稱為；若；若 F 為為，則稱為，則稱為載荷列陣載荷列陣。 （1-1）（1-2）寫成寫成矩陣形式矩陣形式為為q(e)=ui ,vi , zi ,vj ,uj , zjTui ,vi , zi ,vj ,uj , zjF(e)=Fxi ,Fyi ,Mzi ,Fxj ,Fyj ,MzjTFxi ,Fyi ,Mzi ,Fxj ,Fyj ,Mzj27顯然，顯然，梁的節(jié)點力梁的

23、節(jié)點力和和節(jié)點位移節(jié)點位移是有聯(lián)系的。在是有聯(lián)系的。在彈性彈性小變形范圍小變形范圍內(nèi)，內(nèi)，這種關系這種關系是線性的，可用是線性的，可用下式下式表示表示 1112131415162122232425263132333435364142434445465152535455566162636 xiyizixjyjzjFkkkkkkFkkkkkkMkkkkkkFkkkkkkFkkkkkkkkkkM46566 iizijjzjuvuvkk ( )( )( )eeeFKq或或（1-3b）（1-3a）28上上式式(1-3b)稱為稱為，或稱為，或稱為單元剛度方單元剛度方程程，它代表了，它代表了單元的載荷單元的

24、載荷與與位移位移之間（或力與變形之間）之間（或力與變形之間）的聯(lián)系；的聯(lián)系；式中，式中，K(e)稱為稱為單元剛度矩陣單元剛度矩陣，它是單元的特性矩陣。，它是單元的特性矩陣。 對于對于圖圖1所示的所示的平面梁單元問題平面梁單元問題，利用材料力學中的桿，利用材料力學中的桿件受力與變形間的關系及疊加原理，可以直接計算出件受力與變形間的關系及疊加原理，可以直接計算出單元單元剛度矩陣剛度矩陣K(e)中的中的各系數(shù)各系數(shù) kst( s, t = i, j ) 的數(shù)值的數(shù)值292. 虛功原理法虛功原理法下面以下面以平面問題平面問題中的中的三角形單元三角形單元為例，說明利用為例，說明利用虛功原理法虛功原理法來

25、建來建立立單元剛度矩陣單元剛度矩陣的步驟。的步驟。如前所述，將一個如前所述，將一個連續(xù)的彈性體連續(xù)的彈性體分割為分割為一定形狀一定形狀和和數(shù)量數(shù)量的的單元單元，從而使從而使連續(xù)體連續(xù)體轉換為轉換為有限個單元有限個單元組成的組成的組合體組合體。單元單元與與單元單元之間僅通之間僅通過過節(jié)點節(jié)點連結，除此之外再無其他連結。也就是說，一個連結，除此之外再無其他連結。也就是說，一個單元單元上的只能上的只能通過通過節(jié)點節(jié)點傳遞到傳遞到相鄰單元相鄰單元。 從分析對象的從分析對象的組合體組合體中任取中任取一個一個三角形單元三角形單元：設設其編號其編號為為 e ，三個節(jié)點三個節(jié)點的的編號編號為為i、j、m，在定

26、義的在定義的坐標系坐標系 xoy 中，中，節(jié)點坐節(jié)點坐標標分別為分別為(x j , y j)、(xi , y i)、(xm, ym)，如，如圖圖2所示。所示。圖圖2三節(jié)點三角形單元三節(jié)點三角形單元30由由彈性力學彈性力學平面問題的特點可知，平面問題的特點可知，有有兩個位移分兩個位移分量量，即，即每個單元每個單元有有6個自由度個自由度，相應有，相應有6個節(jié)點載荷個節(jié)點載荷，寫成，寫成矩陣形式矩陣形式，即即 單元單元節(jié)點載荷矩陣節(jié)點載荷矩陣：F(e)=Fxi ,Fyi ,Fxj ,Fyj ,Fxm ,FymT單元單元節(jié)點位移矩陣節(jié)點位移矩陣：q(e)=ui ,vi ,uj ,vj ,um ,vmT

27、圖圖2三節(jié)點三角形單元三節(jié)點三角形單元31(1)設定位移函數(shù)設定位移函數(shù) 按照按照有限元法有限元法的的基本思想基本思想：首先需設定：首先需設定一種函數(shù)一種函數(shù)來近似表達單元來近似表達單元內(nèi)部的內(nèi)部的實際位移分布實際位移分布，稱為，稱為位移函數(shù)位移函數(shù)，或，或位移模式位移模式。 三節(jié)點三角形單元三節(jié)點三角形單元有有6個自由度個自由度，可以確定，可以確定 6個待定系數(shù)個待定系數(shù)，故，故三角三角形單元形單元的的位移函數(shù)位移函數(shù)為為 （1-4）式式(1-4)為為線性多項式線性多項式，稱為，稱為線性位移函數(shù)線性位移函數(shù)，相應的單元相應的單元稱為稱為線性單元線性單元。u=u(x,y)= 1+ 2x+ 3y

28、v=v(x,y)= 4+ 5x+ 6y32 1234561 0 0 00 0 0 1 uxydsvx y 上上式式(5-5)也可用也可用矩陣形式矩陣形式表示，即表示，即 式中，式中，d為為單元內(nèi)任意點單元內(nèi)任意點的的位移列陣位移列陣。 （1-5）33由于由于節(jié)點節(jié)點 i、j、m 在在單元單元上，上，它們的位移它們的位移自然也就滿足自然也就滿足位移位移函數(shù)式函數(shù)式(1-4)。設設三個節(jié)點三個節(jié)點的的位移值位移值分別為分別為( ui, vi)、( uj, vj )、( um, vm )，將，將節(jié)節(jié)點位移點位移和和節(jié)點坐標節(jié)點坐標代入代入式式(1-4)，得，得 123123123iiijjjmmmu

29、xyuxyuxy456456456iiijjjmmmvxyvxyvxy34( )12345600000000010002000000eijmiijmiijnjijmjijmmijmmaaaubbbvcccuaaavbbbucccv（1-6）111112221iijjijjmmijimjimmmxyxyx yx yx yx yx yx yxy 式中（1-7）由上可知，共有6個方程，可以求出6個待定系數(shù)。解方程，求得各待定系數(shù)和節(jié)點位移之間的表達式為 為三角形單元的面積。其中： 35, , , , , , ijmmjjmiimmijjiijmjmimijimjjimmjiax yx yax yx

30、yax yx ybyybyybyycxxcxxcxx（1-8）將式(1-7)及式(1-8)、式(1-9)代入式(1-6)中，得到 （1-9）（1-10）36式中，矩陣N 稱為單元的形函數(shù)矩陣； 為單元節(jié)點位移列陣。其中， 為單元的形函數(shù)，它們反映單元內(nèi)位移的 分布形態(tài)，是x, y 坐標的連續(xù)函數(shù)，且有 ( ) eq, , ijmNNN222iiiijjjjmmmmNab xc yNab xc yNab xc y（1-11）式(1-10)又可以寫成, , ,iijjmmiii i j miijjmmiii i j muN uN uN uN uvN vN vN vN v（1-12）上式清楚地表示了

31、單元內(nèi)任意點位移可由節(jié)點位移插值求出。 37(2) 利用幾何方程由位移函數(shù)求應變根據(jù)彈性力學的幾何方程 ，線應變 剪切應變 則應變列陣可以寫成 /, / ,xyuxuy / ,xyuyux ( )( ) 0001 0002iixijmeejyijmjxyiijjmmmmuuvxbbbuucccB qvycbcbcbuuuyxv式中，B稱為單元應變矩陣，它是僅與單元幾何尺寸有關的常量矩陣，即 （1-13）38 00010002ijmijmiijjmmbbbBccccbcbcb（1-14） 上述方程(1-13)稱為單元應變方程，它的意義在于： 單元內(nèi)任意點的應變分量亦可用基本未知量即節(jié)點位移分量來

32、表示。 39(3)利用廣義虎克定律求出單元應力方程根據(jù)廣義虎克定律，對于平面應力問題1()1()12(1)xxyyyxxyxyxyEEGE上式（1-15）也可寫成 ( )( )eeD（1-15）（1-16）式中， 為應力列陣； D 稱為彈性力學平面問題的彈性矩陣，并有 ,Txyxy 40 2101011002ED ( )( )( )( )exeeeyxyDDBq則有如下單元應力方程由式(1-18)可求單元內(nèi)任意點的應力分量，它也可用基本未知量即節(jié)點位移分量來表示。（1-17）（1-18）41(4)由虛功原理求單元剛度矩陣 根據(jù)虛功原理，當彈性結構受到外載荷作用處于平衡狀態(tài)時，在任意給出的微小的

33、虛位移上，外力在虛位移上所做的虛功 AF等于結構內(nèi)應力在虛應變上所存儲的虛變形勢能 A ，即FAA設處于平衡狀態(tài)的彈性結構內(nèi)任一單元發(fā)生一個微小的虛位移，則單元各節(jié)點的虛位移 為 ( )eq ( ), , , , , eTiijjmmquvuvuv（1-20）（1-21）（1-19）則單元內(nèi)部必定產(chǎn)生相應的虛應變，故單元內(nèi)任一點的虛應變 為 ( )e ( ), , eTxyxy42顯然，虛應變和虛位移之間關系為 ( )( )eeBq設節(jié)點力為 ( ), , , , , TexiyixjyjxmymFFFFFFF則外力虛功為 ( )( )TeeFAqF（1-24）（1-22）（1-23）單元內(nèi)的

34、虛變形勢能為 ( ) TeVAdv43根據(jù)虛功原理 ( )( )( )TTeeeVqFdv ( )( )eeDDBq ( )( )( )TTTeeeTBqBq因為（1-26）（1-25）代人式(1-25)，則有 ( )( )( )( )TTeeTeeVqFBqDBqdv式中， ，均與坐標 x, y 無關，故可以從積分符號中提出，可得： ( )Teq ( )eq44 ( )( )( )( )TeeeeVFBDB dvqKq其中，單元剛度矩陣 ( )eTVKBDB dv（1-27）式(1-27)稱為單元有限元方程，或稱單元剛度方程，其中 是單元剛度矩陣。 ( )eK（1-28）因為三角形單元是常應

35、變單元，其應變矩陣B 、彈性矩陣D均為常量，而 ，所以式(1-28)可以寫成 Vdvtdxdyt ( )eTKtBDB （1-29）45式中，t 為三角形單元的厚度； 為三角形單元的面積。 ( )( ) eiiijimejijjjmmimjmmKKKKKKKKKK 對于圖2所示的三角形單元，將D 及B代入式(1-28)，可以得到單元剛度為 （1-30）式中: K為66階矩陣，其中每個子矩陣為22階矩陣，由下式給出 211 22114(1) c22 ( , ,)rsrsrsrsrsrsrsrsrsb bc cb cc bEtKc bb ccb br si j m （1-31）46按照力學的一般說

36、法，任何一個實際狀態(tài)的是這個系統(tǒng)從實際狀態(tài)運動到某一參考狀態(tài)(通常取彈性體外載荷為零時狀態(tài)為參考狀態(tài))時。彈性體的總位能 是一個函數(shù)的函數(shù)，即泛函，位移是泛函的容許函數(shù)。從考慮，變形彈性體受外力作用處于平衡狀態(tài)時，在很多可能的變形狀態(tài)中，使總位能最小的就是彈性體的真正變形，這就是最小位能原理。用求能量泛函的極值方法就是能量變分原理。能量變分原理除了可解機械結構位移場問題以外，還擴展到求解熱傳導、電磁場、流體力學等連續(xù)性問題。3. 能量變分原理法47該方法是將假設的場變量的函數(shù)(稱為試函數(shù))引入問題的控制方程式及邊界條件，利用最小二乘法等方法使殘差最小，便得到近似的場變量函數(shù)形式。該方法的優(yōu)點是

37、不需要建立要解決問題的泛函式，所以，即使沒有泛函表達式也能解題。 4. 加權殘數(shù)法48 有限元解的收斂性有限元解的收斂性有限元解是近似解近似解是否收斂于真實解、近似解收斂速度、近似解的穩(wěn)定性 近似解的收斂條件:1.完備性準則（必要條件）試探函數(shù)（插值函數(shù)）的次數(shù)（m）不小于場函數(shù)的最高可導階數(shù)2.協(xié)調性準則（充分條件）試探函數(shù)在m-1次連續(xù)可導。49有限元分析的誤差有限元分析誤差建模誤差計算誤差離散誤差物理離散誤差幾何離散誤差邊界條件誤差單元形狀誤差舍入誤差截斷誤差插值函數(shù)與真實函數(shù)之間的差異1.減小單元特征尺寸，稱為減小單元特征尺寸，稱為h法法2.提高插值函數(shù)的階次，稱為提高插值函數(shù)的階次，

38、稱為p法法單元組合體與求解對象幾何形狀的差異1.網(wǎng)格局部加密網(wǎng)格局部加密2.選用邊或面上帶有節(jié)點的單元選用邊或面上帶有節(jié)點的單元邊界條件的復雜性1.準確測定，完善模型準確測定，完善模型2.細分邊界網(wǎng)格細分邊界網(wǎng)格單元嚴重畸變而退化細分局部網(wǎng)格或者控制調整關鍵區(qū)域的網(wǎng)格細分局部網(wǎng)格或者控制調整關鍵區(qū)域的網(wǎng)格數(shù)據(jù)儲存計算方法、解題性質、解題規(guī)模注意網(wǎng)格的劃分注意網(wǎng)格的劃分 選擇合適的解算方法選擇合適的解算方法 控制解題的規(guī)模控制解題的規(guī)模減少運算次數(shù)，降低解題規(guī)模減少運算次數(shù)，降低解題規(guī)模選擇合適的解算方法，控制解題規(guī)模選擇合適的解算方法，控制解題規(guī)模50材料成形中的非線性問題 1.材料非線性材料

39、本構方程非線性 彈塑性、剛彈性、剛黏塑性、黏彈塑性 2.幾何非線性 3.邊界非線性512.2 有限差分法基礎 一種直接將微分問題轉變成代數(shù)問題的近似數(shù)值解法。 基本思想基本思想 數(shù)值微分法數(shù)值微分法 是把連續(xù)的定解區(qū)域劃分成差分網(wǎng)格，用有限個節(jié)點代替原連續(xù)求解域。 把原方程和定解條件中的微商用差商來近似 把原微分方程和定解條件近似地用代數(shù)方程組代替，即有限差分方程52差分網(wǎng)格通常為矩形在邊界不規(guī)則或者形狀復雜時精度降低有限元網(wǎng)格有限差分網(wǎng)格53差分概念自變量x的解析函數(shù) y =f (x)，則有：dx,dy自變量和函數(shù)的微分 函數(shù)對自變量的一階導數(shù) 函數(shù)對自變量的一階差商00()( )limli

40、mxxdyyf xxf xdxxx yxdydx差商差商54差分方向 向前差分 向后差分 中心差分()( )yf xxf x ( )()yf xf xx 11()()22yf xxf xx 55差商的截斷誤差差商的截斷誤差 將函數(shù)將函數(shù)f (x+x)按按Taylor級數(shù)展開級數(shù)展開 向前向前 向后向后 中心中心234()()()( )( )( )( ) 0() )1!2!3!xxxf xxf xf xfxfxx23()( )()( )( )( ) 0() )2!3!f xxf xxxf xfxfxxx 23( )()()( )( )( )0() )2!3!f xf xxxxfxfxfxxx23

41、()()()22( )( )0() )2! 3!xxf xf xxfxfxxx56 二階中心差商 通常采用向前差商的向后差商 截斷誤差與(x)2 同一數(shù)量級 一階向前差商 一階向后差商 一階計算精度 一階中心差商 二階中心差商 二階計算精度2222232()2 ( )()()2()( )( )0() )4!yf xxf xf xxxxyxfxfxxx57 我們在彈性體上，用相隔等間距h而平行于坐標軸的兩組平行線織成正方形網(wǎng)格，x=y=h，如圖。 設f=f(x,y)為彈性體內(nèi)的某一個連續(xù)函數(shù)。該函數(shù)在平行于x軸的一根網(wǎng)線上，如在-上,它只隨x坐標的改變而變化。在鄰近結點處，函數(shù)f可展為泰勒級數(shù)如

42、下：.)(! 31)(! 21)(3003320022000 xxxfxxxfxxxfff58 我們將只考慮離開結點充分近的那些結點，即（x-x0）充分小。于是可不計（x-x0）的三次及更高次冪的各項，則上式簡寫為：在結點，x=x0-h, 在結點1, x=x0+h，代入(b) 得：)(.)(! 21)(20022000bxxxfxxxfff)(.20222003cxfhxfhff)(.20222001dxfhxfhff59聯(lián)立(c),(d),解得差分公式： )11(2310hffxf)21 (22031022hfffxf同理，在網(wǎng)線-上可得到差分公式)41 (22)31 (2042022420

43、hfffyfhffyf60從而可導出其它的差分公式如下：)()(461)()(241)()(461)()(41121042040448765432104022411931040447586022fffffhyffffffffffhyxffffffhxfffffhyxf61 相隔2h的兩結點處的函數(shù)值來表示中間結點處的一階導數(shù)值，可稱為中點導數(shù)公式。 以相鄰三結點處的函數(shù)值來表示一個端點處的一階導數(shù)值，可稱為端點導數(shù)公式。 中點導數(shù)公式與端點導數(shù)公式相比，精度較高。因為前者反映了結點兩邊的函數(shù)變化，而后者卻只反映了結點一邊的函數(shù)變化。62邊界元法簡介邊界元法簡介 邊界元法（邊界元法（bounda

44、ry element method） 一種結合有限元法和邊界積分法發(fā)展起來的一種一種結合有限元法和邊界積分法發(fā)展起來的一種新數(shù)值方新數(shù)值方法法 只在定義域的只在定義域的邊界邊界上上劃分單元劃分單元，用滿足控制方程的函數(shù)去，用滿足控制方程的函數(shù)去逼近逼近邊界條件邊界條件。 適用于應力（薄板）、適用于應力（薄板）、流體力學流體力學、聲場、電磁場等問題、聲場、電磁場等問題63邊界元法基本思想邊界元法基本思想 以微分控制方程的基本解為權函數(shù)，利用以微分控制方程的基本解為權函數(shù)，利用加權余量法加權余量法將區(qū)將區(qū)域積分轉化為邊界積分，并結合求解域邊界的離散，構建域積分轉化為邊界積分，并結合求解域邊界的離散

45、，構建基于邊界單元的代數(shù)方程組，然后進行計算求解基于邊界單元的代數(shù)方程組，然后進行計算求解以定義在邊界上的邊界積分方程為控制方程，通過對邊界元插值離散，以定義在邊界上的邊界積分方程為控制方程，通過對邊界元插值離散，化為代數(shù)方程組求解化為代數(shù)方程組求解 降低了問題的降低了問題的維數(shù)維數(shù)，從而顯著降低了，從而顯著降低了自由度數(shù)自由度數(shù)邊界的離散比區(qū)域的離散方便得多，可用較簡單的單元準確地模擬邊邊界的離散比區(qū)域的離散方便得多，可用較簡單的單元準確地模擬邊界形狀，最終得到階數(shù)較低的線性代數(shù)方程組界形狀，最終得到階數(shù)較低的線性代數(shù)方程組64加權余量法簡介加權余量法簡介一種直接從所需求解的一種直接從所需求

46、解的微分方程微分方程及及邊界條件邊界條件出發(fā)，尋求出發(fā)，尋求邊值問題近似邊值問題近似解解的數(shù)學方法。的數(shù)學方法。從從靜力靜力發(fā)展到了發(fā)展到了動力、穩(wěn)定、材料非線性動力、穩(wěn)定、材料非線性和和幾何非線性幾何非線性等各方面。等各方面。在求解域上建立一個在求解域上建立一個試函數(shù)試函數(shù) 試函數(shù)由完備函數(shù)集的子集構成。已被采用過的試函數(shù)有冪級數(shù)、三角級數(shù)、試函數(shù)由完備函數(shù)集的子集構成。已被采用過的試函數(shù)有冪級數(shù)、三角級數(shù)、樣條函數(shù)、貝賽爾函數(shù)、切比雪夫和勒讓德多項式等等。樣條函數(shù)、貝賽爾函數(shù)、切比雪夫和勒讓德多項式等等。試函數(shù)與真實解之間的偏差，即試函數(shù)與真實解之間的偏差，即余量余量（內(nèi)部和邊界）（內(nèi)部和

47、邊界）引入權函數(shù)，定義消除余量的條件引入權函數(shù)，定義消除余量的條件加權余量法就是一種定義近似解與真解之間余量，并設法使其最小的加權余量法就是一種定義近似解與真解之間余量，并設法使其最小的方法。方法。65設問題的控制微分方程為:在V域內(nèi) 在S邊界上 式中 : L、B分別為微分方程和邊界條件中的微分算子； f、g 為與未知函數(shù)u無關的已知函數(shù)域值； u為問題待求的未知函數(shù)。( )0B ug( )0L uf66當利用加權余量法求近似解時，首先在求解域上建立一個試函數(shù) ，一般具有如下形式：u 1niiiuC NNC式中： 待定系數(shù)，也可稱為廣義坐標；iC取自完備函數(shù)集的線性無關的基函數(shù)。iN由于 一

48、般只是待求函數(shù)u的近似解，因此代入后將得不到滿足，若記：u ( )( )IBRL ufRB ug在V域內(nèi)在S邊界上顯然 反映了試函數(shù)與真實解之間的偏差，它們分別稱做內(nèi)部和邊界余量余量。BIRR、67 若在域V內(nèi)引入內(nèi)部權函數(shù) ，在邊界S上引入邊界權函數(shù)則可建立n個消除余量的條件，一般可表示為：IWBW0(1,2, )IiIBiBVSW R dVW R dSin不同的權函數(shù) 和 反映了不同的消除余量的準則。從上式可以得到求解待定系數(shù)矩陣C的代數(shù)方程組。一經(jīng)解得待定系數(shù)，由式(5.1.3)即可得所需求解邊值問題的近似解。BiWIiW68IiW 由于試函數(shù) 的不同，余量 和 可有如下三種情況，依此加

49、權余量法可分為：1內(nèi)部法試函數(shù)滿足邊界條件，也即 此時消除余量的條件成為:2邊界法試函數(shù)滿足控制方程，也即 此時消除余量的條件為:( )0BRB ug( )0IRL uf0(1,2, )IiIVW R dVin0(1,2, )BiBSW R dSinu IRBR693混合法 試函數(shù)不滿足控制方程和邊界條件，此時消除余量的條件為: 顯然，混合法對于試函數(shù)的選取最方便，但在相同精度條件下，工作量最大。對內(nèi)部法和邊界法必須使基函數(shù)事先滿足一定條件，這對復雜結構分析往往有一定困難，但試函數(shù)一經(jīng)建立，其工作量較小。0(1,2, )IiIBiBVSW R dVW R dSin70邊界元法特點邊界元法特點 1 前處理工作量小前處理工作量小 2 解算規(guī)模小解算規(guī)模小 3 求解奇異性問題時計算精度高求解奇異性問題時計算精度高 4 在在載荷集中載荷集中和和半無限域半無限域等問題上有優(yōu)勢等問題上有優(yōu)勢 相對于有限元法，邊界元法發(fā)展較慢，相對于有限元法，邊界元法發(fā)展較慢，71有限元法解決應力集中問題 在應力分析中對于應力集中區(qū)域必須劃分很多的單元，從而增加了求解方程的階數(shù)，計算費用也就隨之增加 用位移型有限元法求解出的應力的精度低于位移的精度，對于一個比較復雜的問題必須劃分很多單元，相應的數(shù)據(jù)