刘静; 唐昌柯; 师志峰; 王林峰; 徐子旦; 邵毅敏 【期刊名称】《《中南大学学报（自然科学版）》》 【年(卷),期】2019(050)010 【总页数】8页(P2417-2424)
【关键词】圆柱滚子轴承; 局部剥落故障; 过渡区; 接触特性; 振动特性 【作 者】刘静; 唐昌柯; 师志峰; 王林峰; 徐子旦; 邵毅敏
【作者单位】重庆大学 机械传动国家重点实验室 重庆 400044; 重庆大学机械工程学院 重庆 400044 【正文语种】中 文 【中图分类】TH133.33
滚动轴承是重要基础部件之一，广泛应用于工业、能源和航天等行业的关键机械装备。滚动轴承的运行状态和寿命直接影响机械装备的运转状态，30%的旋转机械故障是由轴承故障引起的[1-3]。滚动轴承的故障形式包括锈蚀、断裂、磨损、疲劳剥落和胶合等，许多研究者对滚动轴承滚道表面的局部剥落故障进行研究。当滚动轴承出现局部剥落故障时，会导致轴承刚度发生变化，使得轴承的振动水平增大[4-5]，甚至危害整套机械装备。因此，需要深入研究局部剥落故障诱发的滚动轴承振动特征。对滚子与局部剥落故障周围区域接触特性的深入认识将有助于理解轴承局部剥落故障的振动特征，为轴承的故障诊断和振动特征识别提供有益的参考。
针对滚动轴承部件局部剥落故障与滚子的接触特性问题，许多学者进行了大量研究工作。BRANCH等[6]利用试验和有限元的方法，研究了球轴承滚动体冲击剥落故障边缘时引起的临界应力和塑性应变。PATIL 等[7]建立了带有局部剥落故障的球轴承动力学模型，研究了局部剥落故障尺寸对轴承振动特征的影响规律。牛蔺楷等[8]以刚体动力学碰撞理论为基础，在分析滚动轴承运动学特性的基础上，对滚动体通过局部剥落故障时产生的冲击力进行了分析和计算。樊莉等[9]运用
ANSYS/LS-DYNA软件，建立了含局部剥落故障的滚动轴承有限元模型，研究了局部剥落故障对轴承滚动体与滚道之间接触应力的影响规律。张福星等[10]利用有限元方法，分析了深沟球轴承滚动体和内、外圈不同部位接触应力和变形的分布。张雪等[11]通过计算轴承接触表面及次表面应力等证明了有限元分析的正确性。伍济钢等[12]将轴承的表面损伤故障考虑为冲击函数分别加载到轴承的内外圈及滚动体上，计算了轴承振动和频谱分析。剡昌锋等[13]计算了考虑弹流润滑及滑动作用下滚动轴承系统局部缺陷位移激励的动力学响应。尽管上述研究采用集中质量法和有限元分析方法讨论了局部剥落故障过渡区对滚子与滚道之间接触应力及其诱发的振动特征，但是缺乏滚子通过轴承局部缺陷故障时滚子与滚道之间接触特性的相关研究。针对上述问题，本文作者采用有限元分析方法，建立圆柱滚子轴承滚子与滚道接触二维有限元等效模型，分析单个滚子进入和退出滚道局部剥落故障区域的过程中滚子与滚道之间的接触特性。研究局部剥落故障宽度对滚子与滚道之间的接触区宽度、接触变形、接触应力分布状态与接触刚度的影响规律，讨论局部剥落故障过渡区对轴承振动特性的影响规律，为深入认识局部剥落故障诱发的滚动轴承的振动特征提供参考。 1 模型描述
1.1 圆柱滚子与滚道等效接触模型
图1所示为局部剥落故障出现时圆柱滚子通过局部剥落故障区域的过程。图1中：
Q为作用于滚子上的径向载荷；v为滚子运动方向；O1和O2分别为滚子在通过局部剥落故障过程中的不同位置；δ1和δ2分别为不同位置的滚子与滚道之间的接触变形量；w1和w2分别为不同接触位置处滚子与滚道之间的接触宽度。本文假设滚道表面的局部剥落故障采用矩形轮廓进行表征[4]。滚子进入局部剥落故障区间时，滚子与故障边缘位置发生接触。由于故障边缘的曲率半径不同于正常滚道，故滚子与故障滚道的接触变形和滚子与正常滚道的接触变形量存在差异，滚子与故障滚道的接触宽度和滚子与正常滚道的接触宽度也存在差异。为了准确地描述滚子通过局部剥落故障过程中滚子与滚道的接触状态，忽略滚子在运动方向的速度冲击对接触的影响，建立了圆柱滚子轴承滚子与滚道接触的等效模型，如图2所示。 图1 圆柱滚子与滚道之间的位置关系Fig.1 Geometric relationship between roller and defect
图2 圆柱滚子轴承局部剥落故障等效模型Fig.2 Equivalent model of cylindrical roller bearing with localized spalling defect
圆柱滚子轴承滚子与滚道接触的等效模型[14]，如图2所示。从图2可知：滚道表面分为局部剥落故障区、过渡区与无故障区，矩形缺陷表示轴承局部剥落故障，XP为等效模型中轴承滚子与滚道接触的中心位置与轴承局部剥落故障中心位置的距离；Pmax为等效滚子中心位置的最大赫兹接触应力。等效模型的长度为10b，宽度为7b，其中，b为用赫兹接触模型计算得到的轴承接触半宽。接触半宽b表示为[15]
式中：F为作用在轴承上的径向载荷；R1为滚子半径；R2为内圈滚道半径；l 为滚子等效接触长度；E′为等效弹性模量。E′的表达式为
式中：E1和E2分别为滚子与滚道的弹性模量；ν1和ν2分别为滚子与滚道的泊
松比。与滚道对应的滚子等效半径用R表示，其值取为
1.2 建立有限元模型
根据圆柱滚子轴承滚子与滚道接触等效模型(见图3)，以圆柱滚子轴承N306为例，建立含局部剥落故障的滚子与滚道接触有限元二维模型。该模型的宽度和高度与图2所示等效模型一致。其中，滚子半径为R1，局部剥落故障的深度为d，在接触位置对网格进行加密处理，单元网格边长应不大于接触椭圆的半轴长度[16]，单元长度取为1.5 μm，远离接触位置的单元长度逐渐增大。约束模型底部边缘，将载荷竖直施加在滚子中心，滚子与滚道之间接触的目标单元与接触单元分别为TARGE169和CONTA172[17]。
图3 滚子与滚道接触有限元模型Fig.3 Finite element model for contact between roller and raceway of cylindrical roller bearing 1.3 有限元计算理论
滚子与滚道接触有限元模型的求解过程中，有限元模型的刚度矩阵、位移矩阵和载荷矩阵的方程表达式为[18]
式中：K，u 和F 分别为滚子与滚道之间的接触刚度、位移和载荷矩阵。在获取有限元模型的节点位移之后，根据几何方程和材料的应力-应变关系，可计算获得有限元模型各结点的应变和应力。
为了准确求解滚子与滚道接触有限元模型中的接触变形，必须选择合适的接触算法求解。有限元方法常用的接触算法包括：罚函数法、拉格朗日乘子法以及增强拉格朗日算法[19]。罚函数法是将接触非线性问题转化为材料非线性问题，在接触面之间利用接触弹簧单元建立连接关系，该弹簧刚度被定义为接触刚度。接触体的接触压力与接触变形尺寸之间的关系可表示为[17]
式中：p 为接触压力；Kn为接触刚度；un为接触变形尺寸。为了找到精确的拉格朗日乘子，增强拉格朗日算法对罚函数进行的一系列修正迭代，其接触体的接触压力与接触变形尺寸之间的关系可表示为[17] 式中：
式中：ɛ为输入容差；λi为第i次迭代的拉格朗日乘子分量。然而，拉格朗日乘子法则不需要定义接触刚度，仅需定义颤振控制参数。理论上，拉格朗日乘子法在接触闭合时强制零穿透，在发生黏着接触时强制定义为“零滑移”。但是，拉格朗日乘子法增加了模型的附加自由度，需要额外的迭代来稳定接触条件，会增加计算成本。
1.4 等效模型验证
滚子与滚道之间的接触算法采用增强的拉格朗日接触算法，将有限元等效模型与数值计算方法进行对比验证。对于正常钢制圆柱滚子轴承，常用计算滚子与滚道之间的接触变形的方法为Palmgren经验公式[20]，即：
式中：δ为滚子与滚道之间的接触变形。利用式(8)计算正常滚子与滚道之间的接触变形，将计算得到的结果与有限元方法计算得到的结果进行对比，结果如表1所示。
表1 有限元模型和Palmgren法计算得到的接触变形对比Table 1 Comparison of deformations between roller and raceway obtained from FEM and Palmgren’s method载荷/N 100 200 300 400 534 δ/μm有限元0.415 0.77 1.10 1.42 1.71 Palmgren 0.395 0.737 1.062 1.376 1.681接触变形差异/%4.8
4.3 3.5 3.1 1.7
从表1可知：当轴承径向载荷分布范围为100～534 N时，有限元计算得到的接触变形与Palmgren经验公式计算得到的接触变形之间的差异小于5%，故该有限元等效模型适用于研究圆柱滚子轴承滚子与滚道之间的接触特性。 2 局部剥落故障对轴承接触特性的影响
文献[15]的研究表明：对于圆柱滚子轴承，滚子与滚道之间接触宽度、接触变形和接触应力取决于滚子的径向载荷；接触宽度、接触变形和接触应力还与滚子半径、滚道半径和等效接触长度相关。然而，滚道表面出现局部剥落故障时，滚子与滚道之间的等效接触长度将出现变化，从而影响滚子与滚道之间接触宽度、接触变形和接触应力。
以圆柱滚子N306 研究对象，其几何参数如表2所示。该轴承的滚子直径为11 mm，内圈滚道半径为20.25 mm，滚子等效接触长度为11.4 mm，XP取为-1.5b～1.5b。本次计算选取25个等效接触位置进行研究。
表2 圆柱滚子轴承N306的几何尺寸参数Table 2 Geometric parameters of a CRB N306参数内圈滚道半径(R2)/mm外圈滚道半径/mm滚子半径(R1)/mm滚子等效接触长度(l)/mm滚子数目(Z)径向游隙(Cr)/µm数值20.25 31.25 5.5 11.4 12 1
2.1 局部剥落故障对接触应力和接触变形的影响
利用图3所示的有限元模型求解滚子与局部剥落故障周围区域之间的接触变形、接触宽度和接触应力的分布状态。选取滚子的外载荷为19.5 N，方向定义为竖直向下。选取－1.5b≤Xp≤1.5b，b 为式(1)中F 取19.5 N 时对应的接触宽度，滚子的外载荷保持不变。图4所示为局部剥落故障宽度分别为0.5b，0.75b 和b 时滚子通过局部剥落故障区域滚道的过程中滚子与滚道之间接触变形量、接触应力和接触宽度。从图4(a)可知：当局部剥落故障宽度为b时，随着滚子与局部剥落故障
边缘之间距离的减小，滚子与滚道之间的接触变形从0.563 μm 增大为0.812 μm，后又减小到0.682 μm，呈先增大后减小的趋势。从图4(b)可知：当局部剥落故障宽度为b时，随着滚子与局部剥落故障边缘之间距离的减小，滚子与滚道之间的接触应力从386 MPa 增大到1 644 MPa，呈增大的趋势。从图4(c)可知：当局部剥落故障宽度为b时，随着滚子与局部剥落故障边缘之间距离的减小，滚子与滚道之间的接触宽度从73 μm 减小为37 μm，呈减小趋势。仿真结果显示，随着局部剥落故障宽度的增加，滚子与滚道之间的接触变形量与接触应力也随之增加，而滚子与滚道之间的接触宽度则随之减小。
图4 滚道接触位置的变形、应力和接触宽度Fig.4 Deformation,stress and contact width of contact(a)接触变形；(b)接触应力；(c)接触宽度
图5所示为滚子通过局部剥落故障的过程中滚子与滚道上应力分布云图。从图5可知：滚子与滚道之间的最大接触应力所在的位置发生了变化，当接触点距离局部剥落故障边缘的距离较大时，最大接触应力位置位于接触点附近，在滚子趋近局部剥落故障边缘的过程中，最大接触应力位置逐渐移动到局部剥落故障边缘根部A，当滚子进入局部剥落故障区间时，最大接触应力位置到达局部剥落故障边缘顶部B。 2.2 局部剥落故障对接触刚度的影响
径向载荷与滚子-滚道之间接触变形的关系表达式为[21]
式中：n 为载荷-变形系数。为了简化拟合过程，对式(9)两边分别取对数，可将指数型公式转化为线性拟合问题，从而获得不同接触位置时滚子与滚道之间的接触刚度，其表达式为
取XP=-1.5b 到XP=1.5b 过程中，不同载荷作用下滚子与滚道的接触变形，利用数值拟合的方式计算滚子与滚道之间的接触刚度。当滚子径向载荷从10 N 增大到
25 N 时，XP从-1.5b 增大到1.5b过程中不同局部剥落故障宽度下对应的滚子与内圈滚道之间的接触刚度，如图6所示。从图6可知：当局部剥落故障宽度为b时，滚子与滚道之间的最小接触刚度为1.18×108 N/m1/0.723 1，最大接触刚度为4.98×108 N/m1/0.895 7。随着滚子与局部剥落故障中心位置距离的减小，滚子与滚道之间的接触刚度先减小后增大，且滚子位于局部剥落故障区间时，滚子与滚道之间的接触刚度小于滚子与正常滚道之间的接触刚度。 3 局部剥落故障过渡区对轴承振动特性的影响 3.1 圆柱滚子轴承的动力学方程
图5 滚子与滚道接触处应力分布Fig.5 Stress distribution between roller and raceway(a)接触点远离故障区间；(b)接触点趋近故障区间；(c)接触点位于故障区间
图6 接触刚度变化规律Fig.6 Contact stiffness between roller and raceway 为研究局部剥落故障对轴承振动特性的影响，将轴承滚子考虑为非线性接触弹簧，引入文献[22]中的简化圆柱滚子轴承2 自由度系统的动力学方程，其表达式为
式中：min为内圈与转轴总质量；c 为系统阻尼；K(θ)为滚子与内外圈滚道的综合接触刚度，其数值与图6中的一致；Fx和Fy分别为沿x 和y 方向的径向载荷分量；δi为第i个滚子与内外圈的接触变形。
任意角度θi下的滚子与内外圈滚道的综合接触变形为
式中：Hd为滚子通过轴承滚道局部剥落故障位置的时变位移激励。时变位移激励Hd通常用半正弦函数表示[22]，考虑到轴承局部剥落故障过渡区对时变位移激励的影响，将其表示为
式中：Δθ=wt/rd；wt为过渡区的宽度；rd为局部剥落故障所在的套圈的滚道半径。
3.2 局部剥落故障对轴承振动响应的影响
运用定步长的四阶龙格库塔方法对方程(11)进行求解。时间步长Δt=5 μs，方程中轴承内圈与转轴的总质量min=1.2 kg，系统阻尼系数设定为350 N∙s/m[22-23]，转轴速度为4 000 r/min，施加在转轴上的载荷分别为Fx=400 N 和Fy=0 N。系统初始位移为x=1 μm 和y=1 μm。系统初始速度为0 m/s。圆柱滚子轴承N306的参数如表2所示。
取外圈滚道上局部剥落故障的宽度为0.5b 和b，分别计算圆柱滚子轴承的振动响应。其中，局部剥落故障的深度都为b，过渡区宽度取局部剥落故障宽度的1/2，缺陷相对于X 轴正方向的初始角位置θd0=0°，计算结果如图7所示。从图7可知：局部剥落故障的宽度为0.5b 时，包含局部故障过渡区的轴承X向最大加速度为35.5 m/s2，而不含局部故障过渡区的轴承X向最大加速度为20.64 m/s2；局部剥落故障的宽度为b时，包含局部故障过渡区的轴承X 向最大加速度为78.43 m/s2，而不含局部故障过渡区的轴承X 向最大加速度为36.35 m/s2。对比结果显示，考虑轴承局部剥落故障过渡区时轴承振动加速度大于不考虑轴承局部剥落故障过渡区的轴承振动加速度。
图7 不同剥落故障宽度对应的轴承X向加速度对比Fig.7 Comparisons of accelerations inX direction of bearing with different spalling fault widths(a)局部剥落故障宽度为0.5b时轴承振动响应；(b)局部剥落故障宽度为b时轴承振动响应
图8 不同剥落故障宽度对应的轴承X方向位移的FFT变换对比Fig.8
Comparisons of FFT inX direction of bearing with different spalling fault widths(a)故障宽度为0.5b时X方向位移的FFT变换；(b)故障宽度为b时X方向
位移的FFT变换
外圈滚道上局部剥落故障的宽度为0.5b 和b时，轴承X 方向位移的频谱图如图8所示。从图8可知：2种故障宽度下都能获得相同的外圈通过频率及其倍频成分。2种故障宽度下获得的外圈通过频率均为314.7 Hz，与采用文献[15]的理论计算方法获得的结果一致，验证了仿真结果的正确性。 4 结论
1)当滚道存在局部剥落故障时，随着滚子与局部剥落故障中心位置距离的减小，滚子与滚道之间的接触应力增大，接触宽度减小，接触变形先增大后减小。 2)在滚子与局部剥落故障中心位置距离的减小过程中，滚子与滚道之间的最大接触应力位置由故障边缘根部移动到故障边缘顶部。
3)随着局部剥落故障宽度的增加，滚子与滚道之间的接触宽度减小，接触应力与接触变形量增大。
4) 滚子在通过局部剥落故障位置的过程中，滚子与局部剥落故障区域的接触刚度小于正常接触刚度。
5)在圆柱滚子轴承的动力学系统中，考虑轴承局部剥落故障过渡区时对应的轴承振动加速度大于不考虑局部故障过渡区的轴承振动加速度。
【相关文献】
[1]RUBINI R, MENEGHETTI U.Application of the envelope and wavelet transform analyses for the diagnosis of incipient faults in ball bearings[J].Mechanical Systems and Signal Processing,2001,15(2):287-302.
[2]杨国安.滚动轴承故障诊断实用技术[M].北京:中国石化出版社,2012:10-11.YANG Guoan.Practical technology for fault diagnosis of rolling bearings[M].Beijing:Sinopec Press,2012:10-11.
[3]李兴林,张仰平,曹茂来，等.滚动轴承故障监测诊断技术应用进展[J].工程与试验,2009,49(4):1-
5.LI Xinglin, ZHANG Yangping, CAO Maolai, et al.Development of fault detection and diagnosis technology of rolling bearing[J].Engineering&Test,2009,49(4):1-5.
[4]LIU Jing, SHAO Yimin, LIM T C.Vibration analysis of ball bearings with a localized defect applying piecewise response function[J].Mechanism and Machine Theory, 2012, 56(1):156-169.
[5]LIU Jing, SHAO Yimin.Dynamic modeling for rigid rotor bearing systems with a localized defect considering additional deformations at the sharp edges[J].Journal of Sound and Vibration,2017,398:84-102.
[6]BRANCH N A, ARAKERE N K, SVENDSEN V, et al.Critical stresses and strains at the spall edge of a case hardened bearing due to ball impact[J].International Journal of Fatigue,2013,47(2):268-278.
[7]PATIL M S,MATHEW J,RAJENDRAKUMAR P K,et al.A theoretical model to predict the effect of localized defect on vibrations associated with ball bearing[J].International Journal of Mechanical Sciences,2010,52(9):1193-1201.
[8]牛蔺楷,曹宏瑞,何正嘉.滚动轴承表面损伤建模与冲击力的定量计算[J].振动.测试与诊
断,2013,33(S1):5-8.NIU Linkai, CAO Hongrui, HE Zhengjia.Modeling of local surface defect in rolling bearing and quantitative calculation of impact force[J].Journal of Vibration, Measurement &Diagnosis,2013,33(S1):5-8.
[9]樊莉,谭南林,沈栋平.基于显式动力学的滚动轴承接触应力有限元分析[J].北京交通大学学报, 2006, 30(4): 109-112.FAN Li, TAN Nanlin, SHEN Dongping.FEA on stress field of rolling element bearing based on explicit dynamics[J].Journal of Beijing Jiaotong University,2006,30(4):109-112.
[10]张福星,郑源,汪清，等.基于ANSYS Workbench的深沟球轴承接触应力有限元分析[J].机械设计与制造,2012(10):222-224.ZHANG Fuxing, ZHENG Yuan,WANG Qing, et al.Contact stress FEM analysis of deep groove ball bearing based on ANSYS workbench[J].Machinery Design & Manufacture,2012(10):222-224.
[11]张雪,张钢,周凯峰，等.交叉滚子转盘轴承的有限元分析[J].轴承,2012(4):1-4.ZHANG Xue, ZHANG Gang, ZHOU Kaifeng, et al.FEA of cross-roller slewing bearings[J].Bearing,2012(4):1-4.
[12]伍济钢,杨曾增,蒋勉,等.滚动轴承表面损伤故障动力学建模方法研究[J].机械传
动,2016,40(2):27-33.WU Jigang,YANG Zengzeng, JIANG Mian, et al.Research of dynamics modeling method for rolling bearing with surface damage fault[J].Journal of Mechanical Transmission,2016,40(2):27-33.
[13]剡昌锋,康建雄,苑浩,等.考虑弹流润滑及滑动作用下滚动轴承系统局部缺陷位移激励动力学建模[J].振动与冲击,2018,37(5):56-64.YAN Changfeng, KANG Jianxiong, YUAN Hao, et al.Dynamic modeling for local defect displacement excitation in rolling bearing systems under elasto-hydrodynamic lubrication and slip[J].Journal of Vibration and
Shock,2018,37(5):56-64.
[14]SLACK T, SADEGHI F.Explicit finite element modeling of subsurface initiated spalling in rolling contacts[J].Tribology International,2010,43(9):1693-1702.
[15]HARRIS T A, KOTZALAS M N.Rolling bearing analysisessential concepts of bearing technology[M].5th ed.New York:Taylor and Francis,2007:113.
[16]王大力,孙立明,单服兵,等.ANSYS在求解轴承接触问题中的应用[J].轴承,2002(9):1-4.WANG Dali, SUN Liming, SHAN Fubing, et al.Application of ANSYS software in solution of contact problem of bearing[J].Bearing,2002(9):1-4.
[17]ANSYS theory reference.ANSYS release 11.0[M].Canonsburg:ANSYS,Inc.2007:367. [18]MOAVENI S.Finite element analysis: theory and application with ANSYS[J].Pearson Schweiz Ag, 2007, 9(1):136-41.
[19]姜育松,苏超.工程接触问题数值分析方法[J].水力发电,2010,36(4):75-78.JIANG Yusong, SU Chao.Numerical analysis methods for solving contact problems in engineering[J].Water Power,2010,36(4):75-78.
[20]PALMGREN A.Ball and roller bearing engineering[M].Philadelphia:Philadelphia SKF Industries Inc,1959:116.
[21]HERTZ H.On the contact of elastic solids[J].Journal of Reine and Angewandlte Mathematik,1881,92:156-171.
[22]刘静,邵毅敏,秦晓猛，等.基于非理想Hertz线接触特性的圆柱滚子轴承局部故障动力学建模[J].机械工程学报,2014,50(1):91-97.LIU Jing, SHAO Yimin, QIN Xiaomeng, et al.Dynamic modeling on localized defect of cylindrical roller bearing based on non-hertz line contact characteristics[J].Journal of Mechanical Engineering,2014,50(1):91-97.
[23]刘静,吴昊,邵毅敏,等.考虑内圈挡边表面波纹度的圆锥滚子轴承振动特征研究[J].机械工程学报,2018,54(8):26-34.LIU Jing, WU Hao, SHAO Yimin, et al.Investigation for vibrations of tapered roller bearing considering the surface waviness on the rib of the inner race[J].Journal of Mechanical Engineering,2018,54(8):26-34.
因篇幅问题不能全部显示，请点此查看更多更全内容