屈服準(zhǔn)則和硬化模型對(duì)5052鋁板回彈仿真的影響

余海燕，沈嘉怡，王 友

(同濟(jì)大學(xué) 汽車(chē)學(xué)院，上海 201804)

汽車(chē)質(zhì)量每減少10%可使油耗下降6%～8%，汽車(chē)輕量化對(duì)降低油耗非常關(guān)鍵.鋁合金具有低密度、高比強(qiáng)度，在滿(mǎn)足性能的前提下可比鋼件質(zhì)量減少60%.我國(guó)的鋁資源豐富，且汽車(chē)用鋁的90%都可以循環(huán)利用，因此鋁合金在汽車(chē)上的推廣應(yīng)用是必然趨勢(shì)［1］.同時(shí)，由于鋁合金彈性模量只有鋼板的1/3，所以鋁合金具有比鋼板更加嚴(yán)重的回彈問(wèn)題［2-3］.回彈的準(zhǔn)確預(yù)測(cè)及其有效控制是推廣鋁合金應(yīng)用面臨的關(guān)鍵問(wèn)題.

從塑性力學(xué)角度，一個(gè)完整的材料模型應(yīng)該包括3個(gè)方面:描述材料應(yīng)力應(yīng)變關(guān)系的流動(dòng)方程、判斷材料是否進(jìn)入塑性屈服的屈服準(zhǔn)則、定義后繼屈服面變化軌跡的硬化模型［4-5］.現(xiàn)有的流動(dòng)方程可分為關(guān)聯(lián)和非關(guān)聯(lián)流動(dòng)方程［6］.關(guān)聯(lián)流動(dòng)理論認(rèn)為，塑性勢(shì)函數(shù)與屈服準(zhǔn)則采取相同的函數(shù).目前，在金屬薄板的成形仿真中多采用關(guān)聯(lián)流動(dòng)方程.常用的屈服準(zhǔn)則有包括 Hill’48［7］、Hill’79［8］在內(nèi)的Hill系列以及Barlat’89［9］、Barlat’91［10-11］等系列的屈服準(zhǔn)則.硬化模型主要分為各向同性硬化模型、隨動(dòng)硬化模型和混合硬化模型［12-13］.不同的屈服準(zhǔn)則和硬化模型的組合加上關(guān)聯(lián)流動(dòng)方程就組成了不同的材料模型［14-15］.

5052鋁板屬于Al-Mg系鋁合金，由于其具有優(yōu)異的防銹、高抗疲勞強(qiáng)度、良好的塑性與耐腐蝕性等優(yōu)點(diǎn)在汽車(chē)上應(yīng)用越來(lái)越廣泛.由于鋁合金屬于面心立方晶體結(jié)構(gòu)，且存在嚴(yán)重的織構(gòu)，這些使得5052鋁合金的彈塑性力學(xué)行為區(qū)別于其他系列的鋁合金，選擇適合于5052鋁合金的材料模型對(duì)提高沖壓成形和回彈仿真精度尤為重要［16］.因此，本文對(duì) LS-Dyna軟件中 4個(gè)材料模型MAT_122、MAT_36、MAT_125、MAT_226［17］所采用的屈服準(zhǔn)則和硬化模型進(jìn)行比較分析，并采用這4個(gè)材料模型進(jìn)行U形件的回彈仿真，通過(guò)U形件的沖壓成形實(shí)測(cè)回彈，分析屈服準(zhǔn)則和硬化模型對(duì)回彈精度的影響.

1 材料模型

1.1 MAT_36模型

LS-Dyna軟件中，MAT_36材料模型采用的是Barlat三參數(shù)屈服準(zhǔn)則和Hollomon各向同性硬化模型.Barlat三參數(shù)屈服準(zhǔn)則是由Barlat和Lian［9］在1989年提出來(lái)，故又稱(chēng)為Barlat’89屈服準(zhǔn)則.該準(zhǔn)則適合于平面應(yīng)力條件下面內(nèi)各向異性的材料，其屈服函數(shù)如式(1)所示:

式中:σx、σy及σxy是平面應(yīng)力分量;m 為材料常數(shù)，鋁合金屬于面心立方晶體結(jié)構(gòu)，其m值為8;參數(shù) a、c、h 及 p 可由各向異性指數(shù) r0、r45和 r90計(jì)算得到.

各向同性硬化理論認(rèn)為，后繼屈服面的尺寸在改變，而屈服面位置不變，即塑性應(yīng)力只與塑性應(yīng)變累積有關(guān).Hollomon各向同性硬化模型［18］可表示為

式中:σ和εp為流動(dòng)應(yīng)力和塑性應(yīng)變，K為材料強(qiáng)度系數(shù)，n為應(yīng)變硬化指數(shù).

1.2 MAT_122模型

LS-Dyna中，MAT_122材料模型采用了考慮面內(nèi)各向異性的Hill’48屈服準(zhǔn)則和Hollomon硬化模型，適用于殼單元和各向異性材料的成形仿真.Hill’48屈服準(zhǔn)則如式(4)所示:

式中，F(xiàn)、G、H、L、M 和 N 等參數(shù)可由 3個(gè)方向的屈服強(qiáng)度和剪切強(qiáng)度來(lái)確定［14］.

1.3 MAT_125模型

LS-Dyna軟件中的MAT_125材料模型是由Hill’48屈服準(zhǔn)則與 Yoshida-Uemori(以下簡(jiǎn)稱(chēng)Y-U模型)［19-20］非線性隨動(dòng)硬化模型結(jié)合而來(lái)，適用于殼單元和體單元分析.與各向同性硬化模型相比，隨動(dòng)硬化模型認(rèn)為后繼屈服面的位置隨著塑性應(yīng)變而改變.描述后繼屈服面中心點(diǎn)的軌跡矢量就稱(chēng)為背應(yīng)力.目前常用的隨動(dòng)硬化模型有Chaboche系列模型［21-23］和Y-U模型，其中YU模型由日本學(xué)者Yoshida和Uemori合作提出.Y-U模型屬于雙屈服面模型，該模型采用屈服面和邊界面來(lái)描述金屬材料的塑性屈服和后繼屈服行為.邊界面與屈服面一樣具有隨動(dòng)硬化特性，邊界面的存在使正向與反向應(yīng)變路徑產(chǎn)生非各向同性硬化，由此可以描述硬化遲滯與后續(xù)整體硬化現(xiàn)象.該模型中背應(yīng)力的演變遵從:

式中:dα和dβ分別是屈服面中心和邊界面中心;dθ為屈服面相對(duì)邊界面的移動(dòng)矢量是θ的等效塑性勢(shì)，其物理意義是邊界面尺寸與屈服面尺寸的差;B是邊界面尺寸的初始值;Y是初始屈服應(yīng)力;R是各向同性硬化模量，可由式(9)確定:

1.4 MAT_226模型

LS-Dyna中，MAT_226材料模型采用了Barlat’89屈服準(zhǔn)則與Y-U隨動(dòng)硬化模型.屈服準(zhǔn)則表達(dá)式與式(1)相同，硬化模型由式(5)～(9)確定.

2 回彈仿真和U形件回彈試驗(yàn)

2.1 回彈仿真

U形件模具尺寸如圖1所示.仿真中需要的鋁合金5052的材料參數(shù)由單向拉伸試驗(yàn)獲得.板坯為長(zhǎng)300 mm、寬40 mm的矩形板條，長(zhǎng)度方向與軋制方向相同.考慮模型和邊界條件的對(duì)稱(chēng)性，仿真采用1/2模型.板料網(wǎng)格尺寸為2 mm，在凹模圓角對(duì)應(yīng)的部位對(duì)網(wǎng)格進(jìn)行了細(xì)分.沖壓速度為1 000 mm/s，壓邊力為10 kN，和試驗(yàn)壓邊力一致.為分析前述4種材料模型的回彈仿真精度，分別采用 MAT_36、MAT_122、MAT_125 和 MAT_226模型對(duì)U形件的回彈進(jìn)行仿真.所有仿真均采用16號(hào)全積分單元，板厚方向7個(gè)積分點(diǎn)，接觸類(lèi)型采用雙面接觸模型FORMING_ONE_WAY_S_S，摩擦系數(shù)0.125.顯式算法模擬沖壓成形，隱式算法進(jìn)行回彈仿真.仿真中采用的基本材料參數(shù)密度、彈性模量、泊松比、屈服強(qiáng)度均采用試驗(yàn)獲得的數(shù)據(jù).Y-U硬化模型參數(shù)采用參考文獻(xiàn)中的數(shù)值:Y=80 MPa、B=140 MPa、Rsat=180 MPa、b=22 MPa、C=900［19-20］.

圖1 U形件沖壓模具尺寸

2.2 單向拉伸試驗(yàn)

本文研究對(duì)象5052鋁板由諾貝里斯公司生產(chǎn).為獲得回彈仿真所需的基本力學(xué)性能參數(shù)，對(duì)5052鋁板進(jìn)行了單向拉伸試驗(yàn)，板料厚度1.0 mm.按照ASTM E8標(biāo)準(zhǔn)準(zhǔn)備試樣尺寸，標(biāo)距50 mm.采用線切割沿板料軋制方向成0°、45°和90°方向加工試樣.拉伸速度為3 mm/min.

拉斷后的試樣如圖2所示，其斷口截面方向與試樣的拉伸方向斜交.

圖2 拉斷后的單向拉伸試樣

圖3所示是試驗(yàn)測(cè)得的工程應(yīng)力-應(yīng)變曲線.采用數(shù)據(jù)擬合的方法得到彈性模量為70.6 GPa.由于圖3所示曲線沒(méi)有明顯的屈服平臺(tái)，取0.2%工程應(yīng)變下的工程應(yīng)力為屈服強(qiáng)度，大小為165.6 MPa.抗拉強(qiáng)度為215.7 MPa，平均均勻延伸率約為10%.根據(jù)標(biāo)準(zhǔn)ASTM E646求得應(yīng)變硬化指數(shù)為0.12.

(2)生石灰活性可影響合成硬硅鈣石纖維質(zhì)量，石灰活性越高，生成的硬硅鈣石纖維直徑越小，體積密度也越小。當(dāng)煅燒溫度為1 000 ℃時(shí)，生成的生石灰與平均粒徑23 μm的晶質(zhì)石英粉反應(yīng)，合成的硬硅鈣石纖維直徑約為82 nm，體積密度為70.4 kg/m3。

圖3 5052鋁合金板工程應(yīng)力-應(yīng)變曲線

2.3 沖壓試驗(yàn)和回彈變形的測(cè)量

U形件沖壓試驗(yàn)?zāi)＞叱叽绾兔鞒叽缗c仿真完全相同.沖壓行程為 45 mm，沖壓速度為4 mm/s，沖壓之前用潤(rùn)滑油對(duì)毛坯進(jìn)行雙面潤(rùn)滑，壓邊力為10 kN.

圖4 試驗(yàn)U形件

圖4所示為已經(jīng)發(fā)生回彈的U形件，可以看到U形件的底部、側(cè)壁和法蘭部分都發(fā)生了比較明顯的回彈.U形件底部產(chǎn)生了向外凸起，同時(shí)側(cè)壁向內(nèi)凸起形成了具有一定曲率的圓弧，法蘭部分下垂.形狀的回復(fù)使得側(cè)壁和底部的夾角由直角變?yōu)殁g角，法蘭與側(cè)壁的夾角變?yōu)殇J角.

沖壓結(jié)束后將U形件輪廓印制在白紙上，掃描獲得U形件輪廓，將其導(dǎo)入到CAD軟件中提取輪廓.然后參考 Numisheet’2011會(huì)議確定的回彈測(cè)量方法［24］在CAD軟件中對(duì)U形件輪廓的回彈角和法蘭的下垂量進(jìn)行測(cè)量［19］.如圖5所示，以U形件底部中心點(diǎn)的水平切線為基準(zhǔn)線，將基準(zhǔn)線向上平移10 mm得到其和U形件輪廓的交點(diǎn)A，再以A點(diǎn)為圓心做半徑為20 mm的圓弧和U形件輪廓交于B點(diǎn)，連接A、B兩點(diǎn)并沿正反方向延長(zhǎng)，AB延長(zhǎng)線與底部基準(zhǔn)線的夾角定義為θ1，AB延長(zhǎng)線和法蘭延長(zhǎng)線CD的夾角定義為θ2，D點(diǎn)在回彈前后的高度方向的改變量為下垂量Δz.

為盡可能減少試驗(yàn)偶然誤差的影響，試驗(yàn)共測(cè)量了4件，取平均值作為最終回彈試驗(yàn)結(jié)果.

圖5 回彈測(cè)量方法示意圖

3 結(jié)果分析與討論

3.1 回彈前后Von Mises應(yīng)力

U形件的沖壓成形屬于二維拉彎變形，主要變形區(qū)域是凹模圓角區(qū)、凸模圓角區(qū)和側(cè)壁區(qū).材料流經(jīng)圓角區(qū)要經(jīng)歷彎曲-反彎曲變形，側(cè)壁區(qū)主要是單向拉伸變形，底部與法蘭區(qū)變形很小.這類(lèi)變形極易使材料經(jīng)受循環(huán)硬化，循環(huán)硬化下材料極易表現(xiàn)包辛格效應(yīng)［25］.另外，沖壓件的回彈就是彈性卸載引起的形狀回復(fù)，因此，U形件圓角區(qū)和側(cè)壁區(qū)的等效應(yīng)力在回彈前后的變化情況體現(xiàn)了相應(yīng)材料模型在描述材料彈塑性行為準(zhǔn)確性的一個(gè)方面.下面以MAT_125材料模型為例說(shuō)明回彈前后的Von Mises等效應(yīng)力變化情況，如圖6所示.其他材料模型計(jì)算的等效應(yīng)力變化列于表1中.

圖6 MAT_125模型仿真Von Mises應(yīng)力

表1 回彈前后的等效應(yīng)力

由圖6(a)可知，回彈前最大等效應(yīng)力出現(xiàn)在凹模圓角區(qū)和凸模圓角區(qū)上方的區(qū)域，側(cè)壁區(qū)域應(yīng)力也比較大，這3個(gè)區(qū)域是U形件的主要變形區(qū).在圖6(b)所示的回彈后等效應(yīng)力有顯著下降，最大等效應(yīng)力由回彈前的296.1 MPa降為回彈后的176.6 MPa.尤其是上下2個(gè)圓角部位，回彈前最大應(yīng)力在這2個(gè)部位，而回彈后最大應(yīng)力出現(xiàn)圖6(b)中箭頭所示的側(cè)壁區(qū).這說(shuō)明回彈過(guò)程中圓角區(qū)的應(yīng)力得到了充分釋放.

如前所述，圓角區(qū)經(jīng)歷了彎曲-反彎曲的循環(huán)變形，也就是材料經(jīng)歷了拉伸到壓縮再到拉伸的變形循環(huán)，多數(shù)材料在循環(huán)變形下會(huì)表現(xiàn)出隨動(dòng)硬化現(xiàn)象，如包辛格效應(yīng).本文所用4個(gè)材料模型中，MAT_125和MAT_226采用的是Y-U隨動(dòng)硬化模型，而其他2個(gè)模型采用的是各向同性硬化模型.各向同性硬化理論認(rèn)為，材料在拉壓狀態(tài)下具有相同的塑性行為，無(wú)法描述包辛格效應(yīng)，因此仿真所用材料模型是否具備描述這個(gè)現(xiàn)象的能力極大影響了回彈仿真精度.

3.2 回彈仿真結(jié)果分析

圖7所示為采用4個(gè)材料模型仿真所得U形件輪廓.與卸載前輪廓相比，MAT_226模型預(yù)測(cè)的回彈最大，其次是MAT_125和MAT_36，MAT_122最小.

圖7 U形件回彈仿真所得零件輪廓

按照2.3小節(jié)中的回彈測(cè)量方法對(duì)角θ1和θ2以及法蘭的下垂量Δz進(jìn)行了測(cè)量，測(cè)量結(jié)果如表2 所示.此處角度 θ1和 θ2與圖 5 中含義相同，θ1和θ2與90°的差值就為相應(yīng)的回彈量，Δz是回彈前與回彈后法蘭端部在沖壓方向z向的位移量.因此，較大的θ1和較小的θ2以及較大的Δz意味著較大的回彈.由表2可得:試驗(yàn)測(cè)得的θ1與θ2的回彈量分別為 12.9°和 7.2°，采用 MAT_36、MAT_122、MAT_125和MAT_226模型仿真所得θ1的回彈量分別為 7.2 °、3.9 °、10.7 °和 10.8 °，相應(yīng)的 θ2的回彈量分別為 5.1 °、2.3 °、7.6 °和 7.4 °.預(yù)測(cè)結(jié)果與試驗(yàn)值相比有如下特點(diǎn):

1)4個(gè)材料模型仿真的回彈結(jié)果均比試驗(yàn)值小，其中MAT_226模型預(yù)測(cè)的回彈角度與試驗(yàn)值最接近，MAT_125模型次之，MAT_122模型預(yù)測(cè)結(jié)果與試驗(yàn)值相差較大.

2)各向同性硬化模型比隨動(dòng)硬化模型預(yù)測(cè)的回彈小.MAT_122和MAT_36模型均建立在各向同性硬化的基礎(chǔ)上，這2個(gè)模型預(yù)測(cè)的回彈均比基于Y-U隨動(dòng)硬化模型的 MAT_125模型和MAT_226模型預(yù)測(cè)結(jié)果小.

3)硬化模型相同時(shí)，Barlat’89屈服準(zhǔn)則比Hill’48屈服準(zhǔn)則具有更高的回彈預(yù)測(cè)精度.如MAT_36模型預(yù)測(cè)的回彈角及法蘭下垂量均比MAT_122模型預(yù)測(cè)結(jié)果更接近試驗(yàn)值，MAT_226模型預(yù)測(cè)結(jié)果比MAT_125模型預(yù)測(cè)結(jié)果更接近試驗(yàn)結(jié)果.這說(shuō)明對(duì)本文5052鋁合金板而言，Barlat’89具有比Hill’48更高的精度.

4)硬化模型對(duì)回彈預(yù)測(cè)精度的影響程度要大于屈服準(zhǔn)則.MAT_226與MAT_125模型均是基于Y-U隨動(dòng)硬化模型，而MAT_36與MAT_122模型均采用了Hollomon各向同性硬化模型.表2中MAT_226模型與MAT_125模型的回彈預(yù)測(cè)精度顯著高于MAT_36和MAT_122模型.

表2 試驗(yàn)與仿真回彈測(cè)量結(jié)果

4 結(jié)論

1)對(duì)本文5052鋁合金而言，基于Yoshida-Uemori隨動(dòng)硬化和Barlat’89屈服準(zhǔn)則的材料模型預(yù)測(cè)的回彈結(jié)果與試驗(yàn)值最為接近;而基于各向同性硬化和Hill’48屈服準(zhǔn)則的材料模型回彈預(yù)測(cè)精度最低.

2)隨動(dòng)硬化模型比各向同性硬化模型具有更高的回彈預(yù)測(cè)精度;

3)硬化模型對(duì)回彈預(yù)測(cè)精度的影響比屈服準(zhǔn)則的影響大.

［1］KIM H S，KOC M.Numericalinvestigations on springback characteristics of aluminum sheet metal alloys in warm forming conditions［J］.Journal of Materials Processing Technology，2008，204:370-383.

［2］陳超，桂楓，陳明安.鋁合金板材彎曲成形性能［J］.鍛壓技術(shù)，2013，38(1):25-30.CHEN Chao，GUI Feng，CHEN Mingan.Bending formability of aluminum alloy sheet［J］.Forming ＆Stamping Technology，2013，38(1):25-30.

［3］王強(qiáng).鋁合金車(chē)身覆蓋件沖壓成形回彈仿真方法研究［J］.農(nóng)業(yè)裝備與車(chē)輛工程，2012，50(4):50-54.WANG Qiang.Research of springback simulation method in aluminum alloy auto body panel stamping［J］.Agricultural Equipment＆ Vehicle Engineering，2012，50(4):50-54.

［4］肖煜中，陳軍.金屬板料沖壓數(shù)值模擬中的宏觀硬化模型研究現(xiàn)狀［J］.塑性工程學(xué)報(bào)，2009，16(4):51-58.XIAO Yuzhong，CHEN Jun.A review of research in macroscopic hardening models in numerical simulation of sheet metal forming［J］.Journal of Plasticity Engineering，2009，16(4):51-58.

［5］LEE M G，KIM C，PAVLINA E J，et al.Advances in sheet forming-materials modeling， numerical simulation and press technologies［J］.J Manufacturing Sci Eng，2011，133(6):10011-10012.

［6］SAFAEI M，ZANG S L，LEE M G，et al.Evaluation of anisotropic constitutive models:mixed anisotropic hardening and non-associated flow rule approach［J］.Int J Mech Sci，2013，73:53-68.

［7］HILL R.A theory of the yielding and plastic flow of anisotropic metals［J］.Proceedings of the Royal Society of London.SeriesA，Mathematicaland Physical Sciences，1948，193:281-297

［8］HILL R.Constitutive modeling of orthotropic plasticity in sheet metals［J］.J Mech Phys Solids，1989，38:405-417.

［9］BARLAT F，LIAN J.Plastic behavior and stretchability of sheet metals.Part I:a yield function for orthotropic sheet under plane stress condition［J］.Internation-al Journal of Plasticity，1989，5:51-66.

［10］BARLAT F，LEGE D J，BREM J C.A six-component yield function for anisotropic materials［J］.International Journal of Plasticity，1991，7:693-712.

［11］BARLAT F，BREM J C，YOON J W，et al.Plane stress yield function for aluminum alloy sheets-Part 1:theory［J］.International Journal of Plasticity，2003，19:1-23.

［12］CAO J，LEE W，CHENG H S，et al.Experimental and numerical investigation of combined isotropickinematic hardening behavior of sheet metals［J］.International Journal of Plasticity，2009，25:942-972.

［13］NOBUTADA O，WANG J D.On modeling of kinematic hardeningforratcheting behavior ［J］.Nuclear Engineering and Design，1995，153:205-212.

［14］EGGERTSEN P A，MATTIASSON K.On constitutive modeling for springback analysis［J］.Int J Mech Sci，2010，52:804-818.

［15］WAGONER R H，LIM H，LEE M G.Advanced issues in springback ［J］.Int J Plast，2013，45:3-20.

［16］ZHU Y X，LIU Y L，YANG H，et al.Development and application of the material constitutive model in springback prediction of cold-bending［J］.Materials＆ Design，2012，42:245-258.

［17］Hallquist J O.LS-DYNA keyword user's manual［M］.California:Livermore Software Technology Corporation，2007，970.

［18］GRONOSTAJSKI Z.Theconstitutiveequationsfor FEM analysis［J］.Journal of Materials Processing Technology，2000，106:40-44.

［19］YOSHIDA F，UEMORI T.A model of large-strain cyclic plasticity describing the Bauschinger effect and workhardening stagnation［J］.International Journal of Plasticity，2002，18:661-686.

［20］YOSHIDA F，UEMORI T.A model of large-strain cyclic plasticity and its application to springback simulation ［J］.Int J Mech Sci，2003，45:1687-1702.

［21］CHABOCHE J L.Thermodynamic formulation of constitutive equations and application to the viscoplasticity and viscoelasticity of metals and polymers［J］.Int J Solids and Structures，1997，34(18)，2239-2254.

［22］CHABOCHE J L.A review of some plasticity and viscoplasticity constitutive theories［J］.International Journal of Plasticity，2008，24(10):1642-1693.

［23］CHABOCHE J L，GAUBERT A，KANOUT?P，et al.Viscoplastic constitutive equations of combustion chamber materials including cyclic hardening and dynamic strain aging ［J］.International Journal of Plasticity，2013，46:1-22.

［24］CHUNG K，KUWABARA T，VERMA R，etal.Numisheet 2011 Benchmark 4:Pre-strain effect on spring-back of 2D draw bending［C］//The 8th International Conference and Workshop on Numerical Simulation of 3D sheet metal forming processes.New York:AIP Publishing，2011:171-175.

［25］莊京彪，劉迪輝，李光耀.基于包辛格效應(yīng)的回彈仿真分析［J］.機(jī)械工程學(xué)報(bào)，2013，49(22):84-90.ZHUANG Jingbiao，LIU Dihui，LI Guangyao.Analysis of springback simulation based on bauschinger effect［J］.Jounal of Mechanical Engineering，2013，49(22):84-90.