高性能齿轮是重要的机械基础件，工况多为高速重载。为提高齿轮疲劳寿命，优化齿轮服役性能，一般在全流程制造末端加入喷丸工序。喷丸强化利用高速弹丸冲击零件表面，使工件受喷表面发生不均匀塑性变形，诱导残余压应力产生。较好的残余应力分布可 提高齿轮的接触疲劳强度、弯曲疲劳强度、抗胶合能力，从而改善零件服役性能。
采用喷丸工艺提高齿轮疲劳强度，需充分发挥喷丸强化因素影响效果，即最大程度提高残余压应力值。因此，了解喷丸后齿面残余应力分布与工艺参数关联规律，优化喷丸工艺参数是提高齿轮服役性能的关键。生产中利用试验摸索进行工艺编制，需耗费大量的时间和精力，开展喷丸工艺仿真研究可以节约资源，为工业制造提供指导。
现有喷丸残余应力预测模型主要基于有限元方法。Miao 等建立了与实际喷丸过程相近的随机多弹丸有限元模型，因此被广泛应用。基于此，李源等利用 Python 语言编写了 ABAQUS 子程序，模拟多弹丸随机入射 Q235 钢板的过程，结果表明多弹丸模型模拟的残余压应力峰值和试验结果很接近。Bagherifard 等利用试验对仿真结果进行了验证，并研究了部分工艺参数(如喷丸速度和弹丸尺寸等)与喷丸后残余应力的关系。目前喷丸工艺仿真模型主要针对平面或圆柱面试样，无需考虑弹丸之间的碰撞。但对于螺旋锥齿轮而言，其复杂曲面使弹丸间碰撞概率大大提高，弹丸碰撞成为仿真模型中必须考虑的因素。基于离散元的方法建立模型，可计算弹丸从喷嘴口到齿面的运动过程。结合离散元和有限元方法计算螺旋锥齿轮喷丸工艺强化效果，提高了仿真结果计算的准确性，受到广泛使用。
本文通过建立螺旋锥齿轮的喷丸工艺仿真模型，计算喷丸后轮齿齿面表层残余应力分布状态。对比喷丸试验结果与仿真结果，以验证仿真计算准确性。基于该仿真模型，研究喷丸时间、弹丸速度、弹丸直径与喷丸后齿面残余压应力场的关联规律。
一、喷丸工艺试验
试验样件：喷丸样件为某型号螺旋锥齿轮，材料为 AISI 9310 高强度合金钢，齿部经渗碳淬火处理。齿面为喷丸强化的重要区域，其表面至 0.56 mm 深度内，硬度均高于 60HRC。
试验参数：试验采用 MP1000Ti 型号数控喷丸机开展。弹丸类型采用 ASH110(55~62HRC)，喷丸强度为 0.178~0.228 mmA，覆盖率 200%、喷嘴角度 17°、空气压力 0.25 MPa、移动速度 70 mm/min、喷丸时间 144 s。对齿根、凸面、凹面分别进行强化，喷嘴与目标表面垂直，且距离 150 mm，齿轮转台转速 30 r/min，弹丸流量 5 kg/min。
残余应力测量：对齿面进行残余应力测量，测量点位置如图 1 所示，a、b、c 均位于齿面节锥线上，依次为齿宽的四分点，a 靠近轮齿大端。为得到齿轮表层残余应力场，沿齿面法向对轮齿进行电解抛光，电解抛光深度由白光干涉仪测量得到。
图 1 残余应力测量点(a、b、c)位置示意图
利用加拿大 Proto 公司的 X 射线衍射仪测量齿面残余应力。管电压 25 kV，管电流 5 mA，X 射线管 Cr_K– Alpha，光圈直径 1 mm，波长 2.291 A，曝光时间 3 s，曝光次数 7，最大 β 角 20°。测试照片如图 2 所示。电解抛光所用仪器为 8818–V3 型号电解抛光仪，设置电压 40 V，流量 8，抛光时间 3 s，抛光电流 2.8~3 A。
图 2 残余应力测试示意图
二、仿真模型建立
由于螺旋锥齿轮曲面复杂，弹丸在齿面附近易发生相互碰撞，从而影响弹丸实际冲击速度及次数。为提高计算准确度，本文基于离散元及有限元耦合的方法建立喷丸工艺仿真模型，计算螺旋锥齿轮喷丸后齿面表层残余应力分布，计算流程如图 3 所示。建立螺旋锥齿轮喷丸离散元模型，获取齿面弹丸冲击速度矢量、单位面积冲击次数。将离散元计算结果作为输入，结合喷丸前齿面三维粗糙度、残余应力场、齿轮及弹丸材料模型，建立有限元模型，计算喷丸后轮齿表层残余应力场，分析喷丸时间、弹丸速度、弹丸直径等对喷丸后齿面残余应力场的影响。
图 3 离散元与有限元耦合计算模型流程图
离散元模型
模型建立：利用 EDEM 仿真软件设置离散元模型，模拟弹丸从喷嘴到齿面的过程。
根据 1.2 节工艺参数设置模型参数，几何模型如图 4 所示。为提取轮齿表面目标位置冲击信息，利用 Hypermesh 软件对轮齿表面进行网格划分。在模型中，导入齿面单元作为齿轮模型，设置喷嘴分别垂直于齿轮的齿根、凸面、凹面，且距离目标表面 150 mm，并沿齿宽方向进行移动。弹丸初始速度由经验公式(1)确定。
式中，d 为弹丸直径，mm;qm 为丸流量，kg/min;P 为喷嘴气压，MPa。
图 4 DEM 几何模型示意图
数据提取与处理：利用 EDEM 软件导出计算结果，并利用单元编号判断冲击点位置，筛选出齿面节锥线上四等分点处的冲击，进行数据处理。
(1)冲击速度矢量。
由于齿轮为曲面，各位置与弹丸束的夹角不同，且有限元模型中靶板为齿面上的局部区域。因此，需将弹丸的绝对速度转换为齿面单元与弹丸的相对速度。
如图 5 所示，曲面 OABC 为齿面上某一区域，划分为网格后简化为四边形 OABC。以 O 为坐标原点，OA 为 x′ 轴，OC 为 y′ 轴，根据空间坐标系右手法则确定 z′ 轴，建立笛卡尔坐标系。其坐标轴的单位向量在大地坐标系中表示为
。在大地坐标系中，弹丸 i 的速度矢量表示为 vi(xi ， yi ，zi ，0)，在新坐标系中表示为 vi ′(xi ′，yi ′，zi ′，0)。根据转轴变换原理得
图 5 弹丸速度场坐标系变化示意图
(2)冲击次数。
当实际喷丸时间为 t1 时，利用式(3)计算目标位置的单位面积冲击次数 n。
式中，n1 为该单元与弹丸发生的冲击次数；S1 为该单元面积；tdem 为 EDEM 模型计算时间。
有限元模型
基于 ABAQUS/CAE 商业有限元软件建立有限元模型，计算喷丸后轮齿齿面的残余应力场。
网格划分及边界条件： 弹丸直径为 0.3 mm，单元类型为 C3D8R。靶板 尺寸如图 6 所示，除无限单元体网格外，靶板尺寸为 1 mm×1 mm×0.5 mm，单元类型为 C3D8R，在该区域内设置预定义场，以喷丸前齿面残余应力实测值作为模型的初始残余应力。其中区域Ⅰ用于限制弹丸冲击中心点位置和提取残余应力计算结果。区域Ⅰ和区域Ⅱ为 网格细化部分，单元尺寸均为 10 μm×10 μm×10 μm。区域Ⅲ为有限元网格由细到粗的过渡区。区域Ⅳ的网格类型为无限单元体 CIN3D8，用于消除应力波在靶板边界的反射。
图 6 冲击靶板几何示意图
模型中对靶板底面进行完全固定。弹丸与靶板之间的接触关系为Surface to Surface，法向行为定义为“硬”接触，切向行为定义为罚摩擦，摩擦系数设置为 0.2。
材料模型：弹丸设置为弹塑性体，采用各向同性本构模型，弹丸参数设置为杨氏模量 E=210 GPa，泊松比 μ=0.3，密度 ρ=7800 kg/m3 ，屈服强度 σs=1400 MPa。靶体材料为 9310 渗碳淬火钢，材料参数设置为杨氏模量 E=210 GPa，泊松比 μ=0.3，密度 ρ=7800 kg/m3 ，塑性应力应变曲线采用 Johnson–Cook 模型，其基本模型如式(4)和(5)所示， 参数由 Hopkinson 压杆试验获得，分别为 A=1234.38，B=881，C=0.018，n=0.238，m=0.686。为防止应力波振荡，设置靶板材料阻尼 α=6×106 s-1。
式中，A、B、n、C、m 为材料物理特性参数，由压杆试验结果拟合得到;A 为参考应变率下的屈服强度;B 为应变硬化系数;n 为应变硬化指数;C 为应变率硬化系数;m 为温度软化指数; 
为试验应变率;
为参考应变率;T 为试验温度;Tmelt 为材料的熔点(9310 合金钢的熔点为 1510 ℃)；Troom 为参考温度，值为 25 ℃。
弹丸冲击信息：弹丸冲击角度、速度、次数由 2.1.2 节计算结果确定。根据 ASH110 型号弹丸尺寸标准，设置弹丸为球形，直径为 0.3 mm。考虑实际情况下弹丸冲击位置随机分布，仿真中利用 Python 语言编程，采用 Random 函数产生随机位置，弹丸冲击中心点均在图 6 所示的区域Ⅰ内。
三、结果与讨论
仿真与试验结果对比：为验证仿真模型准确性，在所建立的喷丸工艺仿真模型中，设置与试验相同的工艺参数。图 7 为喷丸后凸面 b 点实测和仿真三维形貌。根据式(6)计算三维粗糙度 Sa 值，实测值为 0.35 μm，仿真计算结果 0.386 μm，误差 10.3%。其中，在 x、y 方向，仿真结果分辨率低于实测结果，原因在于实测数据的采样间距为 0.5 μm，有限元模型的靶板网格尺寸为 10 μm×10 μm。
式中，n 为数据点个数；Zi 为第 i 个节点的高度值。
图 7 喷丸后齿面三维形貌
对比喷丸后实测与仿真计算的轮齿表面齿高方向残余应力值，如表 1 所示。其中，a、b、c 点位置与图 1 标注一致。结果表明，轮齿各点的表面残余应力均在 –800~–880 MPa 之间，计算误差在 6% 以内。图 8 展示了凸面、凹面 b 点的实测与仿真计算的残余应力沿深度方向的分布，其中 σx 为齿宽方向残余应力;σy 为齿高方向残余应力。喷丸处理后，轮齿凸面、凹面残余应力变化趋势一致，凹面残余压应力值略大于凸面。且齿宽方向残余应力和齿高方向分布相近，表面残余应力约为 –800~–850 MPa，最大残余压应力值在 1200~1300 MPa 之间，深度约为 20~30 μm。统计轮齿表面残余压应力和最大残余压应力值，如表 2 所示，实测值与仿真结果误差均小于 10%，证明该仿真模型可准确预测螺旋锥齿轮喷丸后的轮齿表层残余应力场。
表 1 试验与仿真计算的轮齿表面齿高方向残余应力对比
图 8 凸面、凹面 b 点的试验与仿真计算残余应力场对比
喷丸时间影响：喷丸覆盖率受丸流量和喷丸时间影响，丸流量和喷丸时间的增加都会提高弹丸冲击次数，增加喷丸覆盖率。本文主要研究喷丸时间对残余应力场的影响，结果如图 9 所示。基于 Miao 等的方法，对喷丸覆盖率进行了计算。
图9 喷丸时间对残余应力结果的影响
在弹丸初始速度为 40 m/s 的条件下，喷丸时间达到 72 s 时，弹丸冲击次数为 360 个/mm2 ，喷丸覆盖率达到 98%，即满覆盖率。当喷丸时间为 144 s 时，覆盖率为 200%。随着喷丸时间的增加：(1)靶板表面残余压应力升高，当喷丸时间达到 72 s 后，逐渐稳定在 800 MPa 左右；(2)10 μm 处残余压应力呈现逐渐升高的趋势，并在喷丸时间达到 96 s 后稳定，约为 1000 MPa；(3)20 μm 处残余压应力持续提高，在 72 s 后，20 μm 处残余压应力高于 10 μm 处，表明最大残余压应力深度由 10 μm 提高至 20 μm。这表明，提高喷丸时间，会有效增加材料次表层残余压应力值，提高最大残余压应力层深。
表 2 凸面、凹面 b 点试验与仿真计算残余应力结果对比
弹丸速度影响：由式(3)可知，弹丸初始速度由喷嘴气压、丸流量、弹丸直径确定，其中喷嘴气压为主要因素。喷嘴气压的提升，会使弹丸速度明显增加，提高喷丸强度。
参考 Nordin 等的试验弹丸速度范围设置弹丸初始速度，残余应力计算结果如图 10 所示，云图如图 11 所示。结果表明，弹丸速度由 30 m/s 提高至 40 m/s、50 m/s 时，喷丸诱导的残余压应力提高，表面残余压应力值变化较小，最大残余压应力值分别为 –1104.7 MPa、–1144.9 MPa、–1167.3 MPa，最大残余压应力深度 分别为 20 μm、25 μm、30 μm。结果表明，随着弹丸速度的提高，靶板表面残余压应力值变化不大，但最大残余压应力值、最大残余压应力深度会不断提高。
弹丸直径影响：计算不同弹丸直径冲击下的残余应力场如图 12 和 13 所示。其中弹丸直径根据 SAE 标准设置，采用 ASH70、ASH110、ASH170 类型弹丸，表 3 展示了各类型弹丸与其直径的对应关系。结果表明，弹丸直径提高时，表面残余应力值变化较小，但会明显提高最大残余压应力值和最大残余压应力深度。当弹丸直径由 0.18 mm 提高至 0.30 mm 和 0.42 mm 时，表面残余压应力均在 –760~–840 MPa 之间，最大残余压应力值由 –893.6 MPa 提高至 –1145.0 MPa 和 –1251.5 MPa，最大残余压应力深度由 10 μm 提高至 30 μm 和 40 μm。但弹丸直径提高会增加齿轮表面粗糙度。其原因在于当弹丸直径较高时，弹丸动能增加，与齿面冲击时，材料塑性应变增加，从而增大弹坑半径，会使齿轮表面粗糙度明显提高。
表 3 不同类型弹丸对应直径及残余应力计算结果
四、结论
本文基于离散元与有限元耦合的方法，建立 AISI 9310 钢螺旋锥齿轮的喷丸工艺残余应力场仿真计算模型，探寻喷丸工艺参数与喷丸后齿面残余压应力的关联规律。
(1)本文所建立的仿真模型对齿面残余应力场的预测误差在 10% 以内，计算结果准确。
(2)当喷丸时间为 72 s 时，弹丸冲击个数达到 360 个/mm2 ，靶板覆盖率达到满覆盖率。当喷丸时间继续增加时，靶板表面压应力逐渐稳定在 –800 MPa 左右。
(3)喷丸工艺主要影响轮齿表面至 50 μm 深度内的残余应力场，深度大于 50 μm 时，残余应力场由喷丸前工艺确定。
(4)当喷丸覆盖率为 200% 时，弹丸直径和速度的改变对本文所用零件的表面残余压应力影响较小。
(5)当弹丸初始速度和直径提高时，靶板表面残余压应力值变化较小，靶板表层最大残余压应力值和最大残余压应力深度都得到明显提高，其中最大残余压应力值可提高到 –1251.5 MPa，最大残余压应力深度可提高到 40 μm。
参考文献略.