有限元法基礎(chǔ)與板料成形模擬

計算機(jī)數(shù)字模擬技術(shù)有限元法基礎(chǔ)與板料成形模擬教材董湘懷.《材料成形計算機(jī)模擬》.北京：機(jī)械工業(yè)出版社，2001.

講課：12學(xué)時主要參考書王勖成，邵敏.《有限單元法基本原理與數(shù)值方法》.北京：清華大學(xué)出版社，1996.

R.D.庫克著，程耿東等譯.《有限元分析的概念和應(yīng)用》科學(xué)出版社.板料成形快速仿真軟件FASTAMP網(wǎng)址：軟件下載網(wǎng)址：/technology.htm上機(jī)地點(diǎn)：材料學(xué)院機(jī)房（10學(xué)時*2）上機(jī)實(shí)習(xí)軟件考試與成績有限元部分占總成績的50分.筆試:開卷考試,30分上機(jī):采用FASTAMP軟件或DYNAFORM軟件的沖壓成形分析報告2份,20分有限元法基礎(chǔ)有限元發(fā)展過程有限元應(yīng)用有限元發(fā)展方向有限元法的基本思想基本思想

1）將連續(xù)的求解系統(tǒng)離散為一組由節(jié)點(diǎn)相互聯(lián)在一起的單元組合體2）在每個單元內(nèi)假設(shè)近似函數(shù)來分片表示系統(tǒng)的求解場函數(shù)

有限元法的基本思想有限元法分類1）位移法：基于最小勢能原理或虛功原理

2）力法：基于最小余能原理3）雜交法：基于修正余能原理4）混合法：基于Reissner變分原理有限元法的基本思想位移法基本過程1）離散化過程

3）約束處理過程

2）單元平衡方程組裝過程

5）應(yīng)變、應(yīng)力回代過程

4）方程組求解過程

離散化過程最小勢能原理

彈性體的勢能為彈性體變形后所具有的內(nèi)能

為彈性體所受的外力功

離散化過程

為彈性體的應(yīng)變

為彈性體的應(yīng)力

u為彈性體的可容位移彈性體處于平衡狀態(tài)時，其勢能應(yīng)為最小

0離散化過程單元插值關(guān)系

單元幾何關(guān)系單元本構(gòu)關(guān)系

N為單元形函數(shù)矩陣

L為單元幾何微分算子為單元彈性矩陣

單元節(jié)點(diǎn)自由度向量離散化過程B稱為應(yīng)變矩陣

單元平衡方程或單元剛度方程

k稱為單元剛度矩陣

f稱為單元載荷向量

單元剛度矩陣的特性

對稱性

奇異性主元恒正且對角占優(yōu)離散化過程線彈性問題幾何方程—三維問題

三維問題線彈性問題幾何方程—二維問題

二維問題平面應(yīng)力和平面應(yīng)變狀態(tài)

線彈性問題幾何方程—二維問題

二維問題軸對稱狀態(tài)

線彈性問題幾何方程—一維問題

一維問題線彈性問題本構(gòu)方程—三維問題

三維問題E為彈性模量；為泊松比

線彈性問題本構(gòu)方程—平面應(yīng)力

二維問題平面應(yīng)力狀態(tài)

線彈性問題本構(gòu)方程—平面應(yīng)力

平面應(yīng)力狀態(tài)

線彈性問題本構(gòu)方程—平面應(yīng)變

二維問題平面應(yīng)變狀態(tài)

線彈性問題本構(gòu)方程—平面應(yīng)變

平面應(yīng)變狀態(tài)

線彈性問題本構(gòu)方程—軸對稱

二維問題軸對稱狀態(tài)

線彈性問題本構(gòu)方程—軸對稱

二維問題軸對稱狀態(tài)

線彈性問題本構(gòu)方程—軸對稱

軸對稱狀態(tài)

線彈性問題本構(gòu)方程—一維問題

一維問題常用單元模型

單元模型插值關(guān)系一一對應(yīng)單元類型一維單元、二維單元、三維單元等參單元、超參單元、次參單元常用單元模型一維單元

2節(jié)點(diǎn)線單元3節(jié)點(diǎn)線單元梁單元常用單元模型二維單元3節(jié)點(diǎn)三角形線性單元6節(jié)點(diǎn)三角形二次單元常用單元模型二維單元10節(jié)點(diǎn)三角形三次單元4節(jié)點(diǎn)四邊形雙線性單元常用單元模型二維單元8節(jié)點(diǎn)四邊形二次單元12節(jié)點(diǎn)四邊形三次單元常用單元模型三維單元4節(jié)點(diǎn)四面體線性單元10節(jié)點(diǎn)四面體二次單元常用單元模型三維單元8節(jié)點(diǎn)六面體線性單元20節(jié)點(diǎn)六面體二次單元常用單元模型準(zhǔn)三維空間單元桁架單元一維2節(jié)點(diǎn)線單元+單元局部隨體坐標(biāo)系

為什么要建立單元局部隨體坐標(biāo)系？簡化分析問題的復(fù)雜程度。在局部坐標(biāo)系中，空間桁架的每根桿每變成了一維2節(jié)點(diǎn)線單元常用單元模型準(zhǔn)三維空間單元框架單元三維梁單元+一維2節(jié)點(diǎn)線單元+單元局部隨體坐標(biāo)系

兩端都是剛性聯(lián)結(jié)

可以要承受拉壓、彎曲、扭轉(zhuǎn)3種變形模式

框架單元的特點(diǎn)常用單元模型準(zhǔn)三維空間單元板單元薄板單元中厚板單元彎曲和橫向剪切2種變形模式抵抗板的變形如果板很薄，忽略橫向剪切抗力，認(rèn)為抵抗載荷的主要因素是彎矩常用單元模型準(zhǔn)三維空間單元?dú)卧?/p>

抵抗拉壓變形的二維單元+板單元+單元局部隨體坐標(biāo)系。適合于薄殼單元和中厚殼單元從幾何上分為薄殼單元和中厚殼單元①組合單元常用單元模型準(zhǔn)三維空間單元②殼理論單元

由空間殼理論嚴(yán)格構(gòu)造的殼單元。適合于薄殼單元和中厚殼單元

③退化單元

由三維實(shí)體單元退化成的殼單元。只適合于中厚殼單元

單元模型構(gòu)造

有限元法的基本思想

通過單元分片近似，在每個單元內(nèi)假設(shè)近似函數(shù)來分片表示系統(tǒng)的場函數(shù)

選擇近似函數(shù)簡單、實(shí)用的原則在有限元法中，近似函數(shù)稱為插值函數(shù)

單元模型構(gòu)造插值函數(shù)

一般都采用多項(xiàng)式函數(shù)，主要原因是:

采用多項(xiàng)式插值函數(shù)比較容易推導(dǎo)單元平衡方程，特別是易于進(jìn)行微分和積分運(yùn)算。隨著多項(xiàng)式函數(shù)階次的增加，可以提高有限元法的計算精度。從理論上說，無限提高多項(xiàng)式的階數(shù)，可以求得系統(tǒng)的精確解。單元模型構(gòu)造方法

整體坐標(biāo)系法局部坐標(biāo)系法

Lagrange插值方法Hermite插值方法單元模型構(gòu)造方法2節(jié)點(diǎn)線單元12

oxu1u2x1x2ux1.假設(shè)插值多項(xiàng)式2.利用節(jié)點(diǎn)值求a0和a1單元模型構(gòu)造方法3.代入a0和a1，得插值多項(xiàng)式u(x)4.按u1和u2合并同類項(xiàng)，設(shè)l=x2-x1單元模型構(gòu)造方法關(guān)鍵

如何構(gòu)造插值多項(xiàng)式u?二維問題三維問題，如何構(gòu)造插值多項(xiàng)式?收斂性條件

①在單元內(nèi)，場函數(shù)必須是連續(xù)的；②完備性：插值多項(xiàng)式的階次必須由低到高依次增加，不能出現(xiàn)跳躍現(xiàn)象；

③協(xié)調(diào)性：各單元邊界必須連續(xù)，單元邊界不能出現(xiàn)開裂現(xiàn)象。插值多項(xiàng)式收斂性條件

收斂：當(dāng)單元逐漸縮小時，如果插值多項(xiàng)式滿足收斂性條件，則數(shù)值解將收斂于精確解

插值多項(xiàng)式收斂性條件協(xié)調(diào)單元

滿足插值多項(xiàng)式收斂性條件①和③的單元

完備單元

滿足插值多項(xiàng)式收斂性條件②的單元cr

階連續(xù)性

插值多項(xiàng)式的第r階導(dǎo)數(shù)是連續(xù)的

插值多項(xiàng)式收斂性條件非協(xié)調(diào)單元與部分協(xié)調(diào)單元對于一般固體力學(xué)問題來說，協(xié)調(diào)性要求單元在變形時，相鄰單元之間不應(yīng)引起開裂、重疊或其它不連續(xù)現(xiàn)象。例如，梁、板、殼等單元，在單元邊界不但要求位移是連續(xù)的，而且其一階導(dǎo)數(shù)也必須是連續(xù)的。板、殼單元位移函數(shù)沿單元邊界的法向?qū)?shù)（轉(zhuǎn)角）的連續(xù)性一般比較難實(shí)現(xiàn)，因此出現(xiàn)了許多不完全滿足協(xié)調(diào)性要求的“非協(xié)調(diào)單元”或“部分協(xié)調(diào)單元”，有時它們的精度也很好。

插值多項(xiàng)式選擇條件

插值多項(xiàng)式應(yīng)該盡可能滿足其收斂性條件（收斂性）由插值多項(xiàng)式所確定的場函數(shù)變化應(yīng)該與局部坐標(biāo)系的選擇無關(guān)（各向同性）

假設(shè)的插值多項(xiàng)式系數(shù)的數(shù)量應(yīng)該等于單元的節(jié)點(diǎn)數(shù)（解的唯一性）

選擇條件插值多項(xiàng)式選擇條件深入分析由收斂性條件②可知，插值多項(xiàng)式中必須含有常數(shù)項(xiàng)（剛體位移項(xiàng)），高階項(xiàng)的次數(shù)必須依次增加，不允許有跳躍插值多項(xiàng)式選擇條件由選擇條件②可知，插值多項(xiàng)式函數(shù)在所有自由度方向上要滿足各向同性性，這樣就不會隨局部坐標(biāo)系變化而改變了

深入分析插值多項(xiàng)式選擇條件深入分析選擇條件③是為了能由單元節(jié)點(diǎn)值唯一確定插值多項(xiàng)式

4節(jié)點(diǎn)四邊形的插值多項(xiàng)式應(yīng)該是

插值多項(xiàng)式系數(shù)i(i=0,1,2,3)

也是4個

單元模型構(gòu)造—整體坐標(biāo)系法基本思想

針對彈性體有限元網(wǎng)格建立一個統(tǒng)一的坐標(biāo)系，每個單元的插值多項(xiàng)式都在這個坐標(biāo)系上建立

單元模型構(gòu)造—整體坐標(biāo)系法2節(jié)點(diǎn)線單元12

oxu1u2x1x2ux1.假設(shè)插值多項(xiàng)式2.利用節(jié)點(diǎn)值求a0和a1單元模型構(gòu)造—整體坐標(biāo)系法3.代入a0和a1，得插值多項(xiàng)式u(x)4.按u1和u2合并同類項(xiàng)，設(shè)l=x2-x1單元模型構(gòu)造—整體坐標(biāo)系法N1和N2稱為單元的形函數(shù)；N稱為單元的形函數(shù)矩陣；u稱為單元節(jié)點(diǎn)位移向量。

2節(jié)點(diǎn)線的單元形函數(shù)單元模型構(gòu)造—整體坐標(biāo)系法二維3節(jié)點(diǎn)三角形單元

建立整體坐標(biāo)系oxy

單元模型構(gòu)造—整體坐標(biāo)系法1.假設(shè)插值多項(xiàng)式2.首先，利用節(jié)點(diǎn)值求0、

1和

2二維3節(jié)點(diǎn)三角形單元

單元模型構(gòu)造—整體坐標(biāo)系法A為單元面積單元模型構(gòu)造—整體坐標(biāo)系法3.將0、

1和

2代入插值多項(xiàng)式，按u1、u2、u3合并同類項(xiàng)單元模型構(gòu)造—整體坐標(biāo)系法4.同理可得單元模型構(gòu)造—整體坐標(biāo)系法5.單元插值多項(xiàng)式為單元模型構(gòu)造—整體坐標(biāo)系法6.單元插值多項(xiàng)式寫成矩陣形式（常用）單元模型構(gòu)造—整體坐標(biāo)系法7.單元插值多項(xiàng)式的另一種矩陣形式（不常用）單元模型構(gòu)造—整體坐標(biāo)系法4節(jié)點(diǎn)四面體單元單元模型構(gòu)造—整體坐標(biāo)系法1.假設(shè)插值多項(xiàng)式2.插值多項(xiàng)式為單元模型構(gòu)造—整體坐標(biāo)系法（i=1,2,3,4）循環(huán)輪換腳標(biāo)1、2、3、4，相應(yīng)可以得到a2，b2，c2，d2、a3，b3，c3，d3、a4，b4，c4，d4

單元模型構(gòu)造—整體坐標(biāo)系法3.單元插值多項(xiàng)式寫成矩陣形式（常用）單元模型構(gòu)造—整體坐標(biāo)系法4.單元插值多項(xiàng)式另一種矩陣形式（不常用）單元模型構(gòu)造—整體坐標(biāo)系法從理論上講，整體坐標(biāo)系法可以求任意單元的形函數(shù)，但計算過程太復(fù)雜只能求一維2節(jié)點(diǎn)線單元、二維3節(jié)點(diǎn)三角形單元和三維4節(jié)點(diǎn)四面體單元3種簡單單元的形函數(shù)復(fù)雜的或二次以上的單元必須采用局部坐標(biāo)系法求位移場u是形函數(shù)Ni的線性組合，因此形函數(shù)Ni同樣具有插值多項(xiàng)式的特性單元剛度矩陣—2節(jié)點(diǎn)線單元一維2節(jié)點(diǎn)線單元單元插值關(guān)系

單元幾何關(guān)系單元本構(gòu)關(guān)系

N=[N1N2]

De=E單元剛度矩陣—2節(jié)點(diǎn)線單元單元剛度矩陣A為單元截面積；l為單元長度矩陣B單元剛度矩陣—三角形單元二維3角形單元單元插值關(guān)系

單元剛度矩陣—三角形單元單元幾何關(guān)系單元剛度矩陣—三角形單元單元本構(gòu)關(guān)系

平面應(yīng)力問題單元剛度矩陣—三角形單元矩陣B單元剛度矩陣—三角形單元單元剛度矩陣h為單元厚度k為對稱的6*6常數(shù)矩陣A為單元面積單元模型構(gòu)造—整體坐標(biāo)系法單元形函數(shù)的特性正規(guī)性：單元形函數(shù)之和等于1。

正交性：形函數(shù)在本節(jié)點(diǎn)的值等于1，在其它節(jié)點(diǎn)的值等于0。

單元模型—等參單元等參單元

單元內(nèi)任意一點(diǎn)的位移u與單元節(jié)點(diǎn)位移ue之間的關(guān)系為

一般單元坐標(biāo)的插值關(guān)系也采用與位移插值關(guān)系相同的變換關(guān)系即單元內(nèi)任意一點(diǎn)的坐標(biāo)x與單元節(jié)點(diǎn)坐標(biāo)xe之間的關(guān)系為

單元模型—等參單元等參單元凡是幾何形狀和位移場采用同階同參數(shù)插值關(guān)系來描述的單元，稱為等參單元

前面介紹的所有單元都屬于等參單元

在描述單元的幾何形狀和位移場時，并不一定非采用同階插值關(guān)系

單元模型—等參單元等參單元3節(jié)點(diǎn)三角形等參單元

單元模型—等參單元超參單元如果幾何形狀插值函數(shù)的階數(shù)高于位移場插值函數(shù)的階數(shù)，稱為超參單元

次參單元如果幾何形狀插值函數(shù)的階數(shù)低于位移場插值函數(shù)的階數(shù)，稱為次參單元

單元平衡方程組裝過程

為什么要組裝?

消除內(nèi)力組裝的原則是什么?

單元自由度與結(jié)構(gòu)自由度對應(yīng)單元平衡方程組裝過程

2

F

1

3

①

②

③

U3U4U2U1U5U6結(jié)構(gòu)自由度向量U單元平衡方程組裝過程

2

1

①

U3U4U2U11’u1u2u3u4

3

②

U61’1

U2U12’u1u2U5u3u42’單元平衡方程組裝過程2

1

①

U3U4U2U11’u1u2u3u42’單元平衡方程組裝過程組裝單元①單元平衡方程組裝過程

3

②

U61’1

U2U12’u1u2U5u3u4單元平衡方程組裝過程再組裝單元②總體剛度方程

K稱為總體剛度矩陣

U稱為位移向量

F稱為載荷向量

總體剛度矩陣K的特性

對稱性

奇異性

稀疏性

非零元素帶狀分布約束處理過程

為什么要約束處理?總體平衡方程組是奇異的消除無限制的剛體運(yùn)動

使總體平衡方程組存在唯一一組解約束處理過程—邊界條件邊界條件分類

力(載荷)邊界條件位移邊界條件

集中載荷力

表面分布力

自重力熱交換引起的溫度載荷

固定位移約束

強(qiáng)制位移約束

關(guān)聯(lián)位移約束

約束處理過程—模型簡化xy約束處理過程—模型簡化yxxy約束處理過程—約束方程123456789101112yx約束處理過程—約束處理方法位移約束處理方法

賦0賦1法

乘大數(shù)法

約束處理過程—賦0賦1法基本原理利用初等變換對求解方程組進(jìn)行相同的行列變換，既保證方程組解不會改變，又可以保持方程組系數(shù)矩陣的對稱性。在進(jìn)行初等變換時，只要保證對方程組系數(shù)矩陣做相同的行列變換，就可以保持方程組系數(shù)矩陣的對稱性。

約束處理過程—乘大數(shù)法乘大數(shù)法基本原理利用矩陣的初等變換不改變方程組解的思想。

約束處理過程—兩種方法比較賦0賦1法在約束處理過程中是嚴(yán)格精確的，而乘大數(shù)法是一種近似約束處理方法，它的精度取決于所乘大數(shù)A值兩種方法都可以消除有限元平衡方程的奇異性，得到符合實(shí)際邊界條件的唯一一組解。但兩種方法還是有很大的區(qū)別約束處理過程—兩種方法比較采用乘大數(shù)法約束處理后的有限元平衡方程在求解時可能造成解的失真，大數(shù)A值越大可能解的偏差會越大，而賦0賦1法就不會出現(xiàn)類似的問題，它在約束過程和求解過程都是精確的乘大數(shù)法相對于賦0賦1法在約束處理過程上簡單一些

約束處理過程—兩種方法比較賦0賦1法實(shí)際上是將關(guān)聯(lián)位移約束方程代入到有限元平衡方程中的，是代入法。而乘大數(shù)是將占絕對優(yōu)勢的關(guān)聯(lián)位移約束方程合并到有限元平衡方程中的，是罰方法，計算誤差來自于合并過程，計算精度取決于關(guān)聯(lián)位移約束方程的優(yōu)勢大小商業(yè)軟件中，位移邊界條件的約束處理都采用賦0賦1法，乘大數(shù)很少被采用主要原因是它是一種近似方法，而且大數(shù)的大小也不好確定，有時還會造成求解失敗

方程組求解過程—特點(diǎn)方程組求解有限元計算過程中很重要的一部分，在有限元法的發(fā)展過程中，有限元方程的求解效率一直是其應(yīng)用的最大瓶頸之一有限元方程組的特點(diǎn)：有限元方程組的系數(shù)矩陣具有對稱、稀疏、帶狀分布以及正定、主元占優(yōu)。有效地利用這些特點(diǎn)，以減少系數(shù)矩陣的存貯量，提高方程組求解效率

方程組求解過程—分類比較線性方程組的解法主要分兩大類:

直接解法：以高斯消去法基礎(chǔ)，以等帶寬或變帶寬方式存貯系數(shù)矩陣內(nèi)元素，對于求解規(guī)模比較大的問題，要存貯的元素非常巨大。

迭代解法：只需要存貯系數(shù)矩陣中非零元素，存貯量很小，一般是變帶寬存貯量的20%或更少，有些算法的求解效率也非常高，適合求解大規(guī)模線性方程組。但是這種解法對接近病態(tài)的方程組很難保證收斂性。

方程組求解過程—帶寬定義有限元方程組系數(shù)矩陣是稀疏的、非零元素呈帶狀分布，帶寬就是它的寬度，帶寬的大小是由系統(tǒng)有限元網(wǎng)格的節(jié)點(diǎn)號排序決定的，具體求法是帶寬=(單元最大節(jié)點(diǎn)號之差+1)*節(jié)點(diǎn)自由度數(shù)

帶寬是網(wǎng)格節(jié)點(diǎn)標(biāo)注方法直接決定的，不同標(biāo)注方法帶寬可能相關(guān)很大

方程組求解過程—帶寬帶寬是網(wǎng)格節(jié)點(diǎn)標(biāo)注方法直接決定的，不同標(biāo)注方法帶寬可能相關(guān)很大

方程組求解過程—帶寬所示四邊形網(wǎng)格的三種節(jié)點(diǎn)號標(biāo)注方法，每個節(jié)點(diǎn)是2個自由度結(jié)構(gòu)的帶寬分別是12，18，56，相差很大，其中12和56之間相差近5倍，這就意味著系數(shù)矩陣的存貯量也是相差5倍，因此，對于大規(guī)模復(fù)雜系統(tǒng)的節(jié)點(diǎn)號優(yōu)化是十分必要的

方程組求解過程—系數(shù)矩陣存貯

系數(shù)矩陣存貯如果節(jié)點(diǎn)號排序優(yōu)化的比較好，系數(shù)矩陣的存貯量就會減少很多。根據(jù)系數(shù)矩陣的對稱性，一般都是按半帶寬存貯。系數(shù)矩陣存貯的方法二維等帶寬存貯一維變帶寬存貯

方程組求解過程—二維等帶寬存貯二維等帶寬存貯

方程組求解過程—二維等帶寬存貯二維等帶寬存貯消除了最大帶寬以外的全部零元素，節(jié)省了系數(shù)矩陣元素的存貯量。但是由于取最大帶寬為存貯范圍，因此不能排除在帶寬內(nèi)的大量零元素。當(dāng)系數(shù)矩陣的各行帶寬變化不大時，適合采用二維等帶寬存貯，方程組求解過程中系數(shù)矩陣元素的尋址也比較方便，求解效率較高。當(dāng)出現(xiàn)局部帶寬特別大的情況時，采用二維等帶寬存貯時，將由于局部帶寬過大而使整體系數(shù)矩陣的存貯大大增加。

方程組求解過程—一維變帶寬存貯

一維變帶寬存貯

一維變帶寬存貯方法就是把變化的帶寬內(nèi)的元素按一定的順序存貯在一個一維數(shù)組中。由于它不按最大帶寬存貯，因此比二維等帶寬存貯更節(jié)省內(nèi)存。按照解法可分為按行一維變帶寬存貯和按列一維變帶寬存貯。

按行一維變帶寬存貯

方程組求解過程—一維變帶寬存貯

輔助的尋址數(shù)組M

一維變帶寬存貯是最節(jié)省內(nèi)存的一種方法，但是由于要借助于尋址數(shù)組尋找系數(shù)矩陣元素的位置，相對二維等帶寬存貯方法來說要復(fù)雜一些，而且在程序?qū)崿F(xiàn)時也要復(fù)雜得多，方程組求解過程中也要消耗一些數(shù)組尋址時間。因此，在選用存貯方法時要權(quán)衡二者的利弊，統(tǒng)盤考慮。一般當(dāng)帶寬變化不大，計算機(jī)內(nèi)存允許時，采用二維等帶寬存貯方法是比較合適的。

方程組求解過程—一維變帶寬存貯

方程組求解過程—求解方法方程組求解方法高斯消去法

三角分解法

雅可比（Jacobi）迭代法

高斯-賽德爾（Gauss-Seidel）迭代法

應(yīng)變、應(yīng)力回代過程

單元應(yīng)變和應(yīng)力回代求解

通過求解有限元平衡方程得到有限元節(jié)點(diǎn)位移后，就可以進(jìn)行系統(tǒng)的剛度校核。如果所分析問題要進(jìn)行強(qiáng)度校核，就要回代求解單元的應(yīng)變和應(yīng)力。由插值關(guān)系和幾何關(guān)系可得單元應(yīng)變，再通過本構(gòu)關(guān)系得到單元應(yīng)力數(shù)值積分

為什么要進(jìn)行數(shù)值積分?2節(jié)點(diǎn)線單元、3節(jié)點(diǎn)三角形單元和4節(jié)點(diǎn)四面體單元3種單元的單元剛度矩陣是常數(shù)矩陣，不需要再進(jìn)行數(shù)值積分運(yùn)算。

除了這3種單元外，一般其它單元的剛度矩陣都是積分變量的函數(shù)，要采用數(shù)值積分方法進(jìn)行計算。

數(shù)值積分—主要方法常用的單元面內(nèi)數(shù)值積分方法主要有：Hammer積分

Gauss積分

數(shù)值積分—積分階次選擇

數(shù)值積分的階次選擇

求解單元平衡方程時，絕大多數(shù)情況要采用數(shù)值積分方法，如何選擇數(shù)值積分的階次將直接影響計算精度和計算量。如果積分階次選擇不當(dāng)，有時甚至?xí)?dǎo)致計算失敗

數(shù)值積分—積分階次選擇選擇積分階次的原則主要依據(jù)以下兩點(diǎn)：

積分精度

積分階次n與被積分多項(xiàng)式的階次m有直接關(guān)系。一般來說，有限元應(yīng)用的經(jīng)驗(yàn)公式積分項(xiàng)有兩個應(yīng)變矩陣B相乘，因此m一定是偶數(shù)，則積分階數(shù)n等于0.5、1.5、2.5、……數(shù)值積分—積分階次選擇常用單元的積分階次選擇

一維單元

一般都采用正規(guī)自然坐標(biāo)系法得到的形函數(shù)

在單元平衡方程中雅可比矩陣中雖然也含有自然坐標(biāo)，但是它只是單剛的一個系數(shù)，只對單剛中的每個元素的大小有相同的影響，不會改變單剛的特性

數(shù)值積分—積分階次選擇2節(jié)點(diǎn)線單元只能取高斯積分點(diǎn)n=13節(jié)點(diǎn)線單元可以取n=1或n=24節(jié)點(diǎn)線單元可以取n=2或n=32節(jié)點(diǎn)線單元：m=0，n=0.53節(jié)點(diǎn)線單元：m=2，n=1.54節(jié)點(diǎn)線單元：m=4，n=2.5按經(jīng)驗(yàn)公式計算：實(shí)際應(yīng)用計算：數(shù)值積分—積分階次選擇在有限元法中，把3節(jié)點(diǎn)單元取n=1以及4節(jié)點(diǎn)單元取n=2的積分方案稱為減縮積分，而3節(jié)點(diǎn)單元取n=2以及4節(jié)點(diǎn)單元取n=3的積分方案稱為正常積分

實(shí)際數(shù)值結(jié)果表明，有時減縮積分方案會帶來很大的計算誤差，產(chǎn)生零能模式

正常積分方案的有時計算結(jié)果也會偏小，產(chǎn)生閉鎖現(xiàn)象數(shù)值積分—積分階次選擇造成這些現(xiàn)象的原因有很多，例如，單元形狀、單元相對大小、單元受力狀況、分析問題的類型等等。為了避免零能模式和閉鎖現(xiàn)象的發(fā)生，一般采用減縮積分加阻尼矩陣方法。采用減縮積分方案時，對每個節(jié)點(diǎn)施加一個柔性彈簧，通過彈簧的阻尼增加剛度矩陣的穩(wěn)定性，阻止零能模式的發(fā)生。但是彈簧的剛性系數(shù)越大，計算誤差就越大，因此彈簧系數(shù)的選擇也有一定的困難

板料成形有限元法—分類靜力隱式有限元法靜力顯式有限元法動力顯式有限元法逆算法(一步成形)板料成形有限元法—單元模型單元模型薄膜單元薄殼單元中厚殼單元等效彎曲單元板料成形有限元法—薄膜單元薄膜單元是由二維三角形單元或四邊形單元構(gòu)造的空間板殼單元薄膜單元適用于板料在變形過程中主要以拉伸和壓縮變形為主，局部彎曲變形對整個成形問題不產(chǎn)生大的影響液壓脹形、半球沖頭脹形等一類問題板料成形有限元法—薄膜單元三角形薄膜單元

板料成形有限元法—薄膜單元局部坐標(biāo)系以1節(jié)點(diǎn)為坐標(biāo)原點(diǎn)，x軸與單元1、2邊重合并指向2節(jié)點(diǎn)，z軸與單元法向量nz平行。單位法向量nz為

板料成形有限元法—薄膜單元隨體局部坐標(biāo)系oxyz與空間整體坐標(biāo)系OXYZ之間坐標(biāo)轉(zhuǎn)換矩陣為

板料成形有限元法—薄膜單元單元局部坐標(biāo)自由度向量ue與整體坐標(biāo)自由度向量Ue的變換關(guān)系為

板料成形有限元法—薄膜單元四邊形薄膜單元

板料成形有限元法—薄膜單元局部坐標(biāo)系以1節(jié)點(diǎn)為坐標(biāo)原點(diǎn)，x軸與單元1、2邊重合并指向2節(jié)點(diǎn)，z軸與單元法向量nz平行。單位法向量nz為

板料成形有限元法—薄膜單元隨體局部坐標(biāo)系oxyz與空間整體坐標(biāo)系OXYZ之間坐標(biāo)轉(zhuǎn)換矩陣為

板料成形有限元法—薄膜單元單元局部坐標(biāo)自由度向量ue與整體坐標(biāo)自由度向量Ue的變換關(guān)系為

板料成形有限元法—單元平衡方程單元平衡方程單元插值關(guān)系(局部坐標(biāo)系)

單元幾何關(guān)系(局部坐標(biāo)系)單元本構(gòu)關(guān)系(局部坐標(biāo)系)

最小勢能原理(局部坐標(biāo)系)

坐標(biāo)變換關(guān)系

板料成形有限元法—單元平衡方程代入局部坐標(biāo)系的單元3個關(guān)系式,得代入坐標(biāo)轉(zhuǎn)換關(guān)系式,得板料成形有限元法—單元平衡方程ke

為局部坐標(biāo)系的單元剛度矩陣

fe

為局部坐標(biāo)系的單元載荷向量

Ke

為整體坐標(biāo)系的單元剛度矩陣

Fe

為整體坐標(biāo)系的單元載荷向量

板料成形有限元法—本構(gòu)方程

在金屬塑性大變形有限元分析時經(jīng)常采用流動理論本構(gòu)方程，其它本構(gòu)方程很少采用基于形變理論或非經(jīng)典的角點(diǎn)理論本構(gòu)方程雖然可以比較準(zhǔn)確板料失穩(wěn)后的局部化變形過程，但是板料成形屬于強(qiáng)約束過程，對角點(diǎn)本構(gòu)方程不敏感，而且板料成形也并不十分關(guān)心板料失穩(wěn)后的局部化變形過程目前板料成形有限元模擬的精度也不高，也沒有必要十分強(qiáng)調(diào)本構(gòu)方程的影響

板料成形有限元法—本構(gòu)方程考慮具有光滑屈服面屈服函數(shù)的彈塑性體，假設(shè)溫度對變形速度的影響很小，可以忽略不計。這樣全應(yīng)變率可以分解為彈性應(yīng)變率和塑性應(yīng)變率之和，即

板料成形有限元法—本構(gòu)方程采用Hooke定律，彈性應(yīng)變率

PiolaKirchhoff應(yīng)力的Jaumann速率

用第二表示為

塑性應(yīng)變率用流動法則和屈服函數(shù)f表示為板料成形有限元法—本構(gòu)方程nij是屈服面（應(yīng)力空間f=0的曲面）的單位法向量=0，應(yīng)力處于彈性狀態(tài),應(yīng)力點(diǎn)位于屈服面以內(nèi)

=1，應(yīng)力處于塑性加載狀態(tài),應(yīng)力點(diǎn)位于屈服面上h表示當(dāng)前狀態(tài)的加工硬化率

板料成形有限元法—本構(gòu)方程式中應(yīng)變率板料成形有限元法—本構(gòu)方程由

與Cauchy應(yīng)力的Jaumann速率的關(guān)系

如果材料不可壓縮

材料本構(gòu)方程為

板料成形有限元法—J2流動理論J2流動理論

Mises屈服函數(shù)為

表示應(yīng)力的偏量板料成形有限元法—J2流動理論材料的J2流動本構(gòu)關(guān)系

式中

h可以由單向拉伸實(shí)驗(yàn)確定

板料成形有限元法—J2流動理論E為彈性模量Et為單向拉伸真應(yīng)力—對數(shù)應(yīng)變曲線的切線模量K為材料強(qiáng)化系數(shù);n

為材料強(qiáng)化指數(shù)板料成形有限元法—屈服函數(shù)Hill正交各向異性函數(shù)

一般的，若把各向異性主軸作為隨體坐標(biāo)系的x,y,z軸，則Hill屈服函數(shù)可以表示成F、G、H、L、M、N為各向異性參數(shù)，由實(shí)驗(yàn)確定

板料成形有限元法—屈服函數(shù)在平面應(yīng)力狀態(tài)下Hill各向異性函數(shù)

板料成形有限元法—迭代解法非線性方程組迭代解法板材沖壓成形數(shù)值模擬是一個強(qiáng)非線性問題，涉及到幾何、材料和邊界三重非線性。如果采用隱式有限元法就要求解非線性有限元方程組。非線性方程組一般是采用線性化方法，通過一系列線性解逼近非線性解。但是這種方法是有局限性的，而且有時解的漂移誤差很大。因此，一般都采用迭代法求解非線性有限元方程組，即Newton-Raphson法

板料成形有限元法—迭代解法Newton-Raphson法

對具有一階導(dǎo)數(shù)連續(xù)的函數(shù)F(u)在Un處作一階Taylor展開，并用un表Un，則它在Un處的線性近似公式為板料成形有限元法—迭代解法非線性方程組在Un附近的近似方程組是一個線性方程組，即假設(shè)因此Newton-Raphson法的迭代方程為

板料成形有限元法—迭代解法修正的Newton-Raphson法在Newton-Raphson法中，剛度矩陣是與un有關(guān)的，在每個迭代步都要重新計算一次。為了減少計算量，用某一不變的剛度矩陣代替，得到修正的Newton-Raphson法迭代方程

板料成形有限元法—迭代解法擬Newton-Raphson法

Newton-Raphson法和修正的Newton-Raphson法都是用切線剛度矩陣進(jìn)行迭代平衡的。實(shí)際計算表明，這種直接迭代法不僅計算量大，而且常常不收斂。因此又提出一種擬Newton-Raphson法，用割線剛度矩陣進(jìn)行迭代平衡。

板料成形有限元法—收斂準(zhǔn)則

迭代收斂準(zhǔn)則

位移收斂準(zhǔn)則

常用位移收斂準(zhǔn)則為當(dāng)系統(tǒng)含有剛體位移時，會比較大，此時不適合采用位移收斂準(zhǔn)則板料成形有限元法—收斂準(zhǔn)則常用失衡力收斂準(zhǔn)則為為第n迭代步的失衡力失衡力收斂準(zhǔn)則當(dāng)系統(tǒng)處于失穩(wěn)狀態(tài)時，失衡力的微小變化將引起位移增量的很大偏差，此時不能采用失衡力收斂準(zhǔn)則

板料成形有限元法—收斂準(zhǔn)則能量收斂準(zhǔn)則常用能量收斂準(zhǔn)則為位移收斂準(zhǔn)則和失衡力收斂準(zhǔn)則都有一定的缺陷，相對而言能量收斂準(zhǔn)則是比較好的，因?yàn)樗瑫r控制位移增量和失衡力。能量收斂準(zhǔn)則是把每次迭代后的內(nèi)能增量與初始內(nèi)能增量相比較，即板料成形有限元法—主要方法板料成形有限元方法主要有：基于Green應(yīng)變和第二類Kirchhoff應(yīng)力能量共軛的虛功原理，以初始時刻為參考構(gòu)形的全量拉格朗日（TL）有限元方法基于Green應(yīng)變和第二類Kirchhoff應(yīng)力能量共軛的虛功原理，以當(dāng)前時刻為參考構(gòu)形的修正拉格朗日（UL）有限元方法板料成形有限元法—主要方法基于Lagrange（第一類Kirchhoff）應(yīng)力與速度對物質(zhì)坐標(biāo)偏導(dǎo)數(shù)能量共軛的虛功原理，以當(dāng)前時刻為參考構(gòu)形的虛功率增量型有限元方法基于動力學(xué)原理和中心差分方法建立的動力顯式有限元方法基于全量應(yīng)變空間，全量形變理論和虛功原理的有限元逆算法板料成形有限元法—TL法TL法TL法是取初始時刻構(gòu)形作為參考構(gòu)形，在所有的時間步長內(nèi)的計算都參照時刻t0=0構(gòu)形來定義。則t+t時刻的虛功方程為

板料成形有限元法—TL法單元插值關(guān)系

單元本構(gòu)關(guān)系

為時間增量步t內(nèi)第二類Kirchhoff應(yīng)力增量

為時間增量步t內(nèi)Green應(yīng)變增量

為切線本構(gòu)矩陣

板料成形有限元法—TL法單元幾何關(guān)系t時刻

t+t時刻

板料成形有限元法—TL法線性化的TL法有限元平衡方程

KT為t時刻的單元切線剛度矩陣Ks為初應(yīng)力單元剛度矩陣

板料成形有限元法—UL法UL法

UL法是取當(dāng)前時刻構(gòu)形作為參考構(gòu)形，在所有的時間步長內(nèi)的計算都參照當(dāng)前時刻t構(gòu)形來定義。則t+t時刻的虛功方程可以表示為

板料成形有限元法—UL法單元插值關(guān)系

單元本構(gòu)關(guān)系

為時間增量步t內(nèi)第二類Kirchhoff應(yīng)力增量

為時間增量步t內(nèi)Green應(yīng)變增量

為切線本構(gòu)矩陣

板料成形有限元法—UL法單元幾何關(guān)系t時刻

t+t時刻

板料成形有限元法—UL法線性化的UL法有限元平衡方程

Ks為初應(yīng)力單元剛度矩陣

Rs為初應(yīng)力節(jié)點(diǎn)力向量

板料成形有限元法—虛功率增率法虛功原理

由于Lagrange（第一類Kirchhoff）應(yīng)力與速度對物質(zhì)坐標(biāo)偏導(dǎo)數(shù)是能量共軛的，可得時刻現(xiàn)實(shí)構(gòu)形，以初始時刻為參考構(gòu)形的彈塑性大變形虛功方程為

板料成形有限元法—虛功率增率法單元插值關(guān)系

單元本構(gòu)關(guān)系

單元內(nèi)任意一點(diǎn)的速度v與單元節(jié)點(diǎn)速度向量ve的關(guān)系為

板料成形有限元法—虛功率增率法單元幾何關(guān)系變形率張量dij為

為速度梯度張量

以當(dāng)前時刻為參考構(gòu)形時，變形率張量dij等于應(yīng)變率張量

線性化的虛功率增率法有限元平衡方程

板料成形有