基于微態(tài)方法的耦合韌性損傷的彈塑性本構(gòu)模型

， ， Khemais Saanouni，

(1.大連理工大學(xué) 工業(yè)裝備結(jié)構(gòu)分析國(guó)家重點(diǎn)實(shí)驗(yàn)室，大連 116024；2.ICD/LASMIS，University of Technology of Troyes，Troyes 10000)

1 引 言

一般塑性材料在受力變形達(dá)到拉伸極限強(qiáng)度后，會(huì)出現(xiàn)明顯的應(yīng)變軟化直至材料最終斷裂。由材料微缺陷、損傷或溫度等引起的應(yīng)變軟化現(xiàn)象導(dǎo)致了材料失穩(wěn)，進(jìn)而失去承載能力，同時(shí)使得采用局部連續(xù)介質(zhì)損傷力學(xué)模型的準(zhǔn)靜態(tài)問(wèn)題控制方程橢圓性喪失，成為數(shù)學(xué)上病態(tài)的問(wèn)題，且數(shù)值響應(yīng)病態(tài)依賴(lài)于有限元網(wǎng)格，而收斂到不正確的沒(méi)有物理意義的有限元解上[1]。

針對(duì)上述問(wèn)題，已展開(kāi)大量研究，但有效的且根本性的補(bǔ)救措施是通過(guò)在本構(gòu)方程中引入特定尺寸內(nèi)物質(zhì)點(diǎn)間的相互作用，以保證控制方程的橢圓性，即正則化機(jī)制。當(dāng)前引入正則化機(jī)制的途徑主要有四條。一是將經(jīng)典連續(xù)體模型擴(kuò)展為率相關(guān)材料模型[2-4]，該粘性效應(yīng)隱式包含了一種內(nèi)長(zhǎng)度參數(shù)的作用，可以減弱數(shù)值結(jié)果對(duì)網(wǎng)格尺寸的依賴(lài)性。二是采用梯度理論模型[5,6]，最初是由Mindlin[7]提出的應(yīng)變的二階梯度線(xiàn)彈性模型。Triantafyllidis等[8]提出在超彈性材料的應(yīng)變能密度函數(shù)中引入二階變形梯度張量的依賴(lài)項(xiàng)，以保證控制方程的橢圓性。De Borst等[9]提出了廣義梯度塑性模型，并在屈服函數(shù)中引入了等效塑性應(yīng)變的高階空間梯度。三是引入額外自由度以表征材料的非局部行為，Cosserat介質(zhì)[10]引入了旋轉(zhuǎn)自由度并相應(yīng)產(chǎn)生了微曲率，進(jìn)而引入了與微曲率能量共軛的對(duì)偶應(yīng)力和具有特征長(zhǎng)度意義的內(nèi)部長(zhǎng)度參數(shù)。De Borst[11]推導(dǎo)了Mises屈服準(zhǔn)則下的彈塑性Cosserat連續(xù)體模型。李錫夔等[12,13]提出了壓力相關(guān)彈塑性Cosserat連續(xù)體模型，考慮了非關(guān)聯(lián)Drucker-Prager屈服準(zhǔn)則，并對(duì)其本構(gòu)參數(shù)的取值進(jìn)行了分析。Micromorphic介質(zhì)最初由Eringen等[14-16]提出，視每個(gè)質(zhì)點(diǎn)為一個(gè)微小物體，并允許微小物體變形和做剛體轉(zhuǎn)動(dòng)。Forest等[17,18]發(fā)展了包含新的微態(tài)自由度的基于廣義虛功原理和擴(kuò)展的亥姆霍茲自由能的微態(tài)方法。同時(shí)，Eringen的微態(tài)連續(xù)體模型可視為微應(yīng)變模型與Cosserat模型的結(jié)合體，分別對(duì)應(yīng)于微態(tài)方法中微態(tài)變量取為局部應(yīng)變和旋轉(zhuǎn)張量時(shí)的模型。四是采用完全非局部的連續(xù)體理論[19]，其基本思想不再基于局部行為假設(shè)，而是假定整個(gè)域內(nèi)的物質(zhì)點(diǎn)均處于相互作用的狀態(tài)。Ba?ant等[1,20,21]以此建立了針對(duì)控制應(yīng)變硬化所有變量的積分形式的正則化機(jī)制，并應(yīng)用于混凝土的應(yīng)變軟化分析中。

2 本構(gòu)模型

從虛功原理出發(fā)，通過(guò)引入與損傷內(nèi)變量相關(guān)聯(lián)的微態(tài)韌性損傷，以反映損傷演化過(guò)程中的非局部效應(yīng)；在熱力學(xué)定律和克勞修斯不等式的幫助下，通過(guò)亥姆霍茲自由能求解得到各狀態(tài)變量的表達(dá)式；應(yīng)用廣義法向法則根據(jù)彈塑性耗散勢(shì)能獲得耗散內(nèi)變量的演化方程。

2.1 平衡方程

(1)

相應(yīng)地，外力虛功和慣性力虛功可表示為

(2)

應(yīng)用廣義的虛功原理Pint+Pext=Pa,可以得到經(jīng)典的位移平衡方程及邊界條件：

(3)

額外的與新引入的微態(tài)損傷相關(guān)的廣義平衡方程及邊界條件：

(4)

2.2 狀態(tài)方程

微分形式的熱力學(xué)第一定律可擴(kuò)展為

(5)

考慮等溫過(guò)程，結(jié)合式(5)可得到Clausius-Duhem不等式為

(6)

為方便起見(jiàn)，在缺少實(shí)驗(yàn)支撐的前提下，本文假設(shè)微態(tài)損傷變量引起的耗散可以忽略不計(jì)，則通過(guò)將亥姆霍茲自由能函數(shù)代入不等式(6)中，可得狀態(tài)方程為

(7)

與經(jīng)典彈塑性損傷模型一致的機(jī)械耗散項(xiàng)為

(8)

式中R和Xi j分別為與各向同性硬化應(yīng)變和隨動(dòng)硬化應(yīng)變相共軛的應(yīng)力，Y為與損傷共軛的應(yīng)力，也稱(chēng)為損傷能量釋放率。

考慮如下形式的亥姆霍茲自由能：

(9)

通過(guò)狀態(tài)方程(7)，可得各廣義應(yīng)力變量的表達(dá)式為

σ=(1-d)Λεe，R=(1-d)Qr，X=(1-d)Cα

(10)

損傷能量釋放率可進(jìn)一步分解為兩部分之和：

(11)

式中Yloc>0是經(jīng)典彈塑性損傷模型中的損傷能量釋放率，Yn l來(lái)自于損傷和微態(tài)損傷的耦合作用，同時(shí)Y恒為正。

2.3 演化方程

采用式(12)定義的Von Mises屈服函數(shù)和耗散勢(shì)能。

(12)

采用廣義法向非關(guān)聯(lián)流動(dòng)法則得到耗散內(nèi)變量的演化方程為

(13)

(14)

式中Hp為彈塑性硬化模量,

(15)

進(jìn)而可求得應(yīng)變率張量

(16)

其中，四階彈塑性切線(xiàn)模量可表示為

Y(nΛ)?[σ/(1-d)]}

(17)

3 有限元離散

3.1 應(yīng)變空間的平衡方程

利用式(10)的狀態(tài)方程，將式(4)的微態(tài)損傷平衡方程變換到微態(tài)損傷應(yīng)變空間，并忽略其廣義體力與接觸力，可以得到簡(jiǎn)潔形式的偏微分方程和邊界條件：

(18)

式(18)與Engelen等在隱式的梯度損傷模型中提出的亥姆霍茲方程是一致的。

3.2 弱形式的初始邊界值問(wèn)題

采用典型的加權(quán)余量法對(duì)式(3，18)的偏微分方程和Neumann邊界條件進(jìn)行弱化，得到其對(duì)應(yīng)的積分形式的初始邊界值問(wèn)題：

(19)

3.3 有限元空間離散

采用基于Hu-Washizu變分原理的假設(shè)應(yīng)變四邊形單元以減少沙漏和閉鎖現(xiàn)象的影響，則單元離散的弱形式初始邊界值問(wèn)題為

(20)

式中 上標(biāo)e代表參考單元。與位移相關(guān)的阻尼質(zhì)量矩陣和內(nèi)外力向量分別為

(21)

與微態(tài)損傷對(duì)應(yīng)的阻尼質(zhì)量矩陣和內(nèi)外力向量分別為

(22)

由于動(dòng)力顯式求解過(guò)程是條件穩(wěn)定的，須確定穩(wěn)定的最大時(shí)間步長(zhǎng)Δt≤min(Le/Cd),其中，Le為單元的特征長(zhǎng)度，Cd為等效的材料聲速，可通過(guò)Lame常數(shù)近似確定。

3.4 本構(gòu)積分過(guò)程

為了計(jì)算內(nèi)外力向量以求解式(20)的代數(shù)方程組，需要根據(jù)已知tn時(shí)刻的各狀態(tài)變量得到tn +1時(shí)刻的應(yīng)力張量、硬化變量、局部損傷和微態(tài)損傷因子。本文采用完全隱式的彈性預(yù)測(cè)與塑性修正的Simo圖形返回算法更新該非關(guān)聯(lián)彈塑性模型的狀態(tài)變量。

彈性預(yù)測(cè)

(23)

判斷屈服函數(shù)：

(24)

塑性修正

在tn +1時(shí)刻的塑性應(yīng)變與應(yīng)力張量：

(25)

損傷及損傷能量釋放率：

(26)

各項(xiàng)同性硬化應(yīng)變及應(yīng)力：

(27)

隨動(dòng)硬化應(yīng)變及應(yīng)力：

(28)

Δλ>0時(shí)，滿(mǎn)足tn +1時(shí)刻的屈服函數(shù)：

(29)

4 算例與討論

4.1 DP1000材料的拉伸實(shí)驗(yàn)

將對(duì)偶雙相鋼DP1000材料在Zwick 250拉伸機(jī)中以0.1 mm/s拉伸速度進(jìn)行準(zhǔn)靜態(tài)單向拉伸實(shí)驗(yàn)。樣件厚度為1.5 mm，其幾何形狀尺寸如圖1所示。擬合的材料參數(shù)列入表1。

采用局部彈塑性損傷模型計(jì)算三個(gè)不同網(wǎng)格尺寸(0.4 mm，0.8 mm和1.6 mm)下的結(jié)構(gòu)位移-荷載曲線(xiàn),如圖2所示?？梢钥闯?，網(wǎng)格尺寸為1.6 mm的系統(tǒng)響應(yīng)曲線(xiàn)與實(shí)驗(yàn)的位移-荷載曲線(xiàn)吻合，拉伸試件斷裂時(shí)的位移約為9.0 mm。同時(shí)可以觀察到明顯的網(wǎng)格依賴(lài)性，隨著網(wǎng)格尺寸的減小(0.8 mm和0.4 mm)，局部彈塑性損傷模型計(jì)算的系統(tǒng)響應(yīng)曲線(xiàn)下降時(shí)刻越早，斷裂位移相應(yīng)越小(分別對(duì)應(yīng)為8.5 mm和8.2 mm附近)。隨著網(wǎng)格尺寸以1/2的比率遞減，相鄰網(wǎng)格尺寸所對(duì)應(yīng)的斷裂位移間的偏差也越來(lái)越小(約為0.5 mm和0.3 mm)。

圖1 DP1000單向拉伸實(shí)驗(yàn)樣件

Fig.1 Geometry of the DP1000 uniaxial tension specimen

圖2 拉伸實(shí)驗(yàn)局部損傷模型下的位移-荷載曲線(xiàn)

Fig.2 Displacement-force curves of the local damage model of the uniaxial tensile test

圖4分別給出了網(wǎng)格尺寸為0.4 mm的結(jié)構(gòu)拉伸至剛出現(xiàn)破壞單元時(shí)，局部和微態(tài)損傷模型計(jì)算的局部損傷分布圖?？梢钥闯觯B(tài)損傷模型的損傷集中帶比局部損傷模型預(yù)測(cè)的損傷集中帶更寬，損傷分布更加平滑。圖5的曲線(xiàn)分布趨勢(shì)與之類(lèi)似。自峰值處(約18.0 mm)向兩側(cè)延伸方向，與局部損傷模型計(jì)算的有明顯集中現(xiàn)象的局部損傷分布相比，在微態(tài)損傷模型中的局部損傷變量，由于微態(tài)損傷的正則化作用，其分布更加平滑且均有不同程度的抬升(10 mm～25 mm)。

圖3 拉伸實(shí)驗(yàn)微態(tài)損傷模型下的位移-載荷曲線(xiàn)

Fig.3 Displacement-force curves of the micromorphic damage model of the uniaxial tensile test

圖4 0.4 mm單元的損傷分布

Fig.4 Distributions of damage for mesh size is 0.4 mm

4.2 十字形零件的沖壓實(shí)驗(yàn)

圖6為十字形零件的沖壓成形模具的幾何信息，材料為1.5 mm厚的DP1000對(duì)偶雙相鋼，壓邊力為400 kN，摩擦系數(shù)約為0.15?？紤]微裂紋在壓縮狀態(tài)下的閉合效應(yīng)，式(11)的局部損傷能量釋放率修改為式(30)，并取h的值為0.25。

(30)

圖5 0.4 mm單元的損傷變量沿軸向的分布

Fig.5 Distribution of damage along thexaxial for mesh size is 0.4 mm

圖6 十字形零件模具的幾何信息[22]

Fig.6 Geometric information of the cross-section drawing part[22]

圖7 十字形零件成形局部損傷模型的位移-載荷曲線(xiàn)

Fig.7 Displacement-force curves of the local damage model of the cross-section drawing part

圖8 十字形零件成形微態(tài)損傷模型的位移-載荷曲線(xiàn)

Fig.8 Displacement-force curves of the micromorphic damage model of the cross-section drawing part

5 結(jié) 論

本文根據(jù)Forest的微態(tài)方法，提出引入微態(tài)損傷因子及其一階梯度，以考慮應(yīng)變集中的非局部效應(yīng)，建立了新的適用于有限變形框架的熱力學(xué)一致性的微態(tài)損傷彈塑性本構(gòu)模型。通過(guò)在擴(kuò)展的亥姆霍茲自由能中定義局部損傷與微態(tài)損傷的耦合方式，進(jìn)而得到了包含材料內(nèi)部特征尺度的類(lèi)似于隱式梯度理論中補(bǔ)充的亥姆霍茲方程，同時(shí)獲得了新的包含局部損傷和微態(tài)損傷的損傷能量釋放率。由微態(tài)損傷變量對(duì)損傷能量釋放率進(jìn)行正則化作用，進(jìn)而影響了局部損傷的演化過(guò)程，最終達(dá)到了一致性的系統(tǒng)位移-載荷響應(yīng)。對(duì)DP1000材料的單向拉伸和十字形零件的沖壓形實(shí)驗(yàn)進(jìn)行數(shù)值分析，結(jié)果表明，該微態(tài)損傷模型可以有效地避免材料應(yīng)變軟化行為導(dǎo)致的數(shù)值響應(yīng)依賴(lài)于網(wǎng)格的問(wèn)題，并弱化損傷變量的局部集中現(xiàn)象。并且該模型可進(jìn)一步擴(kuò)展包含混合硬化和塑性應(yīng)變等非局部效應(yīng)，以完善對(duì)材料應(yīng)變軟化行為的描述。