噴丸殘余應力場SPH/FEM 耦合建模與仿真

王利平，王建明，裴信超，趙莉莉

(山東大學 機械工程學院，山東 濟南250061)

0 引言

噴丸強化是一種能提高工件表面疲勞壽命的工藝．它通過大量彈丸介質(zhì)以20 ～100 m/s 的速度反復撞擊工件表面，在工件表面形成均勻殘余壓應力層．該殘余壓應力層可抵消工件表面的拉應力，從而有效阻止裂紋的萌生和擴展，顯著提高工件疲勞壽命［1-2］，被廣泛應用于航空、軍工、汽車等領域．因?qū)嶒灧椒ê臅r耗力，且所用的阿爾門試片不能精確反應殘余壓應力的大小和分布． 隨著計算機技術和有限元技術的發(fā)展，已出現(xiàn)許多基于有限元法研究噴丸殘余應力的數(shù)值模型．Mori K 等［3］(1992)通過建立2D 單丸粒對稱模型，對塑性變形工件進行了仿真． S． A． Meguid等［4-5］(1999)通過建立3D 單丸粒動力學模型，研究了彈丸和靶材尺寸對塑性應變及殘余應力的影響．S． A． Meguid 等［5］和G． R． Majzoobi 等［6］利用3D 多丸粒對稱模型對噴丸強化過程進行了系統(tǒng)的分析．凌祥等［7］(2006)建立了模擬噴丸殘余應力場的三維有限元模型．路會龍等［8］(2009)分別建立單、多丸撞擊模型，指出彈丸撞擊速度存在最優(yōu)值． 張洪偉等［9］(2010)利用LS-DYNA 有限元軟件分別建立單、多丸對稱胞元3D 模型，利用單丸模型研究彈丸速度、直徑、入射角、摩擦系數(shù)和重復撞擊次數(shù)對殘余應力及塑性應變的影響;通過分別取4，6，8，9，13 個彈丸建立對稱胞元模型，分析不同覆蓋率對殘余應力場的影響．上述有限元模型僅涉及單丸或數(shù)量有限的多丸，屬微觀模型，與采用大量彈丸反復撞擊工件表面的實際噴丸過程存在明顯差異．

筆者采用SPH 耦合FEM 方法針對高強度鋁合金A12024-T3 通過顯式動力分析軟件LS-DYNA，進行了大量彈丸反復撞擊工件表面的噴丸強化過程數(shù)值模擬;基于兩類不同材料屬性的SPH粒子模擬彈丸流，研究了彈丸速度、噴丸時間等對噴丸殘余應力的影響以及工件表面各典型位置處殘余應力的分布規(guī)律．

1 數(shù)值模型的建立

1．1 兩類SPH 粒子均勻混合模型

筆者采用兩種SPH 粒子分別表示空氣和離散彈丸．其中，空氣對應的SPH 粒子采用常規(guī)氣體狀態(tài)方程，彈丸對應的SPH 粒子采用離散固體的狀態(tài)方程．兩類粒子間的相互作用通過形函數(shù)的特殊定義實現(xiàn)，即在表示粒子形函數(shù)的相鄰粒子支持域中可包含兩類不同性質(zhì)的SPH 粒子．

在SPH 法中，問題域的離散是由有限個具有質(zhì)量并且占有一定空間的粒子組成． 粒子的質(zhì)量可表示為

mi=ΔViρi． (1)

式中:ΔVi為粒子i 占有的體積;ρi為粒子i (i =1，2，…，N)的密度，N 為粒子支持域內(nèi)的粒子數(shù)量．

根據(jù)彈丸尺寸及彈丸流中彈丸與空氣的體積比，分別確定SPH 粒子大小和空氣粒子與彈丸粒子的數(shù)量．為使兩類SPH 粒子均勻混合，使用隨機函數(shù)產(chǎn)生隨機數(shù)，通過MATLAB 編程實現(xiàn)不同屬性的兩類SPH 粒子在彈丸流空間范圍內(nèi)隨機均勻分布．將得到的兩類SPH 粒子隨機分布單元信息寫入LS-DYNA 支持的K 文件中． 圖1 所示為兩類SPH 粒子均勻分布模型，其中左側(cè)粒子代表彈丸，右側(cè)粒子代表空氣．

圖1 兩類SPH 粒子混合模型Fig．1 Hybrid model for two types of SPH particles

1．2 氣體狀態(tài)方程

在LS-DYNA 中，采用MAT_Null 材料建立空氣SPH 粒子，使用狀態(tài)方程EOS_LINEAR_ POLYNOMIAL 表示壓力與密度變化間的關系． EOS_LINEAR_POLYNOMIAL 為線性多項式形式的狀態(tài)方程，可表示為

Ρ=C0+C1μ +C2μ2+C3μ3+(C4μ4+C5μ5+C6μ6)Ε． (2)

式中:μ =ρ/ρ0-1，表示當前密度與參考密度間的比值，ρ0為初始參考密度，在材料模式MAT_Null中定義;C0～C6為各項系數(shù)，對于氣體狀態(tài)方程取

式中:Cp為定壓比熱;Cv為定容比熱．

上述狀態(tài)方程相關參數(shù)使用AUTODYNAVersion6．1 材料庫中的參數(shù)，如表1 所示．

表1 氣體材料和狀態(tài)方程相關參數(shù)Tab．1 Relative parameters of the air material and state equation

1．3 丸粒狀態(tài)方程

丸粒SPH 粒子使用MAT_Null 模擬，考慮彈丸流具有離散介質(zhì)的特性，根據(jù)文獻［10］推導彈丸SPH 粒子的狀態(tài)方程，其為特殊形式的Gruneisen 方程:

式中:Cs為聲音在彈丸材料中的傳播速度，對于鋼丸其值為5．2 km/s;ρ0為彈丸材料初始密度．

故彈丸SPH 粒子也可以使用線性多項式形式的狀態(tài)方程EOS_LINEAR_POLYNOMIAL，此時有

彈丸材料和狀態(tài)方程參數(shù)如下:密度ρ0=7．8 g·cm-3;多項式系數(shù)C1=2．096．

1．4 工件材料本構方程

工件材料為A12024-T3，力學性能參數(shù)如下:彈性模量E=71．3 GPa;密度ρ0=2．78 g·cm-3;泊松比γ=0．3．

由于撞擊過程速度較高且高應變率，需要考慮溫度對應力流動的影響． 工件材料采用Johnson-Cook 黏塑性本構方程模擬，屈服應力σy與等效塑性應變的關系式為

式中:A、B、C、n、m 是通過實驗獲得與應變、應變率、溫度相關的常數(shù)． Dabboussi 等人提出的相關參數(shù)使得仿真結(jié)果與實驗結(jié)果具有較高的一致性［11］，見表為等效塑性應變=等于1 s-1．

式中:T 為任一時刻的瞬時溫度;Tm為熔解溫度，Tr為參考溫度取292K．

表2 Johnson-Cook 相關參數(shù)Tab.2 Relative parameters of Johnson-Cook

1．5 SPH 耦合FEM 的噴丸模型

采用LS-DYNA 軟件建立噴丸模型如圖1 所示．其中彈丸流使用兩類SPH 粒子建模，共包含320 個SPH 粒子;工件使用有限元建模，為減少計算時間，工件尺寸取6 mm×6 mm×4 mm，由于工件中部區(qū)域與彈丸發(fā)生接觸，故對其網(wǎng)格進行細化，由53 760 個solid 單元組成．工件底面采用固定位移邊界，4 個側(cè)面采用非反射邊界以避免應力波在其邊界發(fā)生反射．

通過定義接觸算法耦合FEM 與SPH 模型．在LS-DYNA 中通過關鍵字Contact_Automatic_Nodes_to_Surface 定義有限單元模型上表面與SPH 粒子間的接觸關系，其中將FEM 單元節(jié)點定義為主節(jié)點，SPH 粒子定義為從節(jié)點． SPH 耦合FEM 的過程如圖2 所示．

圖2 SPH 耦合FEM 的過程Fig．2 The process of SPH coupled FEM

2 仿真結(jié)果分析值模型

圖3 為大量彈丸以速度50 m/s 撞擊工件時不同時刻的撞擊效果圖． 基于此模型可模擬大量彈丸撞擊工件的噴丸過程．

圖3 大量彈丸撞擊工件效果圖Fig．3 The process of shots impacting the target

2．1 彈丸速度對殘余應力的影響

由于本文仿真模型中，彈丸流由大量隨機分布的丸粒組成，可通過彈坑附近的殘余壓應力研究大量丸粒撞擊下的噴丸結(jié)果． 圖4 給出彈丸速度分別為25，50，75 和100 m/s 的大量彈丸經(jīng)過120 μs 撞擊后，彈坑底部殘余應力分布圖． 對應上述4 種彈丸速度，其最大殘余應力分別為-632．2，-692．71，-729．0，-738．9 MPa，發(fā)生在距離彈坑表面0．18 ～0．24 mm 的次表層．最大殘余壓應力值和其深度均隨速度的增加而增加．因此在實際生產(chǎn)中，應根據(jù)最大殘余壓應力及應力層深度要求合理選擇噴丸速度．

2．2 噴丸時間對殘余應力的影響

基于本文含大量丸粒的宏觀噴丸模型，可研究殘余應力場隨噴丸時間的變化． 隨著噴丸時間的增加，撞擊工件表面的彈丸數(shù)量、工件表面的有效受噴面積及覆蓋率均隨之增加． 圖5 所示為速度為50 m/s 的彈丸分別經(jīng)過26．6，79．8，102．6，120 μs 時間后，彈坑底部殘余應力分布圖． 圖示彈坑處受多個不同粒子撞擊，其中26．6 μs 對應彈丸首次撞擊后的殘余應力結(jié)果，而79．8 μs 則對應彈丸多次撞擊同一位置后的殘余應力結(jié)果．兩者對應的最大殘余壓應力值分別為-658．86 MPa 和-692．17 MPa，其最大殘余壓應力值有所提高．撞擊時間繼續(xù)增加至102．6 μs 和120 μs時，其對應的殘余應力變化甚微，最大變化幅度小于10%，即噴丸殘余應力存在飽和狀態(tài)．

2．3 不同位置點的殘余應力分布

圖4 噴丸速度對殘余應力的影響Fig．4 Effect of impact velocity on residual stress

圖5 噴丸時間對殘余應力的影響Fig．5 Effect of impact time on residual stress

圖6 典型彈坑位置點Fig．6 Typical positions of craters

圖6 所示為一典型四彈坑分布圖． C、E、G、I分別為4 個彈坑的中心，B、D、H、F 分別為4 個彈坑兩兩重疊區(qū)域，A 點未受到撞擊． ABCD 區(qū)域可視為一典型區(qū)域．圖7 分別表示A、B、C、D 4 個位置點的殘余應力分布情況． C 點為粒子撞擊中心點，其最大殘余壓應力為-653．45 MPa，發(fā)生在深度為0．20 ～0．22 mm 的次表層;B 點和D 點位于兩彈坑重疊處，該處的最大殘余壓應力值為-343．56 MPa，位于深度為0．12 ～0．18 mm 的次表層;A 點的最大殘余壓應力值為-450 MPa，發(fā)生在表面處．

圖7 不同位置殘余應力分布Fig．7 Distribution of the residual stress vs． depth at different positions

3 結(jié)論

針對有限元法難以模擬大量彈丸群集撞擊的噴丸過程，采用SPH 耦合FEM 法進行了噴丸過程數(shù)值建模與分析．仿真結(jié)果如下:

(1)最大殘余壓應力值及所在次表層深度均隨彈丸速度的增加而增加． 故應根據(jù)最大殘余壓應力及應力層的深度要求合理選擇噴丸速度．

(2)開始時，隨著噴丸時間的增加，最大殘余應力值增大;達到一定時間后，繼續(xù)增加噴丸時間，殘余應力變化甚微，最大殘余壓應力的變化幅度小于10%，達到噴丸飽和狀態(tài)．

(3)彈坑中心點處最大殘余壓應力值較高且位于次表層;兩坑交匯處的最大殘余壓應力值較中心點低且其次表層深度較淺;未受撞擊點處的最大殘余壓應力發(fā)生在表面上．

［1］ KAZUYUKI O． Fatigue life enhancement of aluminum alloy for aircraft by Fine Particle Shot Peening(FPSP)［J］． Journal of Materials Processing Technology，2011，211(8):1395 -1399．

［2］ 高玉魁．噴丸強化對23CO14Ni12Cr3MoE 鋼疲勞性能的影響［J］． 材料熱處理學報，2007，28(6):75 -77．

［3］ MORI K，OSAKADA K，MATSUOKA N． Finite element analysis of peening process with plastically deforming shot［J］． Mater Process Technol，1994，45(1 -4):607 -612．

［4］ MEGUID S A，SHAGAL G，STRANART J C． Finite element modelling of shot-peening residual stresses［J］． Journal of Materials Processing Technology，1999，92/93:401 -404．

［5］ MEGUID S A，SHAGAL G，STRANART J C，et al．Three-dimensional dynamic finite element analysis of shot-peening induced residual stresses［J］． Finite Elements in Analysis and Design，1999，31(3):179-191．

［6］ MAJZOOBI G H，AZIZI R，NIA A A． A three-dimensional simulation of shot peening process using multiple shot impacts ［J］． Journal of Materials Processing Technology，2005:164/165:1226 -1234．

［7］ 凌翔，彭薇薇，倪紅芳． 噴丸三維殘余應力場的有限元模擬［J］． 機械工程學報，2006，42 (8):182 -189．

［8］ 路會龍，姚平喜，劉海英．基于ANSYS/LS-DYNA 的受控噴丸工藝過程仿真［J］．機械設計與制造，2009(2):214 -216．

［9］ 張洪偉，張以都，吳瓊．噴丸強化殘余應力場三維數(shù)值分析［J］． 航空動力學報，2010，25(3):603-609．

［10］ GRUJICIC M，PANDURANGAN B，QIAO R，et al．Parameterization of the porous-material model for sand with different levels of water saturation［J］． Soil Dynamics and Earthquake Engineering，2008，28(1):20 -35．

［11］ DABBOUSSI W，NEMES J A． Modeling of ductile fracture using the dynamic punch test ［J］． International Journal of Mechanical Sciences，2005，47(8):1282 -1299．