保持架过梁断裂前后球轴承的动力学行为分析

2024-01-03房晓辉徐自力黄金平

火箭推进 2023年6期

房晓辉,徐自力,王珺,黄金平

(1.西安交通大学机械结构强度与振动国家重点实验室,陕西西安 710049;2.西安航天动力研究所液体火箭发动机技术重点实验室,陕西西安 710100)

0 引言

轴承是保证液体火箭发动机涡轮泵正常工作的关键部件之一,涡轮泵工作在高压力、高转速等极端恶劣工况下[1],并自始至终伴随着振动。涡轮泵球轴承内部的保持架与内、外圈及滚珠之间存在频繁的摩擦和碰撞,碰撞引起的振动使保持架过梁产生较大的剪应力,进而引发疲劳断裂[2],如果发生故障隐患排除不及时,则导致轴承的工作精度下降,轴承的使用寿命和可靠性降低,甚至影响涡轮泵系统的正常运行。因此,进行高速球轴承故障的动力学行为研究,及时有效地预测球轴承保持架断裂故障对防止事故发生具有重要意义。

球轴承特性研究主要分为实验法、力学模型解析法及有限元法3种。在球轴承特性的实验研究方面,Tauqir等针对某发动机轴承疲劳失效引发的飞机坠毁事故,进行了故障调研和实验模拟,指出事故原因在于保持架疲劳断裂[3]。兴长喜等实测了WP6发动机轴承保持架疲劳断裂引起的多起故障,研究了疲劳断裂机理,并提出了预防措施[4]。徐锐等针对某发动机地面试车过程异常问题进行调查,发现原因在于轴承保持架断裂[5]。保持架断裂引发的事故屡屡发生,在实际工程中,保持架断裂是高速运转轴承失效的重要原因之一。

轴承力学模型解析法通过建立滚动体、保持架、内圈的受力平衡方程组,求解非线性方程组,得到轴承的受力、振动响应等动力学行为。Liu等采用时变位移激励函数,应用赫兹接触理论计算轴承接触刚度,研究了缺陷形状和缺陷偏置角对轴承振动响应的影响[6]。Jafari等将轴承之间的接触描述为弹簧—阻尼系统,假设轴承保持架为刚性,研究了轴承外圈剥落缺陷对轴承振动响应的影响[7]。胡爱军等建立了轴承5自由度动力学模型,假设滚珠做纯滚动运动,研究了外圈故障数量、故障间隔及载荷分布对轴承故障特征的影响[8]。王凯等考虑了复合故障、时变位移和滚动体滑动等因素,建立了4自由度复合故障深沟球轴承动力学模型,研究了径向载荷作用下复合故障振动机理[9]。Wen等考虑了离心力和陀螺力矩作用,建立了轴承保持架3自由度模型,并采用了几何约束关系中的局部缺陷效应代替位移激励函数,研究了缺陷尺寸和缺陷位置对轴承动力学行为的影响[10]。在力学模型解析法中,大部分模型采用单一的位移激励或者力激励函数来代替实际故障产生的冲击激励,难以精确反映滚动体经过故障区的实际接触情况和真实的振动情况。同时,力学模型解析法往往采用了刚体假设,或对轴承部件的自由度进行了简化,与实际情况存在一定的偏差。

随着轴承研究及计算机技术的进一步发展,有限元法为轴承复杂模型的求解及动力学行为的精确计算提供了帮助。许多学者应用有限元法进行轴承动力学特性的计算并验证了方法的有效性[11-12]。Singh等运用LS-DYNA建立了外圈滚道缺陷的轴承二维有限元模型,研究了滚珠经过缺陷前后的接触及振动响应的变化[13]。Liu等建立了模拟轴承滚动体和滚道间接触的二维有限元模型,研究了剥落边缘轮廓对轴承塑性变形和接触力的影响[14-15]。马辉等运用ANSYS LS-DYNA分析了轴承故障区域的平滑程度对轴承振动特性的影响[16]。丁东升等模拟外圈滚道多个位置处的缺陷故障,研究了不同缺陷位置处的动力学响应规律[17]。He等运用ANSYS LS-DYNA建立了滚子剥落缺陷的轴承二维有限元模型,研究了滚子剥落缺陷对轴承接触力和滑动行为的影响[18]。有限元法能够真实反映故障轴承的受力和实际运行情况,计算故障轴承的动力学行为更加准确。

目前的研究中,针对内圈、外圈、滚动体等故障轴承的研究较为全面,关于保持架断裂后球轴承的动力学行为的相关研究较少;同时,轴承有限元模型的应用多局限在二维平面,仅能考虑单一径向载荷的影响,而轴承在三维联合载荷(如径向和轴向共同承载)作用下的情况则研究较少。本文针对高速状态下保持架过梁断裂引起的球轴承动力学行为变化的问题,分析了保持架断裂后的受力和运动学关系,在ANSYS LS-DYNA中建立正常轴承和保持架断裂轴承的三维有限元仿真模型,对比了断裂前后保持架和滚珠的瞬时作用力、打滑率、磨损、最大PV值等动力学行为。

1 保持架断裂后球轴承动力学分析

断裂发生在两个滚珠之间的保持架过梁,不考虑断裂后引入的杂质等异物,断裂位置不存在尖角等能够产生应力集中的部位。保持架断裂位置如图1所示。

图1 保持架断裂位置Fig.1 Position of cage fracture

保持架过梁断裂后,球轴承动力学特征变化具体表现为滚珠失去了保持架的均匀隔开作用,工作过程中断裂区的滚珠可能发生相互接触;由于断裂故障引起的保持架部分质量缺失,高速运行状态产生的不平衡作用力不容忽视,保持架不平衡质量me引起的自身不平衡作用力Fme为

(1)

式中:df为保持架内径;dg为保持架外径;nc为保持架实际转动速度。

为准确描述球轴承各部件的受力和运动学关系,以滚珠和保持架接触为例,建立了简化模型,如图2所示。图中,Oi-xiyizi为轴承的惯性坐标系,Ob-xbybzb为滚珠方位坐标系,Oc-xcyczc为固定于保持架坐标系,Of-xfyfzf为保持架兜孔坐标系。

图2 滚珠与保持架的相互作用和位置关系Fig.2 Interaction and position relationship between ball and cage

保持架兜孔中心到滚珠质心的向量可以表示为

(2)

滚珠和保持架兜孔之间的间隙可以表示为

(3)

由赫兹点接触理论得到滚珠对保持架的法向压力Fcbn,即

(4)

式中Ko为滚珠和保持架之间的赫兹接触刚度系数。

滚珠对保持架的摩擦力Fcbf按照Coulomb摩擦定律计算,即

Fcbf=ucFcbn

(5)

其中

uc=FD+(Fs-FD)e-λvrel

(6)

式中:uc为接触摩擦因数;FD为动摩擦因数;Fs为静摩擦因数;λ为指数衰减系数;vrel为接触表面的相对速度。

Mcb=|rpf|Fcb

(7)

式中rpf为保持架兜孔中心到接触点p的向量。

采用上述方法,可以得到引导外圈对保持架的法向压力Fcrn、切向摩擦力Fcrf和作用力矩Mcr。

断裂保持架上的作用力主要有滚珠和外圈对其作用力及自身的不平衡力。断裂保持架受力如图3所示。

图3 断裂保持架受力Fig.3 Forces on the fracture cage

断裂保持架的平移运动微分方程为

(8)

式中:Fcrx、Fcry和Fcrz为外圈对保持架作用力的分量;Fmex、Fmey和Fmez为不平衡力的分量;Fcbxk、Fcbyk和Fcbzk为第k个滚珠对保持架作用力的分量;Nb为滚珠的个数。

断裂保持架的旋转运动微分方程为

(9)

式中:Mcrx、Mcry和Mcrz为外圈对保持架作用力矩的分量;Mcbxk、Mcbyk和Mcbzk为第k个滚珠对保持架作用力矩的分量。

第4步，输入新的生产设备采购量情况样本或者已有生产设备采购量在不同情况下的样本，到支持向量机模型中进行生产设备采购量预测。

基于滚珠和保持架的受力和运动学关系,采用显式算法进行运动微分方程的求解,可以实现轴承瞬时作用力、打滑率等动力学行为的分析。

2 球轴承动力学计算的有限元法

角接触球轴承能够承受联合载荷作用,为考虑轴向载荷对轴承部件动力学行为的影响,采用三维单元来建立轴承有限元模型。内、外圈、滚珠和保持架均采用SOLID164六面体单元,内圈最内侧的单元采用SHELL163面单元,轴承所有部件都为弹性体,可以考虑轴承的弹性变形。采用映射网格划分,得到规则的六面体网格。

接触类型设置为自动面—面接触,定义滚珠为接触面,内、外圈滚道及保持架为目标面。LS-DYNA采用罚函数法计算接触力,在时间点ti+1=ti+Δt判断是否发生接触,如果存在接触则将接触力离散到接触单元节点上,并通过调整接触刚度罚因子fs控制穿透量。设置全局刚度系数IQH、体积黏度系数Q1及线性体积黏度系数Q2等接触参数避免沙漏现象,提高计算精度。为模拟轴承实际工作状况,轴向载荷加载到内圈侧面,径向载荷加载到内圈下表面,转速加载到内圈上。

在ANSYS LS-DYNA中进行有限元前处理、网格划分及接触、载荷、仿真时间等参数设置,并通过LS-DYNA求解器求解,LS-DYNA的离散化结构运动方程为

(10)

其时间积分采用显式差分算法求解,与隐式算法相比,显式算法用更少的时间增量求解高度非线性问题,不需要迭代和收敛准则,也不需要求解切线刚度矩阵,可以节省大量时间。轴承部件的瞬时作用力、打滑率等动力学行为的计算流程如图4所示。

图4 动力学行为的计算流程Fig.4 The computational process of dynamic behavior

3 结果分析

3.1 研究对象及数值验证

本文研究对象为某涡轮泵球轴承。轴承内、外圈及滚珠材料为9Cr18MoV,其密度为7 750 kg/m3,弹性模量为206 GPa,泊松比为0.3;保持架材料为聚四氟乙烯,密度为2 344.5 kg/m3,弹性模量为1.65 GPa,泊松比为0.4。各部件的动摩擦因数FD为0.05,静摩擦因数Fs为0.1。轴承几何参数见表1。

表1 轴承几何参数Tab.1 Geometric parameters of bearing

对于无故障轴承,在外圈固定约束时,可采用理论公式由内圈转速得到保持架内径处的速度vC及滚珠与内圈接触点的速度vB,即

(11)

(12)

式中Ne为内圈转速。

对轴承内圈施加5 000 r/min的转速时,保持架和滚珠的理论速度求解值分别为9.96 m/s和21.20 m/s。分别采用3种不同网格密度的有限元模型进行计算,并与理论结果进行对比,如表2所示。

表2 轴承部件转速求解结果Tab.2 Solution results of bearing component speed

由表2可知,采用网格数为377 208的有限元模型,已能达到较高精度。综合考虑计算效率和数值解的准确性,本文后续计算均采用该网格进行计算,无故障球轴承有限元网格如图5所示。

图5 无故障球轴承有限元网格Fig.5 Finite element mesh of fault free ball bearing

在无故障球轴承有限元网格尺寸的基础上,将两个滚珠之间的保持架过梁进行断裂处理来模拟断裂故障,仿真计算不考虑断裂后引入的杂质等异物,保持架过梁断裂后的有限元网格如图6所示。

图6 保持架断裂后的有限元网格Fig.6 Finite element mesh after cage fracture

在保持架断裂后,通过本文建立的轴承保持架断裂模型,计算了滚珠在运行至不同周向方位角时滚珠与外圈的接触力,并与Harris-Jones轴承模型[19]的理论计算结果进行对比,如图7所示。

图7 不同位置处滚珠与外圈的接触力Fig.7 Contact force between ball and outer ring at different positions

由图7可知,本文计算结果与Harris-Jones理论模型计算结果趋势一致,最大误差为1.32%,表明了保持架断裂模型在模拟保持架断裂后轴承动力学行为的正确性。

3.2 滚珠的动力学行为

保持架的作用是其过梁将滚珠沿圆周方向均匀隔开,在轴向载荷为1.5 kN、径向载荷为500 N、转速为12 000 r/min的工况下,保持架过梁断裂后,断裂位置处相邻的两个滚珠之间的作用力如图8所示。

图8 断裂区滚珠之间的作用力Fig.8 Force between ball and ball in fracture area

碰撞瞬间,轴承部件的速度云图如图9所示。

图9 碰撞瞬间轴承部件的速度云图Fig.9 Velocity nephogram of bearing components at the moment of impact

由图8和图9可知,当过梁断裂后,滚珠失去了保持架的均匀隔开作用,会发生两个滚珠的相互碰撞,瞬时碰撞力高达750 N,碰撞点的速度为36.98 m/s。

对比了无故障轴承和保持架断裂轴承在断裂区滚珠的质心速度,如图10所示。

图10 保持架断裂区和无故障轴承的滚珠质心速度Fig.10 Ball centroid velocity of cage fracture area and normal bearing

由图10可知,无故障轴承的滚珠质心速度变化小,波动幅度仅为0.24 m/s。保持架断裂后,断裂区的滚珠由于彼此碰撞,破坏了原来的运行规律,滚珠运动具有较大随机性,滚珠质心速度变化大且无规律,波动幅度达4.17 m/s。

滚动体的打滑率为

(13)

保持架断裂区和无故障轴承的滚珠打滑率如图11所示。

图11 保持架断裂区和无故障轴承的滚珠打滑率Fig.11 Ballslip rate of cage fracture area and normal bearing

由图11可知,无故障轴承的滚珠打滑率稳定在2%左右,而断裂区滚珠打滑率在0.00%～14.06%之间。这是因为断裂区的滚珠发生相互碰撞后,高速状态下由于离心力的作用,滚珠与外圈的接触作用增强,与内圈的接触作用减弱,内圈对滚珠的拖动力远小于滚珠之间的相互作用力,滚珠原有的运动规律遭到破坏,滚珠更容易与内圈发生打滑。

绘制断裂区滚珠表面某点的运行轨迹如图12所示。图中A、B、C分别为滚珠运行轨迹上不同时刻的空间位置。在AB段,由于滚珠未发生碰撞,滚珠在不断地自转,运行轨迹正常。在B点,滚珠之间发生碰撞,在BC段,滚珠出现多次打滑,运行轨迹发生变化,滚珠运行轨迹的变化也会影响轴承其他部件的运行规律,使得轴承的性能下降。

图12 断裂区滚珠表面的运行轨迹Fig.12 Trajectory of ball bearing surface in fracture area

3.3 保持架的动力学行为

在轴向载荷为1.5 kN、径向载荷为500 N、转速为12 000 r/min的工况下,对该球轴承在保持架断裂前后保持架与外圈之间的作用力进行计算,绘制的作用力变化曲线如图13所示。由图13可知,无故障保持架与外圈之间的平均作用力为83.22 N,瞬时作用力最大值为775.73 N;断裂保持架与外圈之间的平均作用力为227.72 N,增大了173.64%,瞬时作用力最大值为1.006 97 kN,增大了29.81%。保持架断裂后,保持架与外圈之间的作用力远大于断裂前的作用力。这是因为断裂区滚珠发生碰撞后,内圈对其拖动力减弱,导致保持架与滚珠之间产生较大的速度差,滚珠以较大的力将保持架推向外圈引导边,并且由于保持架断裂后滚珠存在较大的运动间隙,速度差持续时间比较长,保持架与外圈接触更频繁。

图13 保持架与外圈之间的作用力Fig.13 Force between cage and outer ring

保持架打滑是导致轴承早期失效的主要因素之一[20]。保持架与外圈间产生剧烈滑动可能导致接触表面滑蹭损伤,过大的摩擦热量可能导致表面材料软化脱落形成麻坑、接触表面几何形貌发生变化、产生裂纹、噪音等[21]。断裂前后保持架的打滑率如图14所示。

图14 断裂保持架和无故障保持架的打滑率Fig.14 Slip rate of fracture cage and normal cage

由图14可知,保持架断裂后打滑率明显增大,最大打滑率由5.32%上升至12.56%。这是因为保持架转动依靠滚珠的推动作用,当保持架断裂后,保持架与外圈之间的作用力增大,使得保持架被紧紧推向外圈引导边,当保持架受到的推动作用小于受到的各种阻力,保持架实际转速便会降低,发生打滑。

轴承的PV值等于接触应力与滑动速度的乘积,表征了轴承的发热量。在接触区上的不同点具有不同的PV值,由赫兹接触理论,选取接触区的中心位置,计算外圈作用下保持架的最大PV值。绘制保持架断裂前后的最大PV值如图15所示,其中F为接触应力,v为滑动速度。

图15 断裂保持架和无故障保持架的最大PV值Fig.15 Maximum PV value of fracture cage and normal cage

由图15可知,保持架断裂后最大PV值显著增大。这是因为保持架断裂后,外圈对保持架的接触作用增强,保持架与外圈频繁接触,保持架打滑率增大,更容易发生接触表面的相对滑动,导致PV值增大。

(14)

式中:Fc(t)为保持架导向面的接触应力;v(t)为保持架导向面的滑动速度;KW为磨损系数;H为材料的硬度。

在轴向载荷为1.5 kN,径向载荷为500 N,内圈转速分别设置为12 000、13 000、14 000、15 000 r/min的工况下,断裂保持架与无故障保持架导向面的磨损比γ随转速的变化如图16所示。由图16可知,保持架断裂后,保持架导向面磨损加剧,在转速12 000 r/min工况下,断裂保持架导向面的磨损为断裂前的2.21倍;随着转速增大,不平衡力的作用增强,断裂保持架的磨损倍率随转速增加而上升,更易引发工作表面变粗糙、摩擦温度升高、振动及噪声加剧等结果。