基于單向拉伸仿真對各向異性薄板本構(gòu)修正的研究

張東凱，李曉廣，余亞克

(河鋼集團邯鋼技術(shù)中心，河北 邯鄲 056038)

傳統(tǒng)板料沖壓成形主要依靠設(shè)計者的經(jīng)驗，即通過不斷修模、試模、調(diào)整間距及壓邊力等試錯方式，來解決成形中出現(xiàn)起皺、破裂等問題。隨著汽車工業(yè)的發(fā)展，板料沖壓成形工藝趨于復雜，傳統(tǒng)方法導致零件生產(chǎn)周期長、效率低，已經(jīng)不能適應成形工藝的發(fā)展。而計算機技術(shù)和有限元方法的發(fā)展給板料成形帶來了新方法，通過計算機仿真模擬沖壓工藝、指導沖壓設(shè)計、減少試制次數(shù)、縮短生產(chǎn)周期，從而降低制造成本，提高生產(chǎn)效率。要提高模擬精度，使仿真模擬更加接近實際成形過程，必須選擇合理的本構(gòu)模型[1]，如屈服準則、硬化模型、損傷斷裂模型等。軋制板料帶有各向異性行為，Hill48屈服準則可以準確預測鋼材的屈服軌跡，胡成武[2]利用Hill48平面應力方程和Hill塑性應變比建立了預測任意方向屈服方程，預測DC04和DC05板料在0°、45°、90°的Hill48屈服值，并將結(jié)果與試驗數(shù)據(jù)對比，對比結(jié)果十分吻合。牛超[3]利用Hill48屈服函數(shù)和Swift等向強化本構(gòu)數(shù)值模擬了0°、45°、90°真應力-應變曲線。趙坤民[4]利用ABAQUS/STANDARD進行拉伸仿真分析，得到模擬頸縮后各向異性塑性變形的合適單元。ManSoo Joun等[5]使用剛體-塑性有限元方法分析了拉伸試驗結(jié)果，迭代后的應力-應變曲線精確地預測了頸縮點，減小了拉伸試驗和分析結(jié)果的之間的載荷差異。王友德[6]將單軸拉伸試驗的載荷-位移曲線作為目標，通過試驗-數(shù)值耦合方法(NE法)，不斷迭代修正輸入的應力-應變曲線，直至載荷-位移曲線達到要求，即得到真實應力-應變曲線，并比較了NE法較傳統(tǒng)理論方法的優(yōu)點。王玉凱[7]利用ABAQUS/EXPLICIT對Q235冷軋鋼進行拉伸斷裂數(shù)值模擬，驗證了Johnson-Cook本構(gòu)模型對鋼材拉伸性能的合理性。陳輝[8]應用TF方法對單軸本構(gòu)模型進行有限元迭代計算使得模擬計算結(jié)果與圓環(huán)試樣壓縮試驗結(jié)果滿足一致性，從而得到304材料的本構(gòu)關(guān)系參數(shù)。姚迪[9]基于TF方法模擬的載荷-位移曲線、漏斗根部直徑-位移曲線和漏斗變形輪廓線等均與試驗結(jié)果吻合良好，試樣表面模擬應變也與DIC測試結(jié)果吻合。蔡力勛[10]基于單軸拉伸真應力-真應變曲線，采用二分法并結(jié)合有限元模擬迭代載荷-位移曲線，得到真實反映材料拉伸全程變形特征的等效真應力-真應變曲線。祁爽[11]基于有限元輔助方法得到了四種工程材料的等效全程單軸本構(gòu)關(guān)系，同時分析了應力三軸度。彭云強[12]基于ANSYS APDL語言和有限元輔助方法獲得了某型鋼的本構(gòu)關(guān)系，并借助ABAQUS中GTN損傷模型預測了裂紋擴展。以上仿真采用的是各向同性屈服準則，未分析各向異性屈服準則的影響。

本文對DC04進行了單向拉伸試驗，得到載荷位移曲線和真應力-塑性應變曲線。基于有限元輔助方法，在ABAQUS/EXPLICIT中合理簡化拉伸模型，耦合邊界條件，引入Hill48各向異性屈服準則，采用動力顯示算法和適當質(zhì)量縮放進行準靜態(tài)單向拉伸仿真，迭代后的載荷-位移曲線與試驗結(jié)果吻合，此時得到DC04板材的合理塑性本構(gòu)，可作為成形仿真的本構(gòu)輸入。最后對比分析了Hill48屈服準則和Mises屈服準則下的仿真結(jié)果，結(jié)果顯示Hill48屈服準則狀態(tài)下的載荷位移曲線與試驗結(jié)果符合更好。

1 試驗材料及結(jié)果

拉伸試驗材料為1 mm的DC04冷軋鋼板，按照GB/T 228.1-2010[13]切取與軋制方向成0°、45°、90°的標準試樣，試樣尺寸見圖1。拉伸試驗在150 kN的電子力學拉伸機上進行，根據(jù)GB/T 5027-2007[14],夾頭移動速度定為2 mm/min。拉伸斷裂試件見圖2，由試驗所得到的載荷位移曲線見圖3，各個方向屈服應力σα和塑性應變比rα見表1。

圖1 單向拉伸試樣形狀及尺寸 (單位：mm)Fig.1 Shape and geometric dimension of uniaxial tensile specimen(unit: mm)

圖2 拉伸后的DC04試件(從上到下依次為0°、45°、90°)Fig.2 The specimens of DC04 after stretching(From top to bottom, 0°, 45°, 90°)

圖3 單向拉伸試驗載荷位移曲線Fig.3 Load-displacement curve of unidirectional tensile test

基于體積不可壓縮理論，將工程應力-工程應變按式(1)換算為真應力-應變曲線，去除彈性應變后得到真應力-塑性應變曲線，見圖4，將其作為單向拉伸仿真初始本構(gòu)輸入。

圖4 真應力-塑性應變曲線 Fig.4 True stress and plastic strain

(1)

式中，σT為真應力，εT為真應變，σE為工程應力，εE為工程應變。

2 有限元數(shù)值模擬

2.1 有限元輔助方法

圖5 迭代修正方法Fig.5 Iterative correction method

2.2 幾何及邊界條件

由于考察的是拉伸件的中部區(qū)域，由圣維南原理，去除了夾持段和圓角過渡段，留下標距段，尺寸為1 mm×20 mm×80 mm，見圖6。

圖6 單向拉伸幾何模型及邊界條件Fig.6 Unidirectional stretching geometric model and boundary conditions

為了模擬邊界條件，將模型兩邊分別設(shè)置兩個參考點RP1(右端)，RP2(左端)，位置為過渡圓角與直段的相交中心，按圖中坐標為(-60,0)，(60,0)，將RP1與右端面耦合，耦合位移轉(zhuǎn)角，同理應用于RP2。

2.3 材料模型

采用動力顯示算法，需要輸入DC04的密度7.85×10-9t/mm3，彈性模量210 000 MPa，泊松比0.3。

塑性本構(gòu)首次輸入試驗真應力-塑性應變曲線，這里以真應力和塑性應變表格形式，將擬合曲線輸入到ABAQUS塑性屬性設(shè)置中，需要注意的是真應力為遞增形式，在循環(huán)迭代本構(gòu)后，重新輸入塑性表格數(shù)據(jù)時需去除遞減真應力段。

由ABAQUS手冊，對于板料(平面應力狀態(tài))，Hill48各向異性屈服準則中的應變比R11，R22，R33，R12，R13，R23用各向異性系數(shù)r值表示為

R11=R13=R23=1

(2)

(3)

(4)

(5)

按公式2～5，DC04的R11，R22，R33，R12，R13，R23值見表2。由于引入各向異性，需要定義材料主方向，這里定義幾何模型縱向為材料主方向。

表2 Hill48屈服準則的應變比參數(shù)

2.4 分析步及網(wǎng)格設(shè)置

由于采用準靜態(tài)單向拉伸仿真，所以分析步采用動力顯示分析，時間間隔為0.1 s，對幾何模型進行適當?shù)馁|(zhì)量縮放，后續(xù)分析表明質(zhì)量縮放因子取1000，滿足準靜態(tài)仿真要求。

由于要循環(huán)修正真應力和塑性應變，所以在ABAQUS場輸出中需勾選MISES，PEEQ作為修正結(jié)果，考慮計算效率，時間均勻間隔200。

為了分析計算結(jié)果，得到載荷位移曲線，需要在歷史輸出中勾選參考點RP1的支反力RF1及縱向位移U1，同理，時間均勻間隔200。

整體網(wǎng)格尺寸為1，厚度方向布置了4個種子，將厚度方向分割了3段，中間部分對縱向邊線增加種子個數(shù)，實現(xiàn)網(wǎng)格細化，見圖7。在仿真中總結(jié)經(jīng)驗，調(diào)節(jié)單元類型，其中各向異性采用非協(xié)調(diào)模式。

圖7 網(wǎng)格模型Fig.7 Grid model

3 仿真結(jié)果

3.1 準靜態(tài)分析

在ABAQUS/EXPLICIT中，能量平衡方程可以表示為

constant=ETOT=EI+EV+EFD+EKE-EW

(6)

其中，EI為內(nèi)能(包括彈性應變能、塑性應變能和與沙漏控制相關(guān)的偽能)，EV為粘性機制耗散的能量，EFD為摩擦耗散的能量，EKE為動能，EW為外力功，ETOT為系統(tǒng)的總能量。

由于不存在粘彈性材料、離散的減震器及材料阻尼，所以粘性機制耗散的能量EV很小。單向拉伸仿真沒有摩擦，故EFD為零。只剩下EI+EKE-EW，因此對于單軸拉伸仿真，若變形材料的動能EKE占內(nèi)能EI的比例很小(一般為1%～5%)，那么拉伸試件的外力功EW約等于試件的內(nèi)能EI，那么仿真就是準靜態(tài)的。

從圖8中可以看出，整個模型的動能ALLKE(紅色曲線)對時間的積分只占整個模型內(nèi)能ALLIE(藍色曲線,基本與X軸重合)對時間積分的很小部分，所以是準靜態(tài)的。

圖8 準靜態(tài)拉伸能量歷程 Fig.8 Quasi-static tensile energy history

3.2 Hill48屈服準則下的載荷位移曲線分析

單向拉伸仿真得到等效Mises應力、等效塑性應變PEEQ云圖見圖9、10，局部出現(xiàn)頸縮，頸縮區(qū)域中間單元的等效Mises應力與等效塑性應變PEEQ即為新的塑性本構(gòu)。按照圖5的修正方法，對材料本構(gòu)進行了兩次迭代，每次迭代仿真后的載荷位移曲線及與試驗數(shù)據(jù)的標準差見圖11。

圖9 等效應力云圖(單位：Mpa) Fig.9 Equivalent stress cloud(unit: MPa)

從圖11中可以看出，將單向拉伸試驗數(shù)據(jù)以真應力-塑性應變的形式輸入，經(jīng)過單向拉伸仿真后的載荷位移曲線與試驗數(shù)據(jù)的標準差為959.01，差異主要在損傷斷裂階段，說明單向拉伸試驗中頸縮部分的本構(gòu)數(shù)據(jù)無法準確描述損傷斷裂段，因此需要對本構(gòu)進行修正。采用迭代修正方法對塑性本構(gòu)進行第一次修正后的仿真載荷位移曲線與試驗數(shù)據(jù)標準差為215.82，誤差下降約3.4倍，經(jīng)過仿真后的載荷位移曲線在彈性和塑性強化階段與試驗數(shù)據(jù)吻合，損傷斷裂段開始接近試驗曲線，損傷斷裂的本構(gòu)得到有效修正。經(jīng)過第二次迭代修正后的仿真載荷位移曲線與試驗數(shù)據(jù)的標準差為208.89，與試驗數(shù)據(jù)的差距進一步變小，標準差相對于第一次修正下降明顯變緩，此時的真應力-塑性應變曲線(見表3)即為目標本構(gòu)，可應用于汽車零部件的沖壓成形仿真。

圖10 等效塑性應變云圖 Fig.10 Equivalent plastic strain cloud

圖11 迭代仿真后的載荷位移曲線 Fig.11 Load-displacement curve after iterative simulation based

表3 迭代后的真應力-塑性應變本構(gòu)

3.3 兩種屈服準則下的載荷位移曲線對比

從圖12中可以看出，Hill48屈服準則狀態(tài)下的載荷位移曲線和Mises屈服準則狀態(tài)下的結(jié)果在彈性階段和塑性強化階段一致，在損傷斷裂段出現(xiàn)差異,且與試驗數(shù)據(jù)的標準差分別為959.01和527.87，說明屈服準則對仿真損傷斷裂段影響較大。

圖12 不同屈服準則下的仿真結(jié)果Fig.12 Simulation results under different yield criteria

從圖13中可以看出，塑性本構(gòu)經(jīng)過第一次迭代后，Hill48屈服準則狀態(tài)下的載荷位移曲線和Mises屈服準則狀態(tài)下的仿真結(jié)果與試驗數(shù)據(jù)的標準差分別為215.82和330.00，Hill48屈服準則下的仿真結(jié)果標準差更小。

圖13 第一次迭代仿真后的結(jié)果Fig.13 Results after the first iteration of the simulation

圖14 第二次迭代仿真后的結(jié)果Fig.14 Results after the second iteration of the simulation

4 結(jié)論

1)按照單向拉伸試驗標準，得到DC04(厚度1 mm)薄板的0°、45°、90°載荷位移曲線，真應力-塑性應變曲線，各向屈服應力和塑性應變比。

2)合理簡化了單向拉伸仿真模型及顯示動力計算，仿真結(jié)果表明模型有效，且經(jīng)分析達到了準靜態(tài)要求。

4)通過對比Hill48屈服準則和Mises屈服準狀態(tài)下的仿真結(jié)果，發(fā)現(xiàn)不同的屈服準則對損傷斷裂階段影響不同，Hill48各向異性屈服下的載荷位移曲線與試驗結(jié)果的一致性較各向同性Mises屈服準則更好。