CN110674585A - 高转速滚动轴承-转子系统局部损伤故障模拟方法 - Google Patents

Abstract

本发明涉及高转速滚动轴承‑转子系统局部损伤故障模拟方法，包括以下步骤：S1、建立带有局部损伤故障的轴承的模型；S2、获取轴承转子系统的初值条件；S3、将载荷带入高转速轴承拟静力学模型中；S4、根据内外接触角和接触变形的拟合公式，计算各个滚动体与轴承套圈之间的时变卷吸速度和接触载荷；S5、计算左右侧滚动体与轴承套圈之间考虑振动位移而导致的法向接触位移量；S6、建立轴承‑转子系统六自由度动力学模型；S7、使用变步长龙格库塔法，对轴承‑转子系统进行求解；S8、获取带有局部损伤故障的轴承转子系统振动响应时域信号以及频域信号。本发明能够更加准确的模拟带有轴承损伤故障的转子轴承系统实时运转情况，并输出相应的振动响应。

Description

高转速滚动轴承-转子系统局部损伤故障模拟方法

技术领域

本发明涉及高转速轴承-转子系统领域，更具体地说，涉及一种高转速滚动轴承-转子系统局部损伤故障模拟方法。

背景技术

滚动轴承-转子系统在旋转机械中广泛应用，高转速滚动轴承-转子系统的轴承转速大于等于10000r/min，它的运行状态对整个机械系统的精度、可靠性和寿命等性能有着重要影响。而滚动轴承损伤故障将会使系统正常运行过程中产生额外的异常振动和噪声，现阶段对故障的诊断大多处于使用传感器收集轴承故障信号，然后对信号进行分析，该方法实验成本高，周期长。因此，建立高转速轴承-转子系统损伤故障动力学模型，定量分析损伤类型及其相应振动响应，能够准确预测和识别滚动轴承的早期故障，成本低，时间短。

从现有文献中可以看出，轴承损伤故障对转子轴承系统振动的影响一直受到国内外学者的广泛关注。Liu等(Jing Liu,Yimin Shao,Dynamic modeling for rigid rotorbearing systems with a localized defect considering additional deformationsat the sharp edges[J],J.Sound Vib.398(2017)84-102)分析了在缺陷位置处额外的时变位移以及时变接触刚度，并将其用公式定量表示，带入到十二自由度动力学模型中，但是该模型不适合高速情况，没有考虑高速情况下滚动体的离心力等的复杂情况。剡昌锋等(剡昌锋，康建雄，苑浩，等.考虑弹流润滑及滑动作用下滚动轴承系统局部缺陷位移激励动力学建模[J].振动与冲击，2018，37(5):56-64)该模型综合考虑了入口油膜刚度，赫兹接触刚度以及赫兹接触区的油膜刚度，同时将损伤表示为一个定量公式，然而，该模型只分析了单个轴承，同时没有考虑高转速下的复杂现象。刘永强等(刘永强，王宝森，杨绍普.含外圈故障的高速列车轴承转子系统非线性动力学行为分析[J].机械工程学报，2018，54(8):17-25)建立了十二自由度模型，分析在损伤故障情况下，高速转子轴承系统的动力学响应，但是在该模型中没有考虑润滑作用以及摩擦力矩对系统的影响。曹宏瑞等(曹宏瑞，李亚敏，何正嘉，等.一种滚动轴承-转子系统损伤振动响应定量计算方法:中国，103712785[P].2017-01-04)在该模型中，建立轴承损伤故障特征参数，并将其带入到轴承平衡时内部结构的位移方程中，同时和滚动体平衡方程联立求解，得出带有损伤故障的轴承的刚度矩阵，然后建立转子轴承动力学方程并求解分析，该模型虽然同时考虑了轴承损伤故障和高速复杂受力情况，但是没有综合考虑油膜对转子轴承系统的影响。

从上述文献检索中可以看出，现有的转子轴承系统损伤故障动力学模型，没有综合考虑高速情况下滚动体的离心力等的复杂情况，以及弹流润滑油膜刚度、赫兹接触刚度导致的滚动体与轴承套圈之间的时变刚度，也没有分别考虑左右两侧轴承受力和损伤情况，这将无法准确模拟轴承出现损伤故障后所表现出的动力学现象；且出现损伤故障后，滚动体与轴承套圈在损伤接触区域，将会产生附加的随时间变化的额外变形量，同时不同类型的损伤，所产生的随时间变化的变形量是不同的，而上述模型没有综合考虑以上情况，使得转子轴承系统损伤故障模型不够完善。

发明内容

本发明要解决的技术问题在于，提供一种高转速滚动轴承-转子系统局部损伤故障模拟方法，能够更加准确的模拟带有轴承损伤故障的转子轴承系统实时运转情况，并输出相应的振动响应。

本发明解决其技术问题所采用的技术方案是：构造一种高转速滚动轴承-转子系统局部损伤故障模拟方法，包括以下步骤：

S1、建立带有局部损伤故障的轴承的模型；

S2、获取轴承转子系统的结构参数、工况参数、局部损伤参数以及初值条件；结构参数包括轴承参数、材料参数和转子的几何参数；工况参数包括转速和外力矢量；局部损伤参数包括损伤的长度、宽度、深度以及所在位置角，初始条件包括拟静力学模型和轴承转子系统动力学模型的求解初始值；轴承转子系统为外圈固定，内圈随转子旋转；

S3、将两个轴承所受的不同的载荷分别带入高转速轴承拟静力学模型中，得到在不同转速工况下，随时间变化的不同位置角处各个滚动体与轴承套圈之间的内外接触角和接触变形，将其拟合成以位置角为变量，以内外接触角和接触变形为函数值的曲线；

S4、根据内外接触角和接触变形的拟合公式，计算各个滚动体与轴承套圈之间的时变卷吸速度和接触载荷，并将其带入已有的弹流润滑油膜刚度阻尼计算公式中，计算油膜刚度，同时计算滚动体与轴承套圈之间的赫兹接触刚度，然后计算等效接触刚度；

S5、计算左右侧滚动体与轴承套圈之间考虑振动位移而导致的法向接触位移量，同时判断滚动体是否处于损伤故障处，并计算由于损伤而导致的随时间变化的额外位移量，然后计算滚动体与轴承套圈之间的作用力；

S6、建立轴承-转子系统六自由度动力学模型，在该模型中同时考虑润滑剂粘性所产生的摩擦力矩的影响以及滚动体与轴承套圈间的作用力的影响；

S7、使用变步长龙格库塔法，对轴承-转子系统进行求解；

S8、获取带有局部损伤故障的轴承转子系统振动响应时域信号以及频域信号。

上述方案中，所述步骤S3中的拟静力学模型通过以下方程表示：

其中，A_1j＝BDsinα₀+δ_a+R_iθcosθ_j，

其中，

γ′≠D/d_m；

其中，

其中，R_r＝A_2j-(f_o-0.5))Dcosα₀-f_iD cosα_ij+0.5D，R_i＝0.5d_m+(f_i-0.5)Dcosα_o；

式中，F_a为轴承所受轴向力，F_r为轴承所受径向力，K_oj为外圈的载荷位移常数，K_ij为内圈的载荷位移常数，M_gj为钢球陀螺力矩，δ_r为轴承在径向总变形量，δ_a为轴承在轴向的总变形量，θ为轴承的角变形，f_i为内圈的曲率系数，f_o为外圈的曲率系数，δ_ij第j个滚动体与内圈接触变形量，δ_oj为第j个滚动体与外圈接触变形量，B为总曲率，α_ij第j个滚动体与内圈接触角，α_oj为第j个滚动体与外圈接触角，α₀为初始接触角，λ_oj＝2，λ_ij＝0，D为滚动体直径，m为滚动体质量，w为内圈角速度，w_m为滚动体公转速度，w_R为滚动体自转速度，M′为轴承所受力矩，Δ_c为径向游隙减小量，β_j为球姿态角。

上述方案中，所述步骤S4中，根据内外接触角和接触变形，分别计算左右轴承卷吸速度u_ij，u_oj和左右轴承接触载荷w_ij，w_oj，并将其带入弹流润滑油膜刚度阻尼计算公式中；

滚动体与轴承套圈之间的卷吸速度用下式来表示：

u_ij＝0.5d_m[(1-γ)(w-ω_m)+γω_R] (3)

u_oj＝0.5d_m[(1+γ)ω_m+γω_R] (4)

其中，γ＝Dcosα/d_m

油膜刚度通过下式计算：

式中，E_p为等效弹性模量，

R_x为滚动体与轴承套圈之间的沿轴向的等效接触半径，M_(i，o)j＝W_(i，o)j(2U_f(i，o)j)^-3/4，L_(i，o)j＝1000(U_f(i，o)j)^1/4，U_f(i，o)j＝η_oru_(i，o)j/(E_pR_x(i，o))，

η_or参考粘度取0.1，u_(i，o)j为第j个滚动体与轴承套圈间的卷吸速度，w_(i，o)j为第i个滚动体与滚道间接触载荷，

G＝α₁E_p，κ_(i，o)＝1.0339(R_y(i，o)/R_x(i，o))^0.636，

η_o为润滑油粘度；

同时计算滚动体与轴承套圈之间的赫兹接触刚度，赫兹接触刚度公式为：

则滚动体与轴承套圈之间的等效刚度为：

滚动体与轴承套圈之间的时变油膜阻尼为：

其中：

式中，G＝α₁E_p，α₁为粘压系数，在恒定温度下是常数。

上述方案中，所述步骤S5中，根据滚动体的位置角是否处在局部损伤的起始位置角和终止位置角之间来判断滚动体是否在局部损伤处，然后计算由于损伤而导致的随时间变化的额外位移量，通过以下公式计算：

其中，

其中，

其中，

ΔT为凹坑故障的起始角和终止角之间的夹角，用下式来表示：

其中，

在该式中，

其中，H′₂，H′₃的值等于H′_sd，

H为损伤故障的高度；

式中，d为球的直径，L为损伤的长度，B为损伤的宽度，H′₁的值，θ₀为损伤所在位置角，mod()为求余函数，D_o，D_i分别为内外滚道直径。

上述方案中，所述步骤S5中，计算左右侧滚动体与轴承套圈之间考虑振动位移以及损伤处的额外变形量而导致的法向位移量，法向位移量用下式计算：

X_L＝X+l₁sinθ_y，Y_L＝Y-l₁sinθ_x，Z_L＝δ_epsinα_p-Z

X_R＝X-l₂sinθ_y，Y_R＝Y+l₂sinθ_x，Z_R＝δ_epsinα_p+Z

其中，θ_j为滚动体的位置角，用下式表示：

其中，ξ为判断滚动体是否在损伤区域，ξ＝1表示滚动体在损伤区域，ξ＝0表示滚动体在其他地方；

式中，A为内外沟曲率中心之间的距离，α_ij为实时变化的第j个滚动体与内滚道之间的接触角，d_m为轴承节圆直径，θ_j为第j个球的位置角，θ_x，θ_y，θ_z分别为转子绕X，Y，Z轴转动的角度，α_p，δ_ep分别为轴向预紧力导致的接触角和接触位移，l₁，l₂，l分别表示左侧轴承距离转子质心的距离，右侧轴承距离转子质心的距离，ω_c为保持架转速，n为滚动体个数。

上述方案中，滚动体与轴承套圈之间的时变接触力通过下式求得：

F_j(L，R)＝K_j(L，R)*δ_j(L，R)+C_j(L，R)*δ′_j(L，R) (12)

式中，C_j(L，R)为弹流润滑油膜的阻尼，δ′_j(L，R)为接触位移量的导数，L，R分别表示左右两个轴承。

上述方案中，所述步骤S6中轴承-转子系统六自由度动力学模型的动力学方程为：

式中，M为转子的质量，C为系统阻尼主要考虑滚动体赫兹接触变形以及外圈和轴承座的相互作用，I_x，I_y，I_z分别为转子在X，Y，Z方向上的转动惯量，载荷矢量作用点与左侧轴承间的距离，Mf为外载荷作用和润滑剂粘度而产生的摩擦力矩，用下列公式来表示：

M_f＝M_if+M_df (14)

其中，Mif为润滑剂粘性而产生的摩擦力矩，M_df为外载荷引起的摩擦力矩，用下面两个式子来表示：

M_if＝10^-7f₀(v₀N_i)^2/3(d_m)³ (15)

M_df＝f₁(0.9F_zcot(α_ij)-0.1F_x)d_m (16)

式中，f₀为一个与轴承类型相关的系数，v₀为润滑油运动粘度，f₁为与轴承结构和载荷相关的系数。

实施本发明的高转速滚动轴承-转子系统局部损伤故障模拟方法，具有以下有益效果：

本发明采用拟静力学方法求解接触变形和接触角的过程中，考虑了轴承在高速运转时的离心力，使高速轴承模型更加准确，接触角随时间变化，使得滚动体处于不同的位置，与轴承套圈之间的法向接触力也在实时变化，将滚动体与轴承套圈之间的接触等效为弹簧阻尼模型，同时考虑时变弹流润滑油膜刚度和赫兹接触刚度，并将多种类型的损伤模型耦合到动力学模型中，使得该模型能够准确的模拟多种类型损伤故障导致的动力学响应，为准确的预测和识别滚动轴承早期故障提供支持。

附图说明

下面将结合附图及实施例对本发明作进一步说明，附图中：

图1为本发明的流程示意图；

图2为本发明中损伤的简化模型以及其在轴承的位置，黑色点表示滚动体与损伤接触点，(a)损伤类型1，D＞＞min(L，B)，(b)损伤类型2，B＞L，D＞min(L，B)，(c)损伤类型3，B＝L，D＞min(L，B)，(d)损伤类型4，B＜L，D＞min(L，B)，(e)损伤类型5，D＞min(L，B)；

图3为带局部损伤的简化模型图；

图4为可以模拟的轴承转子系统图；

图5为10000r/min时的径向位移对比图；

图6为10000r/min时的径向位移频谱对比图；

图7为10000r/min时内圈有损伤的径向位移频谱图；

图8为15000r/min时有无损伤的径向位移频谱对比图。

具体实施方式

为了对本发明的技术特征、目的和效果有更加清楚的理解，现对照附图详细说明本发明的具体实施方式。

如图1所示，在本发明的高转速滚动轴承-转子系统局部损伤故障模拟方法，包括以下步骤：

1)根据实际情况中轴承的局部凹坑损伤的大小和形状以及滚动体滚过损伤与损伤边缘接触的具体情况，可以将轴承局部损伤类型转化为如图2(II)所示的几种情况，局部损伤故障在轴承上的具体位置如图2(I)所示；将带有局部损伤故障的轴承转化为如图3所示的数学模型。

2)获取轴承转子系统的结构参数、工况参数、局部损伤参数以及初值条件，结构参数主要包括轴承参数、材料参数、转子的几何参数，工况参数包括转速，外力矢量，局部损伤参数包括损伤的长度、宽度、深度以及所在位置角，初始条件包括拟静力学模型和轴承转子系统动力学模型的求解初始值；轴承转子系统为外圈固定，内圈随转子旋转，内圈转速N_i，沿三个方向上的外力矢量分别用F_x，F_y，F_z来表示，其作用在转轴上的相对位置如图4所示。

3)将两个轴承所受的不同的载荷分别带入高转速轴承拟静力学模型中，得到在不同转速工况下，随时间变化的不同位置角处各个滚动体与轴承套圈之间的内外接触角和接触变形，将其拟合成以位置角为变量，以内外接触角和接触变形为函数值的曲线。

考虑高速离心力的拟静力学模型可以用下列式子来表示：

其中，A_1j＝BDsinα₀+δ_a+R_iθcosθ_j,

其中，

γ′＝D/d_m；

其中，

其中，R_r＝A_2j-(f_o-0.5)Dcosα₀-f_iDcosα_ij+0.5D，R_i＝0.5d_m+(f_i-0.5)Dcosα_o；

使用Newton-Raphson对上述非线性方程组联立求解，直到满足精度，可以得到时变接触角和接触变形。

式中，δ_r，δ_a，θ分别为轴承在径向轴向的总变形量，以及角变形，f_i，f_o为内外圈的曲率系数，δ_ij，δ_ij为第j个滚动体与内外圈接触变形量，该值实时变化，B为总曲率，α_ij，α_oj为随时间和转速变化的第j个滚动体与内外滚道之间的接触角，α₀为初始接触角，λ_oj＝2，λ_ij＝0，D为滚动体直径，m为滚动体质量，ω为内圈转速，w_m为滚动体公转速度，w_R为滚动体自转速度。

4)根据内外接触角和接触变形的拟合公式，计算各个滚动体与轴承套圈之间的时变卷吸速度和接触载荷，并将其带入已有的弹流润滑油膜刚度阻尼计算公式中，计算油膜刚度，同时计算滚动体与轴承套圈之间的赫兹接触刚度，然后计算等效接触刚度；

滚动体与轴承套圈之间的卷吸速度在不同位置处是不相同的，因此是时刻变化的，可以用下式来表示：

u_ij＝0.5d_m[(1-γ)(w-ω_m)+γω_R] (3)

u_oj＝0.5d_m[(1+γ)ω_m+γω_R] (4)

其中，γ＝Dcosα/d_m

由于卷吸速度和接触载荷随着时间变化，同时转速不同，其值也不相同，因此弹流润滑油膜的刚度也是时刻变化的，油膜刚度可以用下式计算：

式中，E_p为等效弹性模量，

R_x为滚动体与轴承套圈之间的沿轴向的等效接触半径，M_(i，o)j＝W_(i，o)j(2U_f(i，o)j)^-3/4，L_(i，o)j＝1000(U_f(i，o)j)1/4，U_f(i，o)j＝η_oru_(i，o)j/(E_pR_x(i，o))，

η_or参考粘度取0.1，u_(i，o)j为第j个滚动体与轴承套圈间的卷吸速度，w_(i，o)j为第i个滚动体与滚道间接触载荷，

G＝α₁E_p，κ_(i，o)＝1.0339(R_y(i，o)/R_x(i，o))^0.636，

η_o为润滑油粘度。

同时计算滚动体与轴承套圈之间的赫兹接触刚度，赫兹接触刚度公式为：

则滚动体与轴承套圈之间的等效刚度为：

滚动体与轴承套圈之间的时变油膜阻尼为：

其中：

式中，G＝α₁E_p，α₁为粘压系数，在恒定温度下是一个常数。

5)计算左右侧滚动体与轴承套圈之间考虑振动位移而导致的法向接触位移量，同时判断滚动体是否处于损伤故障处，并计算由于损伤而导致的随时间变化的额外位移量，然后计算滚动体与轴承套圈之间的作用力；

首先，根据滚动体的位置角是否处在局部损伤的起始位置角和终止位置角之间来判断滚动体是否在局部损伤处，然后计算由于损伤而导致的随时间变化的额外位移量，可以用下述公式计算：

其中，

其中，

其中，

ΔT为凹坑故障的起始角和终止角之间的夹角，可以用下式来表示：

其中，

在该式子中，

其中，H′₂，H′₃的值等于H′_sd，

H为损伤故障的高度。

式中，d为球的直径，L为损伤的长度，B为损伤的宽度，H′₁的值，θ₀为损伤所在位置角，mod()为求余函数，D_o，D_i分别为内外滚道直径。

计算左右侧滚动体与轴承套圈之间考虑振动位移以及损伤处的额外变形量而导致的法向位移量，法向位移量用下式计算：

X_L＝X+l₁sinθ_y，Y_L＝Y-l₁sinθ_x，Z_L＝δ_epsinα_p-Z

X_R＝X-l₂sinθ_y，Y_R＝Y+l₂sinθ_x，Z_R＝δ_epsinα_p+Z

其中，θ_j为滚动体的位置角，用下式表示：

其中，ξ为判断滚动体是否在损伤区域，ξ＝1表示滚动体在损伤区域，ξ＝0表示滚动体在其他地方。

式中，A为内外沟曲率中心之间的距离，α_ij为实时变化的第j个滚动体与内滚道之间的接触角，d_m为轴承节圆直径，θ_j为第j个球的位置角，θ_x，θ_y，θ_z分别为转子绕X，Y，Z轴转动的角度，α_p，δ_ep分别为轴向预紧力导致的接触角和接触位移，l₁，l₂，l分别表示左侧轴承距离转子质心的距离，右侧轴承距离转子质心的距离，ω_c为保持架转速，n为滚动体个数。

滚动体与轴承套圈之间的时变接触力可以用下式求得：

F_j(L，R)＝K_j(L，R)*δ_j(L，R)+C_j(L，R)*δ′_j(L，R) (12)

式中，C_j(L，R)为弹流润滑油膜的阻尼，δ′_j(L，R)为接触位移量的导数，L，R分别表示左右两个轴承。

6)建立轴承-转子系统六自由度动力学模型，在该模型中同时考虑润滑剂粘性所产生的摩擦力矩的影响以及滚动体与轴承套圈间的作用力的影响，动力学方程为：

式中，M为转子的质量，C为系统阻尼主要考虑滚动体赫兹接触变形以及外圈和轴承座的相互作用，I_x，I_y，I_z分别为转子在X，Y，Z方向上的转动惯量，载荷矢量作用点与左侧轴承间的距离，M_f为外载荷作用和润滑剂粘度而产生的摩擦力矩，可以用下列公式来表示：

M_f＝M_if+M_df (14)

其中，M_if为润滑剂粘性而产生的摩擦力矩，M_df为外载荷引起的摩擦力矩，可以用下面两个式子来表示：

M_if＝10^-7f₀(v₀N_i)^2/3(d_m)³ (15)

M_df＝f₁(0.9F_zcot(α_ij)-0.1F_x)d_m (16)

式中，f₀为一个与轴承类型相关的系数，v_o为润滑油运动粘度，f₁为与轴承结构和载荷相关的系数。

7)使用变步长龙格库塔法，对轴承-转子系统进行求解，水平位移初值为10^-6，速度、角位移、角速度初值为0，绕Z方向的角速度初值为轴承旋转角速度。

8)得到带有局部凹坑损伤故障的轴承转子系统振动响应时域信号以及频域信号。

为了便于理解，避免遗漏，本实施例涉及到的所有公式参数如下：

X_1j，X_2j为中间过程辅助量，无具体定义；

δ_ij第j个滚动体与内圈接触变形量，为待求未知参量；δ_oj为第j个滚动体与外圈接触变形量，为待求未知参量；

K_oj、K_ij：载荷位移常数；量纲为1的接触位移；∑_ρi、∑_ρo曲率和

M_gj：钢球陀螺力矩；

J：钢球惯性力矩；

w为内圈角速度；

w_mj为第j个滚动体公转速度；

w_Rj为第j个滚动体自转速度；

α_ij第j个滚动体与内圈接触角，为待求未知参量；α_oj为第j个滚动体与外圈接触角，为待求未知参量；

β_j：球姿态角，公式计算；

θ_j：滚动体位置角，ω_c保持架角速度，w内圈即转子的角速度，

m：滚动体质量；n钢球个数；

d_m：节圆半径；

D为滚动体直径；

Δ_c为安装后径向游隙减小量；

B为总曲率；

δ_a为轴承在轴向的总变形量，为待求未知参量；θ为轴承角变形，为待求未知参量；F_a为轴承所受轴向力；F_r为轴承所受径向力；M′为轴承所受力矩；

r_i为轴承内滚道曲率半径；r_o为轴承外滚道曲率半径；

FF_(i，o)第一类完全椭圆积分，EE_(i，o)第二类完全椭圆积分，R_x(i，o)为滚动体与轴承内外套圈之间的沿轴向的等效接触半径；R_y(i，o)为滚动体与轴承内外套圈之间的沿径向的等效接触半径；

u_ij，u_oj：卷吸速度；

E_p：等效样式模量，v₁，v₂为两个接触材料的泊松比；E₁，Ex为两个接触材料的弹性模量；

η_o为润滑油粘度，η_or参考粘度，α₁为润滑油的黏压系数；

k_(i，o)j：赫兹接触刚度；k_S(i，o)j：油膜刚度

w_(i，o)j：接触载荷；

A：滚道沟曲率中心之间的距离；

F_p为轴向预紧力；α_p：预紧力产生的接触角；δ_ep：预紧力产生的轴向位移；

K_j(L，R)等效接触刚度；δ_jL、δ_jR法向位移量；C_j(L，R)为弹流润滑油膜的阻尼；δ′_j(L，R)为法向位移量的导数；

一个点表示一阶导数，两个点表示二阶导数；

M为转子的质量，C为系统阻尼主要考虑滚动体赫兹接触变形以及外圈和轴承座的相互作用；

f₀为一个与轴承类型相关的系数，通过查表可以得到；

v₀为润滑油运动粘度，查润滑油表可以得到；

N_i为轴转速；

f₁为与轴承结构和载荷相关的系数。

下面结合一个实例对本发明进一步做详细说明，本实例并不用于限制本发明。其步骤为：

1)获取轴承转子系统的结构参数、工况参数及初值条件，获取损伤的几何尺寸和位置参数，本次实例模拟的是图2(II)中的第二类损伤，建立的轴承转子系统的模型如图4所示，轴承转子系统的所有参数如表1所示：

表1轴承转子系统参数

2)根据步骤1)中的参数使用拟静力学方法得到时变内外接触角和接触变形的计算公式。

3)根据步骤2)中，内外接触角和接触变形，计算左右轴承卷吸速度u_ij，u_oj和接触载荷w_ij，w_oj，计算滚动体与轴承套圈之间的油膜刚度和赫兹接触刚度，再计算等效刚度以及阻尼。

4)判断滚动体是否位于损伤区域，计算由于损伤导致的额外的位移量，计算左右侧滚动体与轴承套圈总的位移量，根据步骤3)中得到的等效刚度、法向位移量，计算滚动体与轴承套圈之间的接触力。

5)计算摩擦力矩，以及步骤4)中计算所得的接触力，带入到动力学方程组(13)中。

6)用变步长龙格库塔法求解轴承5)中的动力学方程组(13)，并输出动力学响应。

图5和图6分别给出了转速为10000r/min时，左侧有损伤、右侧有损伤、两侧都有且损伤位于重载区域、两侧都有且位于轻载区域的以及无损伤的径向位移时域图和径向位移频谱图，从频谱图中可以看出滚动体通过频率符合理论计算的通过频率；左侧受载较大，损伤位于左侧轴承时幅值最大；两侧轴承都有损伤时，幅值较小，是因为左右两侧振动相互抵消一部分；无损伤时，幅值最小，符合相应理论。图7给出了损伤位于内圈时的频谱图，可以看出其通过频率与理论计算值相同。图8给出了转速为15000r/min时有无损伤的径向位移频谱图，在转速较高时仍然符合理论计算。

上面结合附图对本发明的实施例进行了描述，但是本发明并不局限于上述的具体实施方式，上述的具体实施方式仅仅是示意性的，而不是限制性的，本领域的普通技术人员在本发明的启示下，在不脱离本发明宗旨和权利要求所保护的范围情况下，还可做出很多形式，这些均属于本发明的保护之内。

Claims (7)

1.一种高转速滚动轴承-转子系统局部损伤故障模拟方法，其特征在于，包括以下步骤：

S1、建立带有局部损伤故障的轴承的模型；

S2、获取轴承转子系统的结构参数、工况参数、局部损伤参数以及初值条件；结构参数包括轴承参数、材料参数和转子的几何参数；工况参数包括转速和外力矢量；局部损伤参数包括损伤的长度、宽度、深度以及所在位置角，初始条件包括拟静力学模型和轴承转子系统动力学模型的求解初始值；轴承转子系统为外圈固定，内圈随转子旋转；

S3、将两个轴承所受的不同的载荷分别带入高转速轴承拟静力学模型中，得到在不同转速工况下，随时间变化的不同位置角处各个滚动体与轴承套圈之间的内外接触角和接触变形，将其拟合成以位置角为变量，以内外接触角和接触变形为函数值的曲线；

S4、根据内外接触角和接触变形的拟合公式，计算各个滚动体与轴承套圈之间的时变卷吸速度和接触载荷，并将其带入已有的弹流润滑油膜刚度阻尼计算公式中，计算油膜刚度，同时计算滚动体与轴承套圈之间的赫兹接触刚度，然后计算等效接触刚度；

S5、计算左右侧滚动体与轴承套圈之间考虑振动位移而导致的法向接触位移量，同时判断滚动体是否处于损伤故障处，并计算由于损伤而导致的随时间变化的额外位移量，然后计算滚动体与轴承套圈之间的作用力；

S6、建立轴承-转子系统六自由度动力学模型，在该模型中同时考虑润滑剂粘性所产生的摩擦力矩的影响以及滚动体与轴承套圈间的作用力的影响；

S7、使用变步长龙格库塔法，对轴承-转子系统进行求解；

S8、获取带有局部损伤故障的轴承转子系统振动响应时域信号以及频域信号。

2.根据权利要求1所述的高转速滚动轴承-转子系统局部损伤故障模拟方法，其特征在于，所述步骤S3中的拟静力学模型通过以下方程表示：

其中，A_1j＝BDsinα₀+δ_a+R_iθcosθ_j，

其中，γ′＝D/d_m；

其中，

其中，R_r＝A_2j-(f_o-0.5)Dcosα₀-f_iDcosα_ij+0.5D，R_i＝0.5d_m+(f_i-0.5)Dcosα_o；

式中，F_a为轴承所受轴向力，F_r为轴承所受径向力，K_oj为外圈的载荷位移常数，K_ij为内圈的载荷位移常数，M_gj为钢球陀螺力矩，δ_r为轴承在径向总变形量，δ_a为轴承在轴向的总变形量，θ为轴承的角变形，f_i为内圈的曲率系数，f_o为外圈的曲率系数，δ_ij第j个滚动体与内圈接触变形量，δ_oj为第j个滚动体与外圈接触变形量，B为总曲率，α_ij第j个滚动体与内圈接触角，α_oj为第j个滚动体与外圈接触角，α₀为初始接触角，λ_oj＝2，λ_ij＝0，D为滚动体直径，m为滚动体质量，w为内圈角速度，w_m为滚动体公转速度，w_R为滚动体自转速度，M′为轴承所受力矩，Δ_c为径向游隙减小量，β_j为球姿态角。

3.根据权利要求2所述的高转速滚动轴承-转子系统局部损伤故障模拟方法，其特征在于，所述步骤S4中，根据内外接触角和接触变形，分别计算左右轴承卷吸速度u_ij，u_oj和左右轴承接触载荷w_ij，w_oj，并将其带入弹流润滑油膜刚度阻尼计算公式中；

滚动体与轴承套圈之间的卷吸速度用下式来表示：

u_ij＝0.5d_m[(1-γ)(w-ω_m)+γω_R] (3)

u_oj＝0.5d_m[(1+γ)ω_m+γω_R] (4)

其中，γ＝Dcosα/d_m

油膜刚度通过下式计算：

式中，E_p为等效弹性模量，

R_x为滚动体与轴承套圈之间的沿轴向的等效接触半径，M_(i，o)j＝W_(i，o)j(2U_f(i，o)j)^-3/4，L_(i，o)j＝1000(U_f(i，o)j)^1/4，U_f(i，o)j＝η_oru_(i，o)j/(E_pR_x(i，o))，

η_or参考粘度取0.1，u_(i，o)j为第j个滚动体与轴承套圈间的卷吸速度，w_(i，o)j为第i个滚动体与滚道间接触载荷，

G＝α₁E_p，κ_(i，o)＝1.0339(R_y(i，o)/R_x(i，o))^0.636，

η_o为润滑油粘度；

同时计算滚动体与轴承套圈之间的赫兹接触刚度，赫兹接触刚度公式为：

则滚动体与轴承套圈之间的等效刚度为：

滚动体与轴承套圈之间的时变油膜阻尼为：

其中：

式中，G＝α₁E_p，α₁为粘压系数，在恒定温度下是常数。

4.根据权利要求3所述的高转速滚动轴承-转子系统局部损伤故障模拟方法，其特征在于，所述步骤S5中，根据滚动体的位置角是否处在局部损伤的起始位置角和终止位置角之间来判断滚动体是否在局部损伤处，然后计算由于损伤而导致的随时间变化的额外位移量，通过以下公式计算：

其中，

其中，

其中，

ΔT为凹坑故障的起始角和终止角之间的夹角，用下式来表示：

其中，

在该式中，

其中，H′₂，H′₃的值等于H′_sd，

H为损伤故障的高度；

式中，d为球的直径，L为损伤的长度，B为损伤的宽度，H′₁的值，θ₀为损伤所在位置角，mod()为求余函数，D_o，D_i分别为内外滚道直径。

5.根据权利要求4所述的高转速滚动轴承-转子系统局部损伤故障模拟方法，其特征在于，所述步骤S5中，计算左右侧滚动体与轴承套圈之间考虑振动位移以及损伤处的额外变形量而导致的法向位移量，法向位移量用下式计算：

X_L＝X+l₁sinθ_y，Y_L＝Y-l₁sinθ_x，Z_L＝δ_epsinα_p-Z

X_R＝X-l₂sinθ_y，Y_R＝Y+l₂sinθ_x，Z_R＝δ_epsinα_p+Z

其中，θ_j为滚动体的位置角，用下式表示：

其中，ξ为判断滚动体是否在损伤区域，ξ＝1表示滚动体在损伤区域，ξ＝0表示滚动体在其他地方；

式中，A为内外沟曲率中心之间的距离，α_ij为实时变化的第j个滚动体与内滚道之间的接触角，d_m为轴承节圆直径，θ_j为第j个球的位置角，θ_x，θ_y，θ_z分别为转子绕X，Y，Z轴转动的角度，α_p，δ_ep分别为轴向预紧力导致的接触角和接触位移，l₁，l₂，l分别表示左侧轴承距离转子质心的距离，右侧轴承距离转子质心的距离，ω_c为保持架转速，n为滚动体个数。

6.根据权利要求5所述的高转速滚动轴承-转子系统局部损伤故障模拟方法，其特征在于，滚动体与轴承套圈之间的时变接触力通过下式求得：

F_j(L，R)＝K_j(L，R)*δ_j(L，R)+C_j(L，R)*δ′_j(L，R) (12)

式中，C_j(L，R)为弹流润滑油膜的阻尼，δ′_j(L，R)为接触位移量的导数，L，R分别表示左右两个轴承。

7.根据权利要求6所述的高转速滚动轴承-转子系统局部损伤故障模拟方法，其特征在于，所述步骤S6中轴承-转子系统六自由度动力学模型的动力学方程为：

式中，M为转子的质量，C为系统阻尼主要考虑滚动体赫兹接触变形以及外圈和轴承座的相互作用，I_x，I_y，I_z分别为转子在X，Y，Z方向上的转动惯量，载荷矢量作用点与左侧轴承间的距离，M_f为外载荷作用和润滑剂粘度而产生的摩擦力矩，用下列公式来表示：

M_f＝M_if+M_df (14)

其中，M_if为润滑剂粘性而产生的摩擦力矩，M_df为外载荷引起的摩擦力矩，用下面两个式子来表示：

M_if＝10^-7f₀(v₀N_i)^2/3(d_m)³ (15)

M_df＝f₁(0.9F_zcot(α_ij)-0.1F_x)d_m (16)

式中，f₀为一个与轴承类型相关的系数，v₀为润滑油运动粘度，f₁为与轴承结构和载荷相关的系数。

Images