高速圆柱滚子轴承保持架运行稳定性分析

2019-05-27王自彬邓四二张文虎黄晓敏

振动与冲击 2019年9期

王自彬，邓四二,2，张文虎，黄晓敏

(1.河南科技大学机电工程学院，河南洛阳 471003；2.辽宁重大装备制造协同创新中心，辽宁大连 116024；3.重庆工业职业技术学院机械工程学院，重庆 401120)

高速圆柱滚子轴承以其优越的高速性能被广泛应用于航空发动机主轴。高速圆柱滚子轴承保持架动态性能及可靠性直接影响着轴承的工作性能，其运动不稳定性会引起轴承早期失效[1]。因此，有必要对影响高速圆柱滚子轴承保持架稳定性的关键因素进行分析，这对提高圆柱滚子轴承寿命及稳定性具有重要意义。

Ghaisas等[2]建立了六自由度圆柱滚子轴承动力学模型，通过保持架质心速度标准偏差比评价保持架稳定性，分析了保持架设计参数、内圈转速、滚子直径以及保持架不对称性对保持架运行稳定性的影响；Gupta[3-5]研究了保持架兜孔间隙、轴承动载荷、润滑剂、保持架不平衡量对保持架动态稳定性的影响；Meeks[6]建立了保持架的6自由度动力学模型，对保持架的设计参数进行了优化，改善了保持架的稳定性；Rivera[7]建立了简化的保持架分析模型，并指出滚动体与保持架之间的相互作用力是导致保持架运动不稳定的主要原因；Tada[8]发现减小保持架兜孔间隙和引导间隙可以有效降低由于保持架运动不稳定引起的噪声；Takafumi等[9]考虑了润滑油非牛顿流体特性和滚子与滚道间热效应，分析了转速、载荷对保持架打滑率的影响；Zhang等[10]建立了高速圆柱滚子轴承动力学模型，研究了四种航空润滑油对保持架质心轨迹的影响，并使用庞加莱图分析了保持架运行稳定性；邓四二等[11]建立了高速圆柱滚子轴承动力学微分方程，对保持架动力学特性进行了理论分析，研究了保持架引导方式、间隙比对保持架质心轨迹及其打滑率的影响；刘秀海等[12]建立了高速圆柱滚子轴承动力学模型，研究了轴承转速、载荷、游隙等对保持架质心运动轨迹的影响；姜维等[13]针对电机专用轴承进行了试验研究，通过电机电流和声音判断保持架稳定性；黄迪山[14]针对微型轴承实体保持架，研究其质心轨迹的检测信号时域分析特点，分析了载荷大小及加载方式对其质心轨迹的影响及时变特征，并讨论了盒维数、滤波、保持架变形对涡动描述的影响。

上述对于保持架运动轨迹稳定性的判定多依靠个人经验，存在较强的主观性。分形理论是现代数学的一个新的分支，它与动力系统的混沌理论交叉结合，从非线性复杂系统自身入手，能够描述不规则现象内部所蕴含的信息[15]。而轴承作为一个复杂的非线性系统，存在许多不规则的现象，因此分形理论在轴承相关的研究中得到了广泛地应用[16-19]。鉴于此，本文基于滚动轴承动力学理论，建立高速圆柱滚子轴承的非线性动力学模型，采用分形理论中的盒维数评价保持架质心轨迹混乱程度，分析保持架间隙比、转速、载荷、径向游隙以及滚子个数对保持架运行稳定性的影响，为高速圆柱滚子轴承保持架设计提供理论依据。

1 高速圆柱滚子轴承各元件相互作用力

本文以NU2214型圆柱滚子轴承为研究对象，外圈固定，内圈旋转，保持架采用外引导方式，保持架兜孔形状为矩形，并假设轴承各零件的质心与形心重合，轴承工作表面为理想表面。为了能较好地描述轴承各零件的运动，通常选取一个整体坐标系{O;X,Y,Z}以及若干个局部坐标系，如图1所示：局部坐标系包括保持架质心坐标系{oc;xc,yc,zc}，内圈质心坐标系{oi;xi,yi,zi}，滚子质心坐标系{or;xr,yr,zr}，保持架兜孔坐标系{op;xp,yp,zp}。

坐标变换的一般形式为

图1 圆柱滚子轴承坐标系系统Fig.1 Coordinate systems of cylindrical roller bearing {R}1=[T]({r}2+{d}2)

(1)

式中：{R}1为在整体坐标系中描述的向量；{r}2为在局部坐标系中描述的向量；{d}2为两坐标系间的平移量；[T]为两坐标系间的转动变换。

1.1 滚子与滚道之间作用力

本文采用“改进切片法”[20]计算滚子与滚道之间的相互作用力。将滚子在自身轴线方向上平均分割成N片，每个切片的宽度为w(w=l/N)。图2为采用切片法分析滚子与滚道之间相互作用问题的示意图，图中：θj为第j个滚子的倾斜角；α为内圈在径向力Fr和倾覆力矩M作用下的倾斜角；l为滚子全长；ls为滚子有效接触长度；Cjm为第j个滚子的第m个切片的凸度减少量。

图2 滚子与套圈间的接触模型示意图Fig.2 Schematic diagram of contact model between roller and rings

第j个滚子的第m个切片与内(外)滚道之间的法向接触力为

(2)

第j个滚子与内(外)滚道之间的法向接触力、拖动力分别为

(3)

(4)

第j个滚子与内(外)滚道之间由法向接触力、拖动力分布产生的附加力矩可分别表示为

(5)

(6)

1.2 滚子与保持架之间作用力

滚子与保持架横梁的接触处主要存在法向的接触力和切向的摩擦力，对于这两种力，本文仍采用“改进切片法”进行计算。滚子与保持架接触模型如图3所示，图中：rcj为第j个滚子的质心or与所在保持架兜孔的中心op沿切向方向上的距离；βj为第j个滚子的歪斜角。qcjm、Fcjm分别为第j个滚子的第m个切片与保持架横梁之间的接触力、摩擦力。

图3 滚子与保持架接触模型示意图Fig.3 Schematic diagram of contact model between roller and cage

第j个滚子的第m个切片与保持架横梁之间的接触变形量为

(7)

第j个滚子与保持架横梁之间的法向接触力、切向摩擦力可分别表示为

(8)

(9)

式中：μcjm为第j个滚子的第m个切片与保持架横梁之间的摩擦因数，可采用文献[21]中的边界润滑模型进行计算。

第j个滚子与保持架横梁之间由法向接触力分布产生的附加力矩为

(10)

1.3 保持架受力分析

除滚子对保持架的作用力之外，外圈引导面还会对保持架定心表面产生引导力，油气混合物会对保持架的表面和端面产生阻力。保持架的受力如图4所示，图中：ec为保持架中心的相对偏移量，ec在y、z轴上的分量分别为Δyc、Δzc；Ψc为整体坐标系与保持架坐标系在径向平面内的夹角。

图4 保持架受力示意图Fig.4 Schematic diagram of cage forces

(11)

(12)

式中：η0为常压下润滑油的动力黏度；Lc为保持架定心表面宽度；C1为保持架引导间隙；ε为保持架中心相对偏心量；u1为润滑油的拖动速度。

保持架与引导面的摩擦力矩为

(13)

式中：V1为引导面与定心面相对滑动速度；Rcd为保持架定心表面半径。

油气混合物对保持架的表面和端面产生的阻滞力矩TCDO、TCDS分别为

TCDO=τAcagercage

(14)

(15)

式中：τ为保持架引导面剪应力；Acage为保持架引导面面积；rcage为保持架非引导面半径；ρeff为油气混合物密度，CN的具体计算方法参见文献[10]。

2 高速圆柱滚子轴承动力学模型

2.1 滚子非线性动力学微分方程组

滚子的受力如图5所示，根据滚子的受力情况可得到第j个滚子的非线性动力学微分方程组为

(16)

图5 滚子受力图Fig.5 Schematic diagram of roller forces

2.2 保持架非线性动力学微分方程组

根据图4保持架的受力情况，建立保持架非线性动力学微分方程组

(17)

2.3 内圈非线性动力学微分方程组

内圈主要受到滚子对内圈法向接触力、摩擦力、附加力矩以及来自外界的力和力矩。根据内圈的受力情况，建立内圈非线性动力学微分方程组

(18)

3 盒维数及其计算方法

3.1 盒维数定义

盒维数是分形维数的一种，它反映了复杂形体占有空间的有效性，是复杂形体不规则性的度量[22]。设F是Rn上任意一非空有界子集，Nδ(F)表示最大直径为δ且能覆盖F的集合的最少个数，则根据文献[23]，F的上下盒维数分别定义为

(19)

(20)

如果上下维相等，则F的盒维数可定义为

(21)

3.2 盒维数计算方法

设离散信号y(i)⊂Y，Y是n维欧氏空间Rn上的闭集。用尽可能细的ε网格划分Rn，Nε是集合Y的网格计数。由于式(21)中的极限无法按定义求出，所以在计算时需采用近似的方法。以ε网格为基准，逐步放大到kε网格，其中k∈Z+。令Nkε为离散空间上的集合Y的网格计数，可由式(22)和式(23)计算得到。

min{yk(i-1)+1,…,yk(i-1)+k+1}|

(22)

式中：i=1,2,…,N/k；N为采样点数；k=1,2,…,M,M

网格计数Nkε为

Nkε=P(kε)/(kε)+1

(23)

式中：Nkε>1。

在lgkε～lgNkε图中确定线性较好的一段为无标度区，设无标度区的起点和终点分别为k1、k2，则：

lgNkε=algkε+bk1≤k≤k2

(24)

最后，用最小二乘法确定该直线的的斜率

(25)

盒维数dB为

(26)

3.3 使用MATLAB计算盒维数

本文在MATLAB中计算保持架质心轨迹盒维数，其具体方法如下：①将像素大小为2n×2n的保持架质心轨迹灰度图像导入到MATLAB，并对其进行二值化处理，生成黑白位图。该二值化图像在MATLAB中以矩阵的形式表示，矩阵的大小为2n×2n，且该矩阵中只含有0或1两种元素，0表示像素点为黑色，1表示像素点为白色；②以大小为k×k的网格对该二值化图像进行覆盖，统计包含0元素的网格个数，并记为Nk，通常取k=1,2,4,…,2k；③在对数坐标系中绘出点(lgk,lgNk)，k=1,2,4,…,2k，采用最小二乘法在无标度区内对其进行直线拟合，并输出拟合得到的直线的斜率，斜率的相反数即为该输入图像的盒维数。

为验证该计算方法的有效性与准确性，选取直线、Koch雪花、Sierpinski三角垫、Sierpinski方毯等几种典型的分形几何图形进行计算验证。计算验证结果如表1所示。由计算结果可知计算误差均控制在1%以内，可以证明该计算方法具有足够可靠的计算精度。

4 结果分析

本文采用预估-校正的GSTIFF变步长积分算法对高速圆柱滚子轴承的非线性动力学微分方程组进行求解，使用盒维数分析保持架间隙比、轴承转速、轴承径向载荷、轴承径向游隙以及滚子个数等参数对保持架运行稳定性的影响。已使用该方法对NU209、NU304、NU2214、NU2220等多个型号的圆柱滚子轴承进行了分析，发现分析结果呈现相同的变化趋势，因此本文仅以NU2214型圆柱滚子轴承为例进行分析，其余型号的轴承的分析结果不再赘述。本文研究发现：保持架质心轨迹越复杂和混乱，盒维数越大，保持架运行越不稳定，反之则越稳定。NU2214型圆柱滚子轴承的主要参数见表2。

表1 典型分形图形盒维数计算误差Tab.1 Calculation error in Box Dimension of typical fractal graphics

表2 轴承主要参数Tab.2 Major parameters of bearing

4.1 保持架间隙比对保持架稳定性的影响

当内圈转速n=20 000 r/min，径向载荷F=5 000 N，径向游隙μr=0.02 mm，滚子个数Z=17时，保持架质心轨迹随间隙比c(间隙比c定义为保持架兜孔间隙与引导间隙的比值)的变化如图6所示。根据图6不同间隙比时的保持架质心轨迹，计算出对应的盒维数如图7所示。

由图6、图7可知：当间隙比小于0.8时，随着间隙比的增大，保持架质心轨迹对应的盒维数先增大后减小，但保持架均处于较好的稳定涡动状态；且当c=0.8时，盒维数最小，保持架运行最稳定；当间隙比大于0.8时，保持架质心轨迹对应的盒维数随着间隙比的增大而增大，并且开始处于不稳定的涡动状态，这是因为随着间隙比的增大，引导间隙减小，引导面与保持架之间的不连续碰撞力增大，使保持架的运行处于不稳定状态。