基于FEM-DEM的齒輪鋼隨機(jī)噴丸模型及殘余應(yīng)力仿真研究

嚴(yán)宏志，伊偉彬，朱鵬飛，吳順興

（中南大學(xué) 機(jī)電工程學(xué)院，長沙 410000）

0 引言

噴丸強(qiáng)化工藝廣泛應(yīng)用于航空航天、汽車制造以及工程機(jī)械等行業(yè)中各種尺寸和樣式的工件，尤其對承受交變載荷齒輪、螺旋槳葉、曲軸以及軸承等金屬制件能大幅提高疲勞強(qiáng)度。

噴丸強(qiáng)化數(shù)值仿真是研究噴丸工藝與機(jī)理的重要途徑，相關(guān)研究已從最初的二維模型發(fā)展到三維模型，從單個彈丸噴射模型發(fā)展到考慮了覆蓋率的多個規(guī)則位置彈丸模型。Johnson[1]等最早建立單彈丸模型，分析了噴丸強(qiáng)化工藝的殘余壓應(yīng)力。Schiffner[2]、Guagliano[3]等建立了規(guī)則的多彈丸噴丸模型，并與單彈丸模型做了對比。近年來隨機(jī)彈丸模型的研究得到學(xué)者們關(guān)注：盛湘飛[4]等對鋁合金進(jìn)行了隨機(jī)噴丸模型研究，通過再次導(dǎo)入殘余應(yīng)力的方式將長時間的仿真分段來提高效率。YangF等[5]結(jié)合規(guī)則分布彈丸與隨機(jī)彈丸，建立了周期性模型。S.M.H.Gangaraj[6]對比了規(guī)則彈丸模型和隨機(jī)彈丸模型，發(fā)現(xiàn)在引入殘余應(yīng)力方面結(jié)果較接近，但規(guī)則彈丸在粗糙度和覆蓋率的預(yù)測上相對弱勢。在受噴對象方面，李智等[7]建立了面對曲面的噴丸模型，對比了不同噴射角度對殘余應(yīng)力的影響。王延忠[8]研究了面向TC4鈦合金板件的噴丸仿真模型，通過提取彈丸坐標(biāo)的方法預(yù)測噴丸覆蓋率。不難看出，隨機(jī)彈丸模型和不規(guī)則表面受噴是噴丸強(qiáng)化工藝仿真的發(fā)展方向，但由于要達(dá)到一定覆蓋率時，既需要大量彈丸，同時受噴表面的網(wǎng)格也需要劃分的十分精細(xì)，因此前人的噴丸隨機(jī)模型不可避免的會帶來巨大的計算量和分析時間。

隨機(jī)彈丸模型可以更準(zhǔn)確反映噴丸機(jī)理，但達(dá)到一定覆蓋率時，既需要大量彈丸，同時受噴表面網(wǎng)格需要更精細(xì)，給噴丸隨機(jī)模型的計算帶來巨大了計算量和分析時間。

應(yīng)用離散元模型可以準(zhǔn)確的表征大量隨機(jī)顆粒，模擬大量彈丸對受噴工件打擊的集群效應(yīng)和彈丸粒子之間的相互影響，還可大大減少了計算時間。而離散元模型難以表征連續(xù)體的應(yīng)力應(yīng)變，有限元模型則能準(zhǔn)確反映受噴工件的彈塑性變形以及引入的殘余應(yīng)力，因此建立FEM-DEM耦合模型可較好的解決隨機(jī)噴丸強(qiáng)化仿真問題。但以往的FEM-DEM模型需要應(yīng)用EDEM軟件與有限元軟件的聯(lián)合仿真，前后處理繁瑣。

本文基于ABAQUS軟件加入P3D3離散單元以及離散元模型，應(yīng)用python語言通過對ABAQUS的前處理inp文件重新編寫，建立了FEM-DEM耦合模型來表征噴丸強(qiáng)化的過程，較好的滿足了隨機(jī)噴丸數(shù)值仿真的高效率需求。

1 仿真模型建立

噴丸強(qiáng)化實(shí)質(zhì)上是大量多數(shù)且高速撞擊在受噴表面，在短時間內(nèi)快速引起彈性和塑性變形。本文所建立的FEM-DEM模型可大致分為兩部分，有限元分析部分是基于ABAQUS/Explicit內(nèi)核建立，主要是受噴件齒輪模型幾何模型和材料本構(gòu)模型；另一部分是根據(jù)PARTICLE GENERATOR顆粒流模塊建立的噴射彈丸的離散元模型。

1.1 齒輪幾何模型

以往的受噴模型多簡化為平面板件，本文則采用了株洲齒輪廠應(yīng)用于某型號車橋中的齒輪模型，該齒輪使用20CrMnTi齒輪鋼，根據(jù)表中參數(shù)使用SolidWorks建立三維模型。為了后續(xù)仿真簡化模型，如圖1所示，截取該齒輪模型中的一個齒根及左右兩齒，并導(dǎo)入abaqus。

圖1 截取的齒輪模型

1.2 材料模型

噴丸強(qiáng)化高度非線性且發(fā)生在瞬時，具有大變形和高應(yīng)變速率。Johnson-Cook本構(gòu)模型能夠準(zhǔn)確反映材料在塑性變形過程中應(yīng)變硬化、應(yīng)變速率硬化效應(yīng)和溫度軟化效應(yīng)[9]。

ε為等效應(yīng)變；

T為溫度，℃；

Tm為材料熔化溫度，Tm=1350℃；

T0為溫度參考值，通常取室溫T0=25℃。

王佳斌對20CrMnTi鋼進(jìn)行動態(tài)壓縮實(shí)驗(yàn)，獲得了不同應(yīng)變率下的動態(tài)應(yīng)力-應(yīng)變曲線，獲得了修正后的J-C本構(gòu)模型[10]，故本文參考該模型來描述材料形變特性。其主要參數(shù)及相應(yīng)數(shù)值如表1所示。其中A為屈服應(yīng)力，單位MPa；B為冪指系數(shù)，單位MPa；n為硬化指數(shù)；C為應(yīng)變率敏感性系數(shù)；m為溫度敏感性系數(shù)。

表1 20CrMnTi合金鋼JC本構(gòu)模型系數(shù)

1.3 網(wǎng)格劃分

網(wǎng)格劃分的越精細(xì)結(jié)果就越準(zhǔn)確，但同時也會導(dǎo)致計算量成倍增加。高應(yīng)變速率的模型對網(wǎng)格尺寸精密程度更加敏感[11]，M.Frija建立了不同的最小單元尺寸模型，用仿真的殘余應(yīng)力結(jié)果對比理論計算結(jié)果，結(jié)論為最小單元的尺寸相比彈丸直徑應(yīng)至少小于1/10[12]。為了優(yōu)化分析時間，采用單個積分點(diǎn)的C3D8R減縮積分單元。同時，在受噴區(qū)域表面做網(wǎng)格細(xì)化，使受噴區(qū)域的網(wǎng)格尺寸小于0.08mm。并在垂直表面方向做縱向梯度細(xì)化，以方便后處理時提取數(shù)據(jù)。劃分網(wǎng)格后的模型如圖2所示，共171360個單元。

圖2 受噴齒的網(wǎng)格劃分

1.4 邊界條件與載荷

為防止在仿真過程中齒輪模型被彈丸撞擊產(chǎn)生位移，對截取出來的齒輪模型底面做完全固定。

1.5 離散顆粒流模型

本文所建立的模型面向高強(qiáng)度齒輪，選用高強(qiáng)度CCW-φ0.8mm/G3高硬度鋼絲切丸。

通過對inp文件的編寫，定義彈丸與受噴表面之間的接觸為hard接觸，動摩擦因數(shù)μ為0.2。對于彈丸這種球形模型，根據(jù)Hertz接觸理論來描述的法向接觸力與法向位移的關(guān)系[13]，來定義彈丸之間的接觸，最終離散元模型如圖所示。其基本參數(shù)如表2所示。

表2 鋼絲切丸參數(shù)

圖3 隨機(jī)彈丸離散模型

2 仿真模擬優(yōu)化與結(jié)果分析

2.1 沙漏控制

由于減縮積分單元只有一個積分點(diǎn)，在某些時候單元有形變，而插值計算的應(yīng)變?yōu)?進(jìn)而導(dǎo)致內(nèi)能為零，這種情況稱為沙漏模式。避免使用減縮積分單元出現(xiàn)沙漏模式的根本方法是細(xì)化網(wǎng)格或是屬于用多積分點(diǎn)單元，但這會增加計算時間，降低計算效率。人工引入“沙漏剛度”，可以一定程度上限制沙漏變形。仿真軟件對沙漏剛度取值范圍為0.2～3，當(dāng)沙漏能與總內(nèi)能的比值低于10%時，可視為該模型的沙漏模式得到比較好的控制，即該模型可用[14]。本文分別取沙漏剛度為0.2，0.5，1，2的沙漏控制進(jìn)行分析。

圖4所示為同一個模型下不同的沙漏剛度值產(chǎn)生的應(yīng)變能與內(nèi)能隨時間的變化。當(dāng)沙漏剛度為0.2，0.5，1時，沙漏能ALLAE隨著仿真時間的增加增大，最大值（15.48mJ，18.82mJ，22.38mJ）與總內(nèi)能ALLIE（136.59mJ）的比值大于了10%，這表明引入的沙漏剛度不足，產(chǎn)生了明顯的沙漏現(xiàn)象。當(dāng)沙漏剛度取到2時，沙漏能增幅度很大，但很快趨于穩(wěn)定不變，這表明引入大的沙漏剛度有效控制了沙漏現(xiàn)象。

圖4 不同沙漏剛度參數(shù)引起的沙漏能

2.2 表面覆蓋率計算方法

表面覆蓋率是噴丸工藝的一項(xiàng)基本工藝參數(shù)，它是指彈丸沖擊工件表面后留下的凹坑面積占受噴面積的比率。在以往的仿真模型中，如圖5所示，多采用規(guī)則的多彈丸模型，按照指定的位置沖擊受噴工件，因而可以預(yù)設(shè)覆蓋率。本文基于隨機(jī)模型研究了一種新型表征覆蓋率與彈丸數(shù)量的方法，以單彈丸模型出發(fā)，再進(jìn)一步擴(kuò)展到隨機(jī)彈丸模型。

圖5 傳統(tǒng)噴丸仿真覆蓋率表征方法

2.2.1 單彈丸模型的受噴區(qū)域界定研究

建立簡化的單彈丸有限元模型，仿真后的塑性變形如圖6所示。提取受噴工件沿深度方向的位移，以離彈坑的距離為橫軸繪制圖7（縱軸參考右軸），則距離彈坑中心第一個z向位移最大時的距離即是彈坑的半徑，可得彈坑直徑為3.38×10-4。

圖6 單彈丸噴射有限元模型

圖7 彈坑縱向變形和等效塑性變形曲線

圖7中的紅色曲線反映了彈丸沖擊下沿深度方向的變形情況，曲線形狀約為彈坑截面，彈坑面的Z向塑性應(yīng)變值范圍為-0.024mm～0.010mm，因此難以用該值界定臨界值來區(qū)分噴與未噴。噴丸工藝中產(chǎn)生的彈丸位置是隨機(jī)的，當(dāng)要保證一定的覆蓋率時，那必然某些微觀表面會被多個彈丸撞擊，產(chǎn)生的塑性變形也會累計受到彈丸沖擊而產(chǎn)生累加的塑性變形，因此引入等效塑性變形（PEEQ）這一無量綱變量，如圖所示繪制了以離彈坑的距離為變量，兩段截面上的等效塑性變形曲線。不難得出，隨著距離彈坑中心原來越近，PEEQ的值逐漸增大，到彈坑中心達(dá)到最大值，PEEQ的值與離心距成正比。上文得出彈坑直徑3.38×10-4，由圖7可知，在彈坑邊緣，Z向位移最大時，PEEQ的值為0.115，因此彈坑內(nèi)的PEEQ≥0.115，彈坑外PEEQ值＜0.115。由此可得，以等效塑性應(yīng)變的值作為界定受噴區(qū)域時，PEEQ≥0.115的單元可視為受噴區(qū)域，反之則視為未噴區(qū)域。

2.2.2 多彈丸模型覆蓋率與彈丸個數(shù)之間的關(guān)系

當(dāng)彈丸直徑、速度與材料等條件一定時，噴丸時間越大，噴射到工件表面的彈丸越多，覆蓋率越高。如圖8所示，本文取部分齒根和齒底直徑6mm的圓作為總受噴區(qū)域，將該部分離散為1015個單元，再提取受噴的齒根區(qū)域1015個單元的PEEQ值，根據(jù)覆蓋率的定義：凹坑面積比總受噴面積，則覆蓋率Cr的值為：

圖8 離散受噴區(qū)域

式中，Np為PEEQ≥0.115的單元個數(shù)；

Nt表示受噴區(qū)域內(nèi)總的單元數(shù)，1015。

圖9～圖12分別是加載100彈丸、220彈丸、340彈丸和460彈丸的受噴面離散點(diǎn)的等效塑性應(yīng)變值，圖中紅線是判別是否為受噴區(qū)域的分界線，則根據(jù)單彈丸模型的結(jié)論，統(tǒng)計表面1015個單元中PEEQ≥0.115的單元數(shù)目，計算得到不同彈丸沖擊下的覆蓋率分別為：100彈丸沖擊下覆蓋率為34.3%，220彈丸沖擊下覆蓋率為72.7%，340彈丸沖擊下覆蓋率為93.6%，460彈丸沖擊下覆蓋率為98.4%。

圖9 100彈丸沖擊下表面離散點(diǎn)PEEQ值

圖10 220彈丸沖擊下表面離散點(diǎn)PEEQ值

圖11 340彈丸沖擊下表面離散點(diǎn)PEEQ值

圖12 460彈丸沖擊下表面離散點(diǎn)PEEQ值

圖13是加載彈丸數(shù)目從100增加到520所對應(yīng)的覆蓋率變化曲線。不難看出，隨著施加的彈丸數(shù)從100增加至520，覆蓋率逐漸從34.3%增加到99.4%，但整體增加幅度不斷減小。當(dāng)彈丸數(shù)達(dá)到460個時，覆蓋率達(dá)到98.4%，由于邊際遞減效應(yīng)，一般情況下認(rèn)為覆蓋率達(dá)到98%時即達(dá)到100%。因此在該模型中，彈丸數(shù)達(dá)到460個時覆蓋率達(dá)到100%。

圖13 不同彈丸數(shù)目對應(yīng)的覆蓋率

2.3 覆蓋率對殘余應(yīng)力的影響

結(jié)合上文結(jié)論，取彈丸直徑為0.8mm，彈丸速度為50m/s，入射角為90°，分別模擬了50%、100%、200%、300%、400%噴丸覆蓋率下的齒根的受噴情況，得到的應(yīng)力云圖如圖14～圖17所示。而根據(jù)模擬結(jié)果，提取出沿層深方向的殘余應(yīng)力曲線，結(jié)果如圖18所示。

圖14 覆蓋率100%表面殘余應(yīng)力云圖

圖15 覆蓋率200%表面殘余應(yīng)力云圖

圖16 覆蓋率300%表面殘余應(yīng)力云圖

圖17 覆蓋率400%表面殘余應(yīng)力云圖

圖18 不同覆蓋率殘余應(yīng)力曲線

由圖14～圖17可知，100%覆蓋率時，彈坑還較明顯，部分區(qū)域應(yīng)力較大，隨著覆蓋率增加，受噴的齒底與齒根區(qū)域表面殘余應(yīng)力趨于均勻。由圖18可以看到，當(dāng)覆蓋率從50%增加到100%時，表面殘余應(yīng)力由-448.22MPa增加到-845.89MPa，最大殘余應(yīng)力由-753.10MPa增加到-967.04MPa，同時最大殘余應(yīng)力層深和引入的殘余應(yīng)力層深也有所增加；當(dāng)覆蓋率增加到200%時，表面殘余應(yīng)力略微減少了45.34MPa，而最大殘余應(yīng)力略微增加了59.36MPa，層深變化不大；覆蓋率繼續(xù)增加，殘余應(yīng)力曲線介于100%和200%之間，這表明覆蓋率超過一定值后，噴丸強(qiáng)化產(chǎn)生硬化現(xiàn)象，隨著覆蓋率增加，難以產(chǎn)生更多的塑性變形，殘余應(yīng)力曲線趨于穩(wěn)定。而噴丸表面完整性的研究同樣表明[15]，高覆蓋率并不強(qiáng)化殘余應(yīng)力，而是從組織強(qiáng)化的角度來進(jìn)一步提升疲勞強(qiáng)度。

4 結(jié)語

本文應(yīng)用ABAQUS軟件建立的FEM-DEM模型，通過對ABAQUS軟件的二次開發(fā)，高效的實(shí)現(xiàn)了面向齒根齒底噴丸強(qiáng)化的模擬。解決了一般的FEM-DEM模型需要應(yīng)用EDEM軟件和有限元軟件聯(lián)合仿真、前處理和后處理繁瑣的問題。

本文研究了在使用減縮積分單元時，不同的沙漏剛度對沙漏模式產(chǎn)生的影響，當(dāng)取沙漏剛度為2時可以有效控制沙漏模式持續(xù)產(chǎn)生。

基于單彈丸模型研究了受噴區(qū)域的判定問題，當(dāng)PEEQ≥0.115時應(yīng)視為受噴區(qū)域；在此基礎(chǔ)上基于隨機(jī)模型研究了當(dāng)噴嘴直徑為4mm時，覆蓋率與彈丸個數(shù)之間的關(guān)系，結(jié)果表明當(dāng)彈丸數(shù)達(dá)到460時覆蓋率可達(dá)到100%。

研究了覆蓋率對殘余應(yīng)力場的影響，當(dāng)覆蓋率低于100%時，隨著覆蓋率增加殘余應(yīng)力會增大；覆蓋率超過200%后持續(xù)的噴丸強(qiáng)化對殘余應(yīng)力場影響甚微。