基于SPH-FEM的超声珩磨空化微射流冲击研究

叶林征， 祝锡晶， 王建青， 郭　策

(中北大学 机械与动力工程学院，太原　030051)

基于SPH-FEM的超声珩磨空化微射流冲击研究

叶林征， 祝锡晶， 王建青， 郭策

(中北大学 机械与动力工程学院，太原030051)

为探究超声珩磨中空化微射流对壁面的冲击作用，考虑壁面弹塑性变形、超声场及珩磨压场，采用SPH-FEM耦合方法建立了空化微射流冲击模型并进行了数值仿真分析，随后进行了试验验证，结果表明：微射流冲击过程中形成的侧向射流速度高于微射流冲击速度，最高可达冲击速度的1.6倍，冲击中后期微射流中部粒子反向运动向上凸起；壁面出现直径约8 μm、深约0.173 μm的微型凹坑，且其边缘处有材料塑性隆起，试验中材料表面出现微型凹坑群；冲击区域内由内而外等效应力先升高后降低，射流边缘附近冲击效应最强。从材料去除角度来看，大量微射流冲击对超声珩磨加工起到积极作用。

流体力学；微射流冲击；光滑粒子流体动力学；超声珩磨；声空化

功率超声珩磨是一种在普通珩磨中加入超声谐振系统的精密加工方法，广泛用于坦克、汽车、摩托车等内燃机气缸内壁加工中，为高硬、高强材料的光整加工提供了新途径[1]。超声珩磨加工中，为冷却及润滑，会向珩磨区中加入大量切削液，珩磨头的几何结构及其往复及轴向旋转运动易于产生空化效应，原理同螺旋桨叶片腐蚀类似，加之超声作用，珩磨区中会产生强烈的超声空化效应。切削液中微小泡核在超声波正负压相的交替作用下，会历经生长、膨胀、压缩、溃灭等一系列行为[2]。空泡的非线性振动及溃灭释放的微射流及冲击波均会对加工壁面产生影响，相关研究表明微射流起主导作用[3]。目前，许多学者在空化领域或射流冲击领域进行研究，但关注空化微射流冲击方面的研究非常少，而研究超声珩磨中空化微射流对壁面的作用具有重要的现实及工程指导意义。

在功率超声珩磨中，靠近壁面的空泡溃灭具有不对称性，远离壁面一侧的空泡自由界面收缩较快，空泡从远离壁面一侧向底部凹陷并最终贯穿，同时产生一束高速微型射流冲击工件壁面。空化微射流冲击壁面属于高速冲击问题，也是一类强非线性流固耦合问题，相似现象还有汽轮机叶片的水滴侵蚀[4]、水下爆炸产生的射流冲击舰船底板[5]、水射流冲击[6]等等，均会涉及液体区域(即微射流)的大变形以及壁面材料的变形损伤，传统的有限元方法(FEM)在处理大变形问题时会出现网格畸变，计算精度不高甚至难以计算。光滑粒子流体动力学(Smoothed Particle Hydrodynamics,SPH)方法，是一种完全拉格朗日无网格粒子方法，用一套粒子代替有限元法中的单元，对积分或微分形式的控制方程组进行离散及近似，避免了计算中由于网格存在而产生的困难，在大变形、瞬态冲击爆炸、网格畸变等问题的处理上具有显著优势，但在处理小变形问题时不如有限元方法精确。为保留两种方法各自的优势，建模过程中可分别建立SPH粒子模型和有限元模型，并进行耦合分析，此即SPH-FEM方法。SPH方法最早由Lucy等[7]提出，主要应用于天体物理学领域，随后广泛地应用于流体力学及连续固体力学相关领域。Libersky等[8]率先在SPH方法中加入材料强度效应，成功地对高速碰撞行为进行了数值模拟计算。随后，Johnson等[9]将SPH方法及有限元方法结合，进行了侵彻贯穿问题的数值计算，取得了较大进展。马利等[10]采用SPH耦合有限元方法模拟了水射流切割现象，为水射流加工工作参数的优化提供了参考数据。

本文通过SPH-FEM耦合方法模拟超声珩磨中空化微射流对工件壁面的冲击过程，建立SPH-FEM耦合模型并进行相关分析，得到了冲击过程中流场速度的变化，并对工件壁面的变形、速度、压强、应力应变等参数进行了仿真分析。仿真结果有助于理解超声空化微射流对壁面的作用机理。

1　SPH基本理论

SPH方法是基于粒子的，所有物理量均由粒子携带。用SPH方法近似偏微分方程组，分为核近似和粒子近似两步[11]。核近似即应用核近似法对任意函数与光滑核函数进行积分，粒子近似即在有限区域内对所有粒子的值进行叠加求和。SPH理论中，变量场Ω内任意连续函数f(x)的核近似表达式为[11]：

〈f(x)〉=∫Ωf(x′)W(x-x′,h)dx′

(1)

式中，h为光滑长度，定义了光滑核函数W(x-x′,h)的影响区域。光滑核函数的选取尤为重要，不同的光滑核函数计算精度不同，本文选用应用广泛的三次B样条光滑核函数。

本研究中SPH粒子区域为液体，用于求解的Navier-Stocks方程经SPH方法离散后为[12]：

(2)

(3)

(4)

2　SPH-FEM耦合模型建立

2.1模型建立

超声珩磨珩磨区中会产生超声空化效应，空泡溃灭可产生速度达几百米每秒的微射流[13]，数量极大的微射流不间断地冲击工件壁面，造成材料的弹塑性变形甚至材料去除。冲击过程中涉及到微射流的大幅度变形以及材料的微小弹塑性变形，为使计算结果精确，采用SPH-FEM耦合方法，选择Abaqus大型商用有限元分析软件对其进行建模分析。建立的微射流冲击模型如图1所示。微射流采用SPH粒子建模，用关键字*Section Controls控制SPH计算相关参数，如人工黏性项、光滑长度、光滑函数类型、SPH区域等等。壁面采用有限单元法划分网格。采用关键字*Contact Inclusions建立SPH粒子模型和有限元模型的接触关系，借以进行耦合分析。

图1　微射流冲击SPH-FEM耦合模型Fig.1 SPH-FEM coupling model of micro-jet impact

2.2模型初始条件

根据超声珩磨加工的具体工况，在建立模型中，忽略了液体的黏性、压缩性、表面张力及介质间的热效应，考虑了壁面的弹塑性变形、超声场及珩磨压场。珩磨液一般为水、煤油或乳化液，粗珩时通常选水，根据加工中壁面表面粗糙度的范围假设其与壁面的摩擦因数为0.1。空化微射流直径约为几微米[14]，选取模型参数为：微射流直径6 μm，长13 μm，以300 m/s的速度垂直冲击壁面，单元类型为PC3D，粒子总数5 544个；壁面为铝板，长40 μm，宽40 μm，高10 μm，单元类型为C3D8R，在受冲击中心区域加密网格，加密单元32 000个，非加密单元12 000个。

液体采用Us-Up状态方程进行控制，即Us=c0+sUp，c0和s分别为液体声速及无量纲参数，定义了冲击速度Us和粒子速度Up间的线性关系,s=1.75[15]。固体选取弹塑性材料且四周固定，在建模区域内施加周期性变化的正弦函数型超声压pa及珩磨压pH[16]，pH大小0.24 MPa，pa=pAsinωt，ω=2πf，pA为声压幅值1.753 MPa，f为超声频率18.6 kHz。材料参数见表1。

表1　材料参数

3　SPH-FEM耦合仿真分析

单次微射流冲击壁面时间极短，只有几十纳秒，选取t=0、10、20、30、40、50 ns时的冲击图像来表示微射

流冲击壁面的过程，如图2所示，本文分析中涉及到的坐标正负也如图中所示。图中可直观地看出粒子的平铺、散落、飞溅等特征。在冲击的中后期，可观察到射流中心轴附近的粒子反方向凸起。仿真分析主要分为微射流区域粒子速度的仿真分析以及工件壁面压强、变形及应力应变等的仿真分析。

3.1微射流粒子速度的数值仿真

在微射流头部圆弧上，从中心点向右依次取5个粒子A、B、C、D、E，如图3所示，粒子编号分别为8、233、4、2、1，其速度大小时历曲线如图4所示。微射流中心粒子，如A、B粒子，冲击壁面后速度急剧降低并始终保持在较小速度最终趋近于0；而射流外侧粒子，如D、E粒子，在到达壁面后，速度先急剧升高，最高可达约477 m/s，随后降低至300 m/s左右保持相对稳定状态并持续较长时间。这是由于射流中心粒子受周围粒子作用，各方向作用力抵消且没有空间运动，故速度急剧降低；而外侧粒子接触壁面时，由于撞击液体内部压强高于外部环境压强，粒子在高压下喷出，形成侧向射流，故外侧粒子侧向流动速度高于微射流冲击速度，最大可达到冲击速度的1.6倍。

图2　不同时刻冲击图Fig.2 Impact figure of different times

图3　微射流粒子点Fig.3 Micro-jet particle point

在射流中心轴线上，自中心点向上依次选取a、b、c、d、e五个粒子，粒子编号依次为8、1 328、2 648、3 968、5 288，其v3时历曲线如图5所示，v3指Z向速度。五个粒子依次撞击壁面随后速度持续振荡降低至0，而射流中部粒子，如b、c、d，在速度降为0后又反向增大。这可能是因为微射流粒子持续冲击壁面，在其交界面上不断有激波产生，并向微射流中传播，带动粒子反方向运动，这也解释了图2中冲击中后期粒子反向凸起的现象。

图4　A、B、C、D、E速度大小时历曲线Fig.4 Velocity curves of A,B,C,D and E over time

图5　a、b、c、d、e v3时历曲线Fig.5 v3 curves of a,b,c,d and e over time

3.2壁面数值仿真

在本仿真分析中，壁面设置为弹塑性材料，受微射流冲击时会发生弹塑性变形，如图6所示。在微射流冲击后，壁面产生了微小的近似圆形凹坑，并在其边缘处有塑性隆起。绘制t=1、2、5、10、50 ns时网格加密区域中心线上的凹坑形貌图，如图7所示。可见凹坑直径约为8 μm，深度约为0.173 μm。凹坑最大深度出现在冲击边缘而不在冲击中心，这与其他学者研究结果一致[17]，证实了SPH-FEM方法分析微射流冲击问题的可行性。

图6　壁面形貌Fig.6 Wall topography

壁面中线上，由中心点开始向右依次选取F、G、H、I四点，见图3，节点编号依次为17 241、17 239、17 237、17 235，其v3时历曲线如图8所示。可见，速度变化主要发生在12 ns之前，随后仅在0附近做微小振动。受冲击质点挤压周围质点从而受到周围介质的反作用力或者应力波反射等因素，有的质点速度会达到正值；有的质点处于塑性隆起区域内，速度会出现沿Z轴正方向增大，如节点I，但最终速度均趋近于0。

F、G、H、I四点的时历压强见图9，壁面质点在受到微射流冲击瞬间压强急剧增大至极大值，但只持续极短时间，随后急剧降低至一相对稳定状态并持续较长时间，最终降低至0。受微射流冲击中心区域压强并非最大，如节点F、G；最大压强出现在微射流冲击边缘附近，如H、I节点，间接论证了凹坑最大深度出现在冲击边缘。为更直观地看出壁面最大压强分布，绘制壁面最大压强分布图，即图10。

图7　不同时刻壁面变形曲线Fig.7 Curves of the wall deformation at different moments

图8　F、G、H、I v3时历曲线Fig.8 v3 curves of F,G,H and I over time

图9　F、G、H、I压强时历曲线Fig.9 Pressure curves of F,G,H and I over time

图10　壁面最大压强分布Fig.10 The deformation of the wallmaximum pressure

采用von mises屈服准则，探求微射流冲击下壁面的等效应力分布。图11为t=1,2,5,7,9,10,15,20 ns时的壁面等效应力等值线分布图，冲击前期等效应力成环形分布，不断向外传播且由内而外逐渐减小，从7 ns开始，等效应力分布逐渐复杂，且冲击中心区域出现等效应力降低现象，最终形成由内而外等效应力先升高后降低的分布层次，且最大等效应力也出现在冲击区域边缘附件。

图11　不同时刻壁面等效应力等值线图Fig.11 Contour map of the wall equivalent stress at different moments

由以上数值仿真分析可知，在超声珩磨中空化效应产生的微射流冲击会造成工件壁面的弹塑性变形。虽然单次微射流冲击的影响很小几乎可以忽略，但是单次微射流冲击的时间极短，仅有几十纳秒，在超声珩磨加工过程中，会有数量极为庞大的微射流不间断地冲击壁面，此时，微射流的作用不可忽略，其冲击及剪切作用甚至会导致材料的去除。超声珩磨的加工效率约为普通珩磨的1.5倍[18]，或许超声珩磨中强烈的空化效应是造成此结果的一个因素，从这方面来讲，微射流冲击对超声珩磨中材料的去除有积极影响。

4　试验验证

在超声珩磨加工过程中，微射流数量极多，且存在时间极短，很难通过试验进行定量分析，为验证微射流对壁面的冲击作用，设计微射流冲击铝箔纸试验进行定性分析。试验原理图见图12，试验仪器主要有250 W的H66MC超声波发生器、功率超声珩磨谐振系统、玻璃水槽、铝箔纸、干燥箱、VHX-600ESO数码显微镜等等。将铝箔纸平整地粘贴在光滑的钢板上，调节超声波发生器的频率使其与功率超声珩磨谐振系统发生共振，将铝箔纸放置在水槽中距油石条1 mm左右位置，经过5 s后取出，放入干燥箱中去除水分，待其完全干燥后，取出放在VHX-600ESO数码显微镜下观察。

图13、14分别为试验前后铝箔试件的表面形貌图，可见试验后铝箔表面出现了微小的近似圆形凹坑，且其尺度均处于微米级。图14(a)中较明显的凹坑可认为是单次或多次微射流冲击同一点形成的，14(b)图中则显示了微小凹坑群，显示了大量微射流随机冲击后壁面的表面形貌。

通过试验直观地观察了微射流对壁面材料的作用，从材料去除的角度，可认为空化微射流对超声珩磨加工有一定的积极作用。

图12　试验示意图Fig.12 Test schematic

图13　试验前铝箔表面形貌Fig.13 Surface morphology of aluminum foil before the test

图14　试验后铝箔表面形貌Fig.14 Surface morphology of aluminum foil after the test

5　结　论

本文应用SPH-FEM耦合方法，将微射流设置为SPH粒子，壁面设置为有限元单元，基于Abaqus对超声空化微射流冲击壁面的过程进行了仿真分析，得到了微射流粒子速度及壁面变形、速度、压强、应力应变等的分布及变化过程，并设计微射流冲击铝箔纸试验进行验证，结论如下：

(1) 运用SPH-FEM耦合方法分析超声珩磨过程中的微射流冲击现象，可形象地描述冲击过程中液固各区域的变化，相关时历过程及数值仿真结果与试验结果一致，说明该方法是可行和有效的。

(2) 在超声珩磨过程中，微射流冲击壁面形成侧向射流，侧向射流的速度大于微射流冲击速度，最大约为微射流冲击速度的1.6倍，在微射流冲击中后期，中部粒子会出现反向运动向上凸起的现象。

(3) 在超声珩磨过程中，微射流冲击壁面瞬间，压强急剧增大至极大值而后迅速降低至一相对稳定状态并持续较长时间。壁面出现直径约8 μm、深约0.173 μm的微型凹坑，且其边缘处有材料隆起。冲击区域内由内而外等效应力先升高后降低，射流边缘附近冲击效应最强。

(4) 在超声珩磨过程中，大量微射流不间断冲击下，会使壁面出现微型凹坑群甚至材料去除，从材料去除角度，空化微射流冲击对超声珩磨加工有一定的积极作用。

[1] 祝锡晶．超声光整加工及表面成型技术[M]．第一版．中国：中国科学文化出版社，2005: 9-13.

[2] 张敏，李晓谦，蒋日鹏，等．超声辐射杆振幅分布与空蚀区域[J]．振动与冲击，2014, 33(13): 59-62.

ZHANG Min, LI Xiao-qian, JIANG Ri-peng, et al. Amplitude distribution and cavitation corrosion field of ultrasonic transducer system [J]. Journal of Vibration and Shock, 2014, 33(13): 59-62.

[3] 陈大融．空化与空蚀研究[J]．中国基础科学，2010, 6: 3-7.CHEN Da-rong. Cavitation and cavitation erosion[J]. China Basic Science, 2010, 6: 3-7.

[4] Li N, Zhou Q, Chen X, et al. Liquid drop impact on solid surface with application to water drop erosion on turbine blades, Part I: Nonlinear wave model and solution of one-dimensional impact [J]. International Journal of Mechanical Sciences, 2008, 50(10): 1526-1542.

[5] 孙龙泉，张忠宇，丁建军，等. 可压缩流体射流冲击平板的压力特性研究[J]. 振动与冲击，2013, 32(3): 178-181.

SUN Long-quan, ZHANG Zhong-yu, DING Jian-jun, et al. Pressure characteristics of panels subjected to collision of compressible jet [J]. Journal of Vibration and Shock, 2013, 32(3): 178-181.

[6] 穆朝民，韩靖．高压水射流冲击煤体的力学特征[J]．爆炸与冲击，2015, 35(3): 442-448.

MU Chao-min, HAN Jing. Mechanical characteristics of high-pressure water jets impinging on coal [J]. Explosion and Shock Waves, 2015, 35(3): 442-448.

[7] Luey L B. Numerical approach to testing of fission hypothesis[J]. Astronomical Journal, 1997, 82(12): 1013-1024.

[8] Libersky L D, Petschek A G. Smoothed particle hydrodynamics with strength of materials[J]. Advances in the Free Lagrange Method, Lecture Notes in Physics, 1990, 395: 248-257.

[9] Johnson G R, Petersen E H, Stryk R A. Incorporation of an SPH option into the EPIC code for a wide range of high velocity impact computations [J]. International Journal of Impact Engineering, 1993, 14(1/2/3/4): 385-394.

[10] 马利，陶伟明，郭乙木，等．SPH耦合有限元方法的水射流弹塑性碰撞模拟[J]．浙江大学学报:工学版，2008,42(2):259-263.

MA Li, TAO Wei-ming, GUO Yi-mu, et al. Elastic/plastic impact simulation of water jet using smoothed particle hydrodynamics and finite element method [J]. Journal of Zhejiang University:Engineering Science, 2008, 42(2): 259-263.

[11] 刘谋斌，宗智，常建忠．光滑粒子动力学方法的发展与应用[J]．力学进展，2011, 41(2): 217-234.

LIU Mou-bin, ZONG Zhi, CHANG Jian-zhong. Developments and applications of smoothed particle hydrodynamics [J]. Advances in Mechanics, 2011, 41(2): 217-234.

[12] 王安文，徐绯，张岳青．SPH方法在液固撞击数值模拟中的应用[J]．计算物理，2012, 29(4): 525-533.

WANG An-wen, XU Fei, ZHANG Yue-qing. SPH method in numerical simulation of liquid-solid impact [J]. Chinese Journal of Computational Physics, 2012, 29(4): 525-533.

[13] Brujan E A, Ikeda T, Matsumoto Y. On the pressure of cavitation bubbles [J]. Experimental Thermal and Fluid Science, 2008, 32(5): 1188-1191.

[14] Hammitt. Mechanical cavitation damage phenomena and corrosion-fatigue [M]. Defense Technical Information Center, 1971.

[15] 李晓杰，张程娇，王小红，等. 水的状态方程对水下爆炸影响的研究[J]. 工程力学，2014, 31(8): 46-52.

LI Xiao-jie, ZHANG Cheng-jiao, WANG Xiao-hong, et al. Numerical study on the effect of equations of state of water on underwater explosions [J]. Engineering Mechanics, 2014, 31(8): 46-52.

[16] 郭策，祝锡晶，王建青，等. 超声珩磨作用下两空化泡动力学特性[J]. 力学学报，2014, 46(6): 879-886.

GUO Ce, ZHU Xi-jing, WANG Jian-qing, et al. Dynamical behaviors of double cavitation bubbles under ultrasonic honing [J]. Chinese Journal of Theoretical and Applied Mechanics, 2014, 46(6): 879-886.

[17] Hsu C Y, Liang C C, Teng T L, et al. A numerical study on high-speed water jet impact [J]. Ocean Engineering, 2013, 72: 98-106.

[18] 祝锡晶. 功率超声振动珩磨技术的基础与应用研究[D]. 南京：南京航空航天大学，2007.

Cavitation micro-jet impact in ultrasonic honing based on SPH-FEM

YE Lin-zheng, ZHU Xi-jing, WANG Jian-qing, GUO Ce

(School of Mechanical and Power Engineering, North University of China, Taiyuan 030051, China)

To explore the effects of cavitation micro-jet impact on a wall in ultrasonic honing, considering elastic-plastic deformations of wall, ultrasonic field and honing pressure field, the SPH-FEM coupled method was used to establish the model of cavitation micro jet impact for numerical simulation, then a test validation was conducted. The results showed that the speed of lateral jet formed in the process of micro-jet impact is higher than the impact velocity of micro-jet, the former is up to 1.6 times of the latter; the central particles of micro-jet have a reverse movement and bulge upward; a micro pit with diameter of about 8 μm and depth of about 0.173 μm appears on the wall, and material plastically uplifts at the edge, micro dimple clusters arise on the material surface in the test; the equivalent stress firstly rises then drops from inside to out side within the impact area; the impact effect is the strongest near the edge of the jet; a large number of micro-jet impacts play a positive role in ultrasonic honing, from the view point of material removal.

fluid mechanics; micro-jet impact; smoothed particle hydrodynamics (SPH); ultrasonic honing; acoustic cavitation

10.13465/j.cnki.jvs.2016.13.012

国家自然科学基金资助项目(51275490)；中北大学研究生科技立项项目(20151211)

2015-10-20修改稿收到日期：2016-01-15

叶林征 男，博士生，1990年4月生

祝锡晶 男，博士，教授，1969年12月生

O427.4；V214.3+2

A