GTN損傷模型的算法研究及試驗(yàn)驗(yàn)證

王效貴,李曉葉

(浙江工業(yè)大學(xué) 機(jī)械工程學(xué)院，浙江 杭州 310014)

GTN損傷模型的算法研究及試驗(yàn)驗(yàn)證

王效貴,李曉葉

(浙江工業(yè)大學(xué) 機(jī)械工程學(xué)院，浙江 杭州 310014)

摘要：基于GTN細(xì)觀損傷本構(gòu)理論，引入了塑性極限載荷模型，將完全隱式積分算法與隱式有限元計(jì)算相結(jié)合，建立了相應(yīng)的數(shù)值積分算法.通過編寫用戶自定義材料子程序(UMAT)，實(shí)現(xiàn)了修正的GTN模型在有限元軟件ABAQUS中的應(yīng)用.通過缺口圓棒試樣損傷與斷裂的數(shù)值模擬，對(duì)修正的GTN模型進(jìn)行了驗(yàn)證，并分析了缺口圓棒試樣的裂紋萌生點(diǎn)位置，裂紋擴(kuò)展情況及應(yīng)力三軸度的變化.數(shù)值預(yù)測(cè)與試驗(yàn)結(jié)果吻合較好，表明修正的GTN模型可以有效的描述材料的損傷和斷裂現(xiàn)象.

關(guān)鍵詞：GTN模型；數(shù)值算法；損傷斷裂；有限元法

Algorithm study and experiment validation of GTN damage model

WANG Xiaogui, LI Xiaoye

(College of Mechanical Engineering, Zhejiang University of Technology, Hangzhou 310014, China)

Abstract:The corresponding numerical integral algorithm is established by adding the plastic limit load model into the GTN mesoscopic damage constitutive model and combining the fully implicit integration algorithm and the implicit finite element. It is implemented in the finite element software ABAQUS by using the user material subroutine (UMAT). The modified GTN model is verified by simulating the damage and fracture of notched round bar. The crack initiation site and propagation path, and the change of stress triaxiality were analyzed. The numerical predications are in good agreement with the experimental observations. It also shows that the modified GTN mesoscopic damage model can effectively reveal the damage behavior and predict the fracture of the material.

Key words:GTN model; numerical algorithm; damage and failure; finite element method

損傷和斷裂是工程材料的常見失效模式.材料中的夾雜物或二相粒子與基體的脫粘或自身開裂，將會(huì)引起孔洞的形核、長(zhǎng)大及聚合，最終形成一定尺度的宏觀裂紋，這是導(dǎo)致材料失效的根本原因.近年來，很多學(xué)者采用不同方法對(duì)損傷斷裂問題展開研究.劉俊新[1]等采用三維重建技術(shù)研究了泥頁(yè)巖蓋層在載荷作用下的裂紋擴(kuò)展與演化規(guī)律；蔡增伸[2]等通過力控制法對(duì)金剛石膜的斷裂韌性進(jìn)行了研究.蘇彬[3]等自行設(shè)計(jì)實(shí)驗(yàn)裝置研究了循環(huán)載荷作用下的微動(dòng)疲勞特性和裂紋位置.

目前應(yīng)用最為廣泛的損傷模型是GTN模型，它是一種以孔洞體積分?jǐn)?shù)為損傷變量反映材料內(nèi)部微小缺陷的損傷模型，基于這一模型國(guó)內(nèi)外學(xué)者取得了諸多研究成果.Vadillo[4]等通過數(shù)值算法研究了GT模型中的損傷參數(shù)q1和q2，發(fā)現(xiàn)這兩個(gè)參數(shù)并不是常數(shù)，而是與應(yīng)力三軸度和初始孔隙率相關(guān)的變量；邢佶慧[5]等通過數(shù)值模擬的方法，分析了初始孔洞體積分?jǐn)?shù)，孔洞形核體積分?jǐn)?shù)及材料失效時(shí)的孔洞破壞體積分?jǐn)?shù)對(duì)斷裂失效的影響.由于原Gurson模型沒有考慮不連續(xù)孔洞的影響，不能預(yù)測(cè)孔洞之間的頸縮現(xiàn)象.因此，筆者基于GTN細(xì)觀損傷理論，引入塑性極限載荷模型，采用Ramberg-Osgood硬化準(zhǔn)則，從數(shù)值積分算法著手，模擬高導(dǎo)無(wú)氧銅(OFHC)材料的損傷演化和斷裂失效過程，并與試驗(yàn)作對(duì)比，驗(yàn)證其可行性.

1GTN模型的算法研究及數(shù)值實(shí)現(xiàn)

GTN模型的屈服函數(shù)[6]表示為

1-(q1f*)2

(1)

式中：σeq為宏觀Mises等效應(yīng)力；σm為宏觀靜水應(yīng)力；σY為流動(dòng)應(yīng)力；q1和q2為損傷參數(shù)；f*為孔洞體積分?jǐn)?shù)函數(shù)，定義為

(2)

式中：f為孔洞體積分?jǐn)?shù)；fc為開始聚合時(shí)的臨界孔洞體積分?jǐn)?shù)；ff為材料斷裂時(shí)的臨界孔洞體積分?jǐn)?shù).

由于原Gurson模型不能預(yù)測(cè)孔洞之間的頸縮現(xiàn)象，因此在GTN模型中引入塑性極限載荷模型來描述這一現(xiàn)象.該模型給出的軸對(duì)稱球形孔洞聚合公式[7-8]為

(3)

式中：Rr為球形孔洞的當(dāng)前半徑；L為單胞的當(dāng)前長(zhǎng)度；σ1為最大主應(yīng)力.

OFHC在塑性屈服后，表現(xiàn)出輕微的硬化現(xiàn)象，為了反映這個(gè)現(xiàn)象，采用Ramberg-Osgood硬化準(zhǔn)則，其曲線關(guān)系為

(4)

式中：σ0為初始屈服應(yīng)力；ε0為初始屈服應(yīng)變；k為硬化系數(shù)；ρ為硬化指數(shù).

根據(jù)Aravas[9]提出的向后Euler完全隱式積分算法，更新應(yīng)力應(yīng)變關(guān)系.n表示增量步，為方便書寫，省略下標(biāo)n+1.假設(shè)給定的應(yīng)變?cè)隽咳珵閺椥詰?yīng)變?cè)隽?，?jì)算彈性試應(yīng)力

σtrial=σn+C：Δε

(5)

式中：C為四階彈性模量；Δε為總的應(yīng)變?cè)隽?

根據(jù)關(guān)聯(lián)流動(dòng)準(zhǔn)則，更新n+1增量步時(shí)的宏觀塑性應(yīng)變?cè)隽繛?/p>

(6)

(7)

根據(jù)Hook定律，結(jié)合式(5，6)，得到n+1增量步時(shí)的宏觀柯西應(yīng)力為

σ=σtrial-KΔDmI-2GΔDeqN

(8)

式中：G為剪切模量；K為體積模量.將式(8)分別映射到I和N方向，得到當(dāng)前時(shí)刻的宏觀靜水應(yīng)力和宏觀等效應(yīng)力為

(9)

因此，宏觀柯西應(yīng)力又可表示為

(10)

Δf=(1-f)ΔDm+

(11)

(12)

將式(4，9，11，12)代入式(1)，得到n+1增量步時(shí)的屈服函數(shù)為

F(ΔDeq,ΔDm)=0

(13)

(14)

根據(jù)式(7，13，14)推導(dǎo)一致切線剛度矩陣T，得到

(15)

式中：Mijkl的表達(dá)式見文獻(xiàn)[9]，Jijkl=(δikδjl+δilδjk)/2.

根據(jù)以上推導(dǎo)的離散公式，下面給出修正的GTN模型在用戶子程序(UMAT)中的實(shí)現(xiàn)流程圖，如圖1所示.

圖1　流程圖Fig.1　Flow chart

2單軸拉伸試驗(yàn)

缺口圓棒試樣的單軸拉伸試驗(yàn)是研究韌性金屬材料損傷斷裂的最好方法，因此，采用如圖2所示的高導(dǎo)無(wú)氧銅光滑和缺口圓棒試樣，進(jìn)行單軸拉伸試驗(yàn).光滑試樣標(biāo)距段的直徑為10mm.缺口試樣有三個(gè)重要幾何參數(shù)，分別為缺口寬度d、缺口半徑r和凈斷面直徑D，實(shí)測(cè)數(shù)據(jù)如表1所示.采用RG4100微機(jī)控制電子萬(wàn)能試驗(yàn)機(jī)，在室溫狀態(tài)下完成單軸拉伸試驗(yàn).試驗(yàn)采用位移控制，引伸計(jì)的標(biāo)距為20mm.光滑試樣的位移加載速率為0.3mm/min.基于不同缺口尺寸，試樣φ6r1，φ6r2，φ8r1，φ8r2的位移加載速率分別為0.05，0.08，0.10，0.15mm/min.試驗(yàn)完畢后，測(cè)量缺口部位的d和D，結(jié)果列于表1.由表1中的斷面收縮率可以看出：凈斷面直徑越大，斷面收縮率越大；同一種凈斷面直徑，缺口半徑越大，斷裂收縮率越大.

圖2　單軸拉伸試樣Fig.2　Uniaxial tensile specimens

試樣編號(hào)實(shí)驗(yàn)前缺口半徑r/mm缺口寬度d/mm凈斷面直徑D/mm試驗(yàn)后缺口寬度d'/mm凈斷面直徑D'/mm斷面收縮率/%?6r1-11.02.106.053.263.8060.55?6r1-21.02.085.953.091)4.1052.52?6r1-31.02.126.073.061)4.0655.26?6r2-12.03.845.955.793.6063.39?6r2-22.03.865.935.733.6362.53?6r2-32.03.835.925.521)3.8058.80?8r1-11.02.157.924.594.4568.43?8r1-21.02.107.944.564.4768.31?8r1-31.02.137.954.524.5067.76?8r2-12.03.237.935.744.3070.60?6r2-22.03.257.945.774.2771.08?6r2-32.03.207.925.794.3470.00

注：1) 為在單軸拉伸試驗(yàn)中材料沒有拉斷.

由光滑試樣單軸拉伸試驗(yàn)得到的真實(shí)應(yīng)力—應(yīng)變曲線，如圖3所示.由彈性階段的斜率，得到OFHC的彈性模量為E=108 GPa.基于國(guó)標(biāo)GB/T 228.1—2010，確定初始屈服應(yīng)變?yōu)棣?=0.026 85，且滿足σ0=Eε0，屈服極限為σy=290 MPa，強(qiáng)度極限為σu=314 MPa.OFHC在塑性屈服后，表現(xiàn)出輕微的硬化現(xiàn)象.采用公式(4)給出的Ramberg-Os-good硬化準(zhǔn)則，曲線擬合后得到硬化系數(shù)為k=0.75和硬化指數(shù)為ρ=27.

圖3　光滑圓棒試樣的真實(shí)應(yīng)力—應(yīng)變曲線Fig.3　True stress-strain curves of smooth round bar specimens

由缺口試樣單軸拉伸試驗(yàn)得到的載荷—位移曲線如圖4所示.對(duì)每種缺口試樣進(jìn)行3次重復(fù)試驗(yàn)，獲得的載荷—位移曲線具有很好的一致性.試驗(yàn)曲線分成三個(gè)階段：第一，試驗(yàn)開始到載荷最大值，這一期間載荷快速增大，而材料的變形程度非常??；第二，載荷最大值到裂紋萌生點(diǎn)(拐點(diǎn))，當(dāng)達(dá)到最大載荷時(shí)孔洞開始聚合，材料承載能力逐漸減弱，同時(shí)缺口半徑處出現(xiàn)頸縮現(xiàn)象，縱向變形快速增大；第三，裂紋萌生點(diǎn)到材料斷裂失效，隨著裂紋的萌生和擴(kuò)展，損傷演化程度加快，載荷急劇降低，材料承載能力下降，最終導(dǎo)致斷裂破壞.圖4可以看出：缺口試樣的凈斷面直徑越小，材料的斷裂失效出現(xiàn)的越早；同一種凈斷面直徑，缺口半徑越小，斷裂破壞產(chǎn)生的越快；凈斷面直徑越大，載荷—位移曲線中峰值越高.

圖4　缺口試樣的載荷—位移曲線Fig.4　Load-Displacement curves of notched specimens

3算例及討論

3.1有限元模型

建四分之一軸對(duì)稱有限元模型，如圖5所示.邊界條件的施加分為3部分：在左端面，施加軸向位移載荷，試樣φ6r1，φ6r2，φ8r1，φ8r2的位移大小分別為0.6，0.75，1.5，1.5 mm；在右端面，施加y向?qū)ΨQ位移約束；在x=0對(duì)稱面上，施加x向?qū)ΨQ位移約束.為了能夠很好地模擬缺口根部區(qū)域較大的塑性應(yīng)力應(yīng)變梯度，將該區(qū)域的網(wǎng)格細(xì)化，最小網(wǎng)格尺寸為0.1 mm.為了提高計(jì)算速度，遠(yuǎn)離缺口區(qū)域的網(wǎng)格尺寸為0.2 mm.采用嵌入GTN模型的有限元軟件ABAQUS，選用八節(jié)點(diǎn)軸對(duì)稱縮減積分單元(CAX8R)，進(jìn)行數(shù)值模擬斷裂預(yù)測(cè).

圖5　缺口圓棒試樣的有限元模型Fig.5　FEA model of notched specimens

3.2GTN模型的參數(shù)確定

由于GTN損傷模型中的體胞是假設(shè)的，因此f0的取值需根據(jù)模擬結(jié)果和試驗(yàn)數(shù)據(jù)對(duì)比而確定.以試樣φ6r1-1為研究對(duì)象，當(dāng)給定的f0能夠使數(shù)值模擬結(jié)果和試驗(yàn)數(shù)據(jù)基本相吻合時(shí)，則認(rèn)為此時(shí)的f0即為要確定的初始孔隙率，最終確定f0=0.001.文獻(xiàn)[4]給出了損傷參數(shù)q1和q2隨應(yīng)力三軸度的變化圖，依據(jù)給定的關(guān)系圖，可以初步確定損傷參數(shù)q1和q2，由于不同材料的實(shí)際損傷情況并不一致，因此根據(jù)數(shù)值預(yù)測(cè)結(jié)果和試驗(yàn)數(shù)據(jù)進(jìn)一步修正，最終得到q1=1.2，q2=0.8.根據(jù)本試驗(yàn)的實(shí)際情況，結(jié)合文獻(xiàn)[5]中介紹的方法，最終確定fN=0.004 5，ff=0.4.

3.3結(jié)果與討論

根據(jù)確定的損傷參數(shù)，將修正的GTN模型嵌入ABAQUS軟件進(jìn)行斷裂預(yù)測(cè)，分析缺口圓棒試樣裂紋萌生點(diǎn)位置和擴(kuò)展情況，以及應(yīng)力三軸度的變化.

3.3.1裂紋萌生點(diǎn)位置及擴(kuò)展路線的預(yù)測(cè)

圖4展示了缺口試樣裂紋萌生點(diǎn)和裂紋擴(kuò)展路線，圖4中的黑點(diǎn)代表裂紋萌生點(diǎn)位置.可以清楚的看到：裂紋萌生最早產(chǎn)生于缺口試樣的中心位置附近，隨后向缺口試樣的尖端方向逐漸發(fā)展，隨著塑性變形的增大，裂紋快速擴(kuò)展并貫通缺口根部，最終造成材料的斷裂失效.這是因?yàn)?，由于塑性變形的增大，孔洞體積分?jǐn)?shù)也隨之增加，當(dāng)達(dá)到塑性極限載荷時(shí)，微孔洞開始聚合，孔洞體積分?jǐn)?shù)逐漸由f變成f*，并急劇增大，當(dāng)f*增加到一定程度時(shí)，缺口試樣的中間部位裂紋開始萌生.隨著裂紋的快速長(zhǎng)大和擴(kuò)展，載荷急劇下降，材料的承載能力降低，當(dāng)f*增加到材料斷裂的臨界值時(shí)損傷失效產(chǎn)生.

從圖4中的試驗(yàn)和數(shù)值模擬曲線看到：裂紋萌生點(diǎn)位置的數(shù)值預(yù)測(cè)與試驗(yàn)結(jié)果較為接近，裂紋擴(kuò)展路線的預(yù)測(cè)結(jié)果與試驗(yàn)曲線趨勢(shì)基本一致.表明修正的GTN模型能夠模擬缺口圓棒試樣的斷裂失效特征，驗(yàn)證了該模型的有效性，該方法的合理性和可行性.

3.3.2應(yīng)力三軸度的變化

應(yīng)力三軸度定義為平均應(yīng)力與等效應(yīng)力的比值，它是影響韌性材料斷裂的一個(gè)重要因素，通常根據(jù)試樣最小橫截面來評(píng)估.圖6給出最大載荷條件下應(yīng)力三軸度的變化過程，圖7給出裂紋萌生點(diǎn)附近應(yīng)力三軸度的變化過程.

以上四幅圖顯示，凈斷面直徑相同，缺口半徑越小，最小橫截面上的應(yīng)力三軸度越大；應(yīng)力三軸度通常在缺口試樣的中間部位最高，在自由邊附近最低.由3.3.1節(jié)知道：裂紋萌生最早發(fā)生于缺口試樣中心部位附近，隨著塑性變形的增大，加劇了缺口試樣中心位置的裂紋萌生擴(kuò)展，從而導(dǎo)致缺口根部的孔洞體積分?jǐn)?shù)快速增大.由GTN模型的屈服函數(shù)可以知道，當(dāng)孔洞體積分?jǐn)?shù)增大到一定程度時(shí)，平均應(yīng)力遠(yuǎn)大于等效應(yīng)力，因此出現(xiàn)了上圖所述規(guī)律.比較圖6,7看到：同一種缺口結(jié)構(gòu)，圖7中的應(yīng)力三軸度要比圖6中的偏大.這是因?yàn)?，?dāng)達(dá)到塑性極限載荷時(shí)微孔洞剛開始聚合，裂紋萌生并沒有完全開始，孔洞體積分?jǐn)?shù)相對(duì)較小，此時(shí)等效應(yīng)力所起的作用比平均應(yīng)力大，應(yīng)力三軸度較小.而在裂紋萌生點(diǎn)附近，由于塑性變形的增大，加劇了裂紋萌生和快速擴(kuò)展，導(dǎo)致此時(shí)的孔洞體積分?jǐn)?shù)比最大載荷工況下更大，平均應(yīng)力的影響遠(yuǎn)遠(yuǎn)大于等效應(yīng)力，所以應(yīng)力三軸度較大.同時(shí)還可以發(fā)現(xiàn)，同一種缺口半徑，凈斷面直徑越小，應(yīng)力三軸度越大.這是因?yàn)閮魯嗝嬷睆皆叫?，此處的?yīng)力集中對(duì)損傷斷裂的影響更加敏感而造成.

圖6　最大載荷下應(yīng)力三軸度沿最小截面的變化Fig.6　Variations in stress triaxiality along minimum section corresponding to maximum load

圖7　裂紋萌生下應(yīng)力三軸度沿最小截面的變化Fig.7　Variations in stress triaxiality along minimum section corresponding to crack initation

圖4中的試驗(yàn)曲線顯示：試樣從早到晚的斷裂情況依次是φ6r1，φ6r2，φ8r1，φ8r2.觀察圖6，7發(fā)現(xiàn)：應(yīng)力三軸度由大到小的順序是φ6r1，φ6r2，φ8r1，φ8r2.而應(yīng)力三軸度越大，孔洞體積分?jǐn)?shù)的增加就越急劇，裂紋萌生和擴(kuò)展速度越快，斷裂破壞產(chǎn)生的就越早.因此GTN模型的模擬結(jié)果和試驗(yàn)結(jié)果吻合良好.

4結(jié)論

采用嵌入修正的GTN模型的有限元軟件ABAQUS進(jìn)行斷裂模擬，得出以下幾點(diǎn)結(jié)論：1) 通過對(duì)缺口圓棒試樣裂紋萌生點(diǎn)位置和擴(kuò)展路徑的分析發(fā)現(xiàn)，裂紋萌生首先產(chǎn)生于缺口試樣的中部，而不是缺口試樣的尖端.一旦裂紋萌生，載荷急劇下降，材料的承載能力快速降低.2) 應(yīng)力三軸度的分析表明，缺口試樣的中間位置應(yīng)力三軸度最大，自由邊附近最低.隨著塑性變形的增加，孔洞體積分?jǐn)?shù)急劇增大，加劇了裂紋的萌生和長(zhǎng)大，從而導(dǎo)致缺口試樣中間部位的塑性應(yīng)變最大，應(yīng)力三軸度也最大，因此斷裂失效最早發(fā)生.3) 缺口試樣的模擬結(jié)果與試驗(yàn)結(jié)果吻合較好，驗(yàn)證了修正的GTN模型具有預(yù)測(cè)材料斷裂行為的可行性.

參考文獻(xiàn)：

[1]劉俊新，楊春和，冒海軍，等．基于CT圖像處理的泥頁(yè)巖裂紋擴(kuò)展與演化研究[J]．浙江工業(yè)大學(xué)學(xué)報(bào)，2015，43(1)：66-72．

[2]蔡增伸，蔣政，王從賢，等．測(cè)定金剛石膜斷裂強(qiáng)度及斷裂韌性的力控制法研究[J]．浙江工業(yè)大學(xué)學(xué)報(bào)，2000，28(3)：216-219．

[3]蘇彬，邢海軍，趙穎娣，等．循環(huán)載荷應(yīng)力比對(duì)微動(dòng)疲勞特性的影響[J]．浙江工業(yè)大學(xué)學(xué)報(bào)，2011，39(4)：445-448．

[4]VADILLO G, FERNANDEZ-SAEZ J. An analysis of Gurson model with parameters dependent on triaxiality based on unitary cell[J]. European Journal of Mechanics. A: Solids，2009，28(3)：417-427．

[5]邢佶慧，史一劍，劉文濤，等．建筑鋼材GTN損傷模型參數(shù)識(shí)別[J]．建筑結(jié)構(gòu)學(xué)報(bào)，2014，35(4)：149-154．

[6]TVERGAARD V, NEEDLEMAN A. Analysis of the cup-cone fracture in a round tensile bar[J]．Acta Metallurgica，1984，32(1)：157-169．

[7]THOMASON P F. Ductile fracture of metals[M]．Oxford：Pergamon Press，1990．

[8]PARDOEN T, HUTCHINSON J W. An extended model for void growth and coalescence[J]．Journal of Mechanics and Physics of Solids，2000，48(12)：2467-2512．

[9]ARAVAS N. On the numerical integration of a class of pressure-dependent plasticity models[J]．International Journal for numerical methods in engineering，1987，24(7)：1397-1416．

[10]CHU C C, NEEDLEMAN A. Void nucleation effects in biaxially stretched sheets[J]．Journal of Engineering Materials and Technology，1980，102(3)：249-256．

(責(zé)任編輯：陳石平)

文章編號(hào)：1006-4303(2015)06-0660-06

中圖分類號(hào)：TG146.1

文獻(xiàn)標(biāo)志碼：A

作者簡(jiǎn)介：王效貴(1973—)，男，山東青州人，教授，研究方向?yàn)榻Y(jié)構(gòu)完整性，E-mail：hpcwxg@zjut.edu.cn.

基金項(xiàng)目：國(guó)家自然科學(xué)基金資助項(xiàng)目(51175469)

收稿日期：2015-05-04