一种无保持架球轴承滚动体打滑的数值分析方法

1.本发明涉及无保持架球轴承技术领域，尤其是一种无保持架球轴承滚动体打滑的数值分析方法。


背景技术：

2.无保持架球轴承是磁悬浮轴承系统中重要的零件之一，是磁悬浮轴承中的辅助轴承，起着保护主轴的作用；无保持架球轴承内圈与主轴连接，无保持架球轴承外圈与轴承座连接，内外圈之间有滚动体，滚动体在无保持架球轴承不工作时主要堆积在轴承的下半圈，滚动体之间相互接触；在紧急情况下(磁悬浮轴承断电、故障)，主轴会下落到无保持架球轴承的内圈并带动其旋转，无保持架球轴承内外圈之间的滚动体会平衡主轴下落的冲击力和主轴旋转产生的力矩，不会使主轴直接掉落在轴承座上，对主轴起到保护作用；无保持架球轴承在工作时滚动体受力时刻发生变化，非常容易出现打滑现象；其中无保持架球轴承出现打滑是指滚动体和轴承的内圈或外圈在圆周方向上发生相对位移；由于不存在保持架分离滚动体，无保持架球轴承的滚动体在打滑时不仅会导致滚动体与内外圈的润滑油膜失效，使滚动体与内外圈直接接触摩擦，也会导致相邻滚动体发生非连续碰撞与高速摩擦，这都会加剧轴承发热和滚动体磨损进而影响轴承的使用寿命；因此，对无保持架球轴承滚动体的打滑分析有着重要意义。目前尚无法通过试验手段获取准确的滚动体滑移速度，进而无法判定滚动体的打滑现象。
3.基于上述问题，提出一种无保持架球轴承滚动体打滑的数值分析方法。


技术实现要素：

4.本发明要解决的技术问题是：为了解决目前尚无法通过试验手段获取准确的滚动体滑移速度，进而无法判定滚动体的打滑现象的问题，本发明提供了一种无保持架球轴承滚动体打滑的数值分析方法，以无保持架球轴承为对象，通过建立无保持架球轴承滚动体的动力学模型并求解，在不同轴承载荷、不同轴承转速，对无保持架球轴承滚动体出现打滑进行分析，进而判定滚动体的打滑现象。
5.本发明解决其技术问题所采用的技术方案是：
6.一种无保持架球轴承滚动体打滑的数值分析方法，包括如下步骤：
7.1)建立滚动体的公转角加速度和自转角加速度的动力学模型，具体为：
[0008][0009]
式中，ψj为滚动体绕轴承回转中心角度，为公转角加速度，由滚动体绕轴承回转
中心角度ψj的二阶导数表示；θj为滚动体绕质心旋转角度，为自转角加速度，由滚动体绕质心旋转角度θj的二阶导数表示；
[0010]fij
为滚动体与内圈的摩擦力；f
oj
为滚动体与外圈的摩擦力；f
bf(j，j 1)
为滚动体与后一个滚动体的摩擦力；f
bf(j-1，j)
为滚动体与前一个滚动体的摩擦力；为滚动体半径；jr为滚动体惯性矩；为轴承内圈沟底半径；为外圈沟底半径；gr为滚动体重力；为轴承节圆半径；f
impact(j，j 1)
为滚动体与后一个滚动体的碰撞力； f
impact(j-1，j)
为滚动体与前一个滚动体的碰撞力；为为为
[0011]
2)求解滚动体的公转角加速度和自转角加速度的动力学模型，求得公转角速度和自转角速度；
[0012]
3)构建滚动体与内圈的的滑移速度公式，即：
[0013][0014]
式中：v
(ij)rel
为滚动体相对轴承内圈的滑移速度；ωi为轴承内圈转速；ω
cj
为滚动体公转角速度；ω
bj
为滚动体自转角速度；
[0015]
4)构建滚动体与外圈的滑移速度公式，即：
[0016][0017]
式中，v
(oj)rel
为滚动体相对于轴承外圈的滑移速度；
[0018]
5)将公转角速度和自转角速度带入公式(3)和(4)中，得到滚动体与轴承内、外圈的滑移速度；
[0019]
6)根据滚动体与轴承内、外圈的滑移速度随滚动体位置角的变化曲线评价分析滚动体的打滑情况。
[0020]
具体的，滚动体与内圈的摩擦力f
ij
具体为：
[0021][0022]
式中，ni为滚动体在和内、外圈接触点的法线方向上受到的内圈的压力。
[0023]
具体的，所述滚动体在和内、外圈接触点的法线方向上受到的内圈的压力ni具体为：
[0024]
[0025]
式中，fr为轴承所受径向载荷，z为滚动体的个数，α为轴承接触角，ε为载荷分布系数，ψ为滚动体位置角，第j个滚动体的位置角为ψj。
[0026]
具体的，所述滚动体与外圈的摩擦力f
oj
具体为：
[0027][0028]
式中，no为外圈对滚动体的反作用力。
[0029]
具体的，所述外圈对滚动体的反作用力no具体为：
[0030][0031]
式中，f
cj
为滚动体离心力，gr为滚动体的自身重力，ψj为为第j个滚动体位置角，f
bf(j,j 1)
＝-μf
impact(j,j 1)
，f
bf(j-1,j)
＝-μf
impact(j-1,j)
，为为为
[0032]
所述滚动体离心力f
cj
为：
[0033][0034]
式中ω
cj
为滚动体公转角速度，m为滚动体的质量。
[0035]
具体的，所述滚动体与后一个滚动体的碰撞力f
impact(j,j 1)
为：
[0036][0037]
式中，δ为滚动体碰撞后的变形量，k为刚度系数，c为阻尼(为刚度的1％)；
[0038]
式中，η为两滚动体综合弹性常数，e1，e2为滚动体材料的弹性模量(在这里滚动体材料一样则e1＝e2)，为与接触点主曲率相关的系数(经验值为1)，ρ为接触点主曲率(在这里为定值4/4.7625)。
[0039]
本发明的有益效果是：本发明提供的一种无保持架球轴承滚动体打滑的数值分析方法，可以得到不同轴承转速、载荷下滚动体与轴承内外圈的滑移速度，如图5、图6所示，当滑移速度不等于0时即表示滚动体出现打滑现象，滚动体处于(2π/3,3π/2)的角位置时，滚动体相对于轴承内、外圈的滑移速度为0，滚动体处于纯滚动状态即没有发生打滑现象；当滚动体处于(0,2π/3)(3π/2,0) 的角位置时，滚动体相对于轴承内、外圈的滑移速度不恒为
0，此时滚动体产生打滑现象；滑移速度小于0时表示滚动体相对于内外圈，往运动的相反方向滑移，滑移速度大于0时表示滚动体相对于轴承内、外圈，往运动的相同方向滑移(由滚动体重力或滚动体碰撞作用产生)。
附图说明
[0040]
图1：滚动体的公转角加速度和自转角加速度的动力学模型求解分析图；
[0041]
图2：滚动体受力与运动示意图；
[0042]
图3：滚动体公转角速度；
[0043]
图4：滚动体自转角速度；
[0044]
图5：极坐标系下不同转速下滚动体与轴承内外圈的滑移速度；
[0045]
图6：极坐标系下不同载荷下滚动体与轴承内外圈的滑移速度。
具体实施方式
[0046]
现在对本发明作进一步详细的说明。
[0047]
首先，分析滚动体的受力和运动，如图2所示，首先分析单个滚动体(第j 个球)的受力，在工作中轴承内圈会受到来自主轴所施加的载荷，受载的内圈会向下对滚动体施加力，因此滚动体在和内、外圈接触点的法线方向上会受到内圈的压力(ni)和外圈对滚动体的反作用力(no)。由于滚动体在轴承中相对于内圈和外圈是运动的，并且滚动体与内外圈是相互接触的，因此滚动体在与内外圈接触点的切线方向会受到来自于内、外圈的摩擦力(fo和fi)。单个滚动体有可能会受到来自相邻滚动体的碰撞力(f
impact
)，当相邻滚动体发生碰撞时，由于滚动体是运动的且在碰撞点的运动速度方向相反，会由此产生滚动体间的摩擦力(f
bf
)。由于滚动体绕质心公转，这会产生离心力(fc)。滚动体的自身重力为(gr)。
[0048]
再分析单个滚动体的运动，轴承外圈固定不动，内圈逆时针旋转驱动所有滚动体运动，滚动体的运动包括：滚动体绕轴承回转中心的公转(ωc)和绕滚动体质心的自转(ωb)。其公转角加速度可由滚动体绕轴承回转中心角度(ψj) 的二阶导数表示，其自转角加速度可由滚动体绕质心旋转角度(θj)的二阶导数表示。
[0049]
一种无保持架球轴承滚动体打滑的数值分析方法，如图1所示，包括如下步骤：
[0050]
1)建立滚动体的公转角加速度和自转角加速度的动力学模型，具体为：
[0051][0052]
式中，ψj为滚动体绕轴承回转中心角度，为公转角加速度，由滚动体绕轴承回转中心角度ψj的二阶导数表示；θj为滚动体绕质心旋转角度，为自转角加速度，由滚动体绕质心旋转角度θj的二阶导数表示；
[0053]fij
为滚动体与内圈的摩擦力；f
oj
为滚动体与外圈的摩擦力；f
bf(j-1，j)
为滚动体与
后一个滚动体的摩擦力；f
bf(j-1，j)
为滚动体与前一个滚动体的摩擦力；为滚动体半径；jr为滚动体惯性矩；为轴承内圈沟底半径；为外圈沟底半径；gr为滚动体重力；为轴承节圆半径；f
impact(j，j 1)
为滚动体与后一个滚动体的碰撞力； f
impact(j-1，j)
为滚动体与前一个滚动体的碰撞力；为为为
[0054]
2)求解滚动体的公转角加速度和自转角加速度的动力学模型，求得公转角速度和自转角速度；
[0055]
3)构建滚动体与内圈的的滑移速度公式，即：
[0056][0057]
式中：v
(ij)rel
为滚动体相对轴承内圈的滑移速度；ωi为轴承内圈转速；ω
cj
为滚动体公转角速度；ω
bj
为滚动体自转角速度；
[0058]
4)构建滚动体与外圈的滑移速度公式，即：
[0059][0060]
式中，v
(oj)rel
为滚动体相对于轴承外圈的滑移速度；
[0061]
5)将公转角速度和自转角速度带入公式(3)和(4)中，得到滚动体与轴承内、外圈的滑移速度；
[0062]
6)根据滚动体与轴承内、外圈的滑移速度随滚动体位置角的变化曲线评价分析滚动体的打滑情况。
[0063]
具体的，滚动体与内圈的摩擦力f
ij
具体为：
[0064][0065]
式中，ni为滚动体在和内、外圈接触点的法线方向上受到的内圈的压力。
[0066]
具体的，所述滚动体在和内、外圈接触点的法线方向上受到的内圈的压力ni具体为：
[0067][0068]
式中，fr为轴承所受径向载荷，在本具体实施例中取值为 1000n，2000n，3000n，4000n，z为滚动体的个数，α为轴承接触角，在本具体实施例中取值为0，ε为载荷分布系数，在本具体实施例中取值为1，ψ为滚动体位置角。
[0069]
具体的，所述滚动体与外圈的摩擦力f
oj
具体为：
[0070][0071]
式中，no为外圈对滚动体的反作用力。
[0072]
具体的，所述外圈对滚动体的反作用力no具体为：
[0073][0074]
式中，f
cj
为滚动体离心力，gr为滚动体的自身重力，ψj为第j个滚动体位置角，f
bf(j,j 1)
＝-μf
impact(j,j 1)
，f
bf(j-1,j)
＝-μf
impact(j-1,j)
，为为为
[0075]
滚动体离心力f
cj
为：
[0076][0077]
式中ω
cj
为滚动体公转角速度，m为滚动体的质量。
[0078]
具体的，滚动体与后一个滚动体的碰撞力f
impact(j,j 1)
为：
[0079][0080]
式中，δ为滚动体碰撞后的变形量，k为刚度系数，c为阻尼(为刚度的1％)；
[0081]
式中，η为两滚动体综合弹性常数，e1，e2为滚动体材料的弹性模量(在这里滚动体材料一样则e1＝e2)，为与接触点主曲率相关的系数(经验值为1)，ρ为接触点主曲率(在这里为定值4/4.7625)。
[0082]
在本具体实施例中，无保持架球轴承模型拟选取的参数值如表1-1所示：
[0083]
表1-1无保持架球轴承模型参数
[0084][0085]
滚动体的公转角加速度和自转角加速度方程程序代码为：
[0086]
function[dy]＝test804(t,y)
[0087]
m＝0.0035519；％滚动体质量
[0088]
g＝0.0348；％滚动体重量
[0089]
dw＝0.009525；％滚动体直径(m)
[0090]
r＝dw/2；滚动体半径(m)
[0091]
jr＝m*dw^2/10；％滚动体转动惯量
[0092]
dm＝0.04605；％节圆半径(m)
[0093]
di＝0.03648；％内滚道沟底半径(m)
[0094]
do＝0.05553；％外滚道沟底半径(m)
[0095]
rm＝dm/2；％节圆半径(m)
[0096]
ri＝di/2；％内滚道沟底半径(m)
[0097]
ro＝do/2；％外滚道沟底半径(m)
[0098]
dy＝zeros(60,1)；
[0099]
dy(1)＝y(2)；
[0100]
dy(2)＝(mui(y(2),y(4))*ni(y(1))*ri*deltavi(y(2),y(4))-muo(y(2),y(4))*deltavo (y(2),y(4))*(ni(y(1)) m*rm*y(2)*y(2)-g*cos(y(1))-0.002*contactball(y(1),y(5))*c os((y(1)-y(5))/2) 0.002*contactball(2*pi y(57),y(1))*cos((2*pi y(57)-y(1))/2) con tactball(y(1),y(5))*sin((y(1)-y(5))/2) contactball(2*pi y(57),y(1))*sin((2*pi y(57)
‑ꢀ
y(1))/2))*ro-g*sin(y(1) pi)*rm contactball(y(1),y(5))*cos((y(1)-y(5))/2)*(rm-r* sin((y(1)-y(5))/2))-contactball(2*pi y(57),y(1))*cos((2*pi y(57)-y(1))/2)*(rm-r*si n((2*pi y(57)-y(1))/2)))/m/rm/rm；％mui代表滚动体与内圈的摩擦系数；muo代表滚动体与外圈的摩擦系数；deltavo是判断滚动体相对于内外圈滑移的方向，和运动方向同向为正，反向为负；contactball表示滚动体间的
碰撞力；pi代表π。
[0101]
dy(3)＝y(4)；
[0102]
dy(4)＝(mui(y(2),y(4))*ni(y(1))*deltavi(y(2),y(4))-contactball(y(1),y(5))*0.00 2-contactball(2*pi y(53),y(1))*0.002)*dw/2/jr (muo(y(2),y(4))*deltavo(y(2),y(4)) *(ni(y(1)) m*rm*y(2)*y(2)-g*cos(y(1)) 0.002*contactball(y(1),y(5))*cos((y(1)-y (5))/2)-0.002*contactball(2*pi y(57),y(1))*cos((2*pi y(57)-y(1))/2) contactball(y( 1),y(5))*sin((y(1)-y(5))/2) contactball(2*pi y(57),y(1))*sin((2*pi y(57)-y(1))/2)))* r/jr；
[0103]
end
[0104]
程序2
[0105]
微分方程求解程序
[0106]
w_i＝100*pi；％轴承内圈转速
[0107]
dw＝0.009525；％滚动体直径
[0108]
r＝dw/2；滚动体半径(m)
[0109]
dm＝0.04605；％节圆半径(m)
[0110]
di＝0.03648；％内滚道沟底半径(m)
[0111]
do＝0.05553；％外滚道沟底半径(m)
[0112]
rm＝dm/2；％节圆半径(m)
[0113]
ri＝di/2；％内滚道沟底半径(m)
[0114]
ro＝do/2；％外滚道沟底半径(m)
[0115]
w_r＝r_i*r_o*w_i/(2*r*rm)；％滚动体自转角速度
[0116]
w_m＝r_i*w_i/(2*rm)；％滚动体公转角速度
[0117]
xspan＝[0,1]；％设置求解时间
[0118]
s＝2*asin(r/rm)；％相邻滚动体接触时的质心到回转中心的夹角
[0119]
y0＝[14*s,w_m,0,w_r,13*s,w_m,0,w_r,12*s,w_m,0,w_r,11*s,w_m,0,w_r,10*s,w _m,0,w_r,9*s,w_m,0,w_r,8*s,w_m,0,w_r,7*s,w_m,0,w_r,6*s,w_m,0,w_r,5*s,w_m,0, w_r,4*s,w_m,0,w_r,3*s,w_m,0,w_r,2*s,w_m,0,w_r,s,w_m,0,w_r,0,w_m,0,w_r]；
[0120]
％微分方程初始值的确定(数字代表滚动体编号从0开始)
[0121]
h＝1e-6；％设置积分步长(每间隔10-6
s计算一次)
[0122]
[t,y]＝adams4(@test804,xspan,y0,h)；％adams算法；test804表示在此处调用程序1。
[0123]
获得滚动体公转角速度如图2所示，滚动体自转角速度如图3所示。
[0124]
程序3
[0125]
不同轴承内圈转速下(1500r/min、2000r/min、3000r/min、4500r/min)，极坐标系下，滚动体滑移速度曲线随滚动体位置角变化的绘制程序
[0126]
w_i1＝50*pi；
[0127]
w_i2＝200/3*pi；
[0128]
w_i3＝100*pi；
[0129]
w_i4＝150*pi；
[0130]
％代表不同转速1500r/min、2000r/min、3000r/min、4500r/min
[0131]
figure('color',[1 1 1])
[0132]
g_1＝y1500(:,17)；
[0133]
x_1＝find(g_1》14.5*2*pi&g_1《15.5*2*pi)；％滚动体转动到第14.5-15.5圈的数据。
[0134]
g_1(x_1)；
[0135]
t＝mod(g_1(x_1),2*pi)；％将滚动体时间坐标系下的角度变化转换成以2π为周期的极坐标下的角度变化。
[0136]
a＝x_1(1:1)；
[0137]
b＝x_1(end:end)；
[0138]
i_1＝(r_m-r)*(w_i1-y1500(a:b,18))-r*(y1500(a:b,20))；％转速1500r/min下滚动体的滑移速度
[0139]
polarplot(t,i_1,'b','linewidth',1.2)
[0140]
hold on
[0141]
g_2＝y2000(:,17)/(2*pi)；
[0142]
x_2＝find(g_2》14.5&g_2《15.5)；
[0143]
g_2(x_2)；
[0144]
c＝x_2(1:1)；
[0145]
d＝x_2(end:end)；
[0146]
i_2＝(r_m-r)*(w_i2-y2000(c:d,18))-r*(y2000(c:d,20))；％转速2000r/min下滚动体的滑移速度
[0147]
polarplot(g_2(x_2)*2*pi,i_2,'r','linewidth',1.2)
[0148]
hold on
[0149]
g_3＝y3000(:,17)/(2*pi)；
[0150]
x_3＝find(g_3》14.5&g_3《15.5)；
[0151]
g_3(x_3)；
[0152]
e＝x_3(1:1)；
[0153]
f＝x_3(end:end)；
[0154]
i_3＝(r_m-r)*(w_i3-y3000(e:f,18))-r*(y3000(e:f,20))；％转速3000r/min下滚动体的滑移速度
[0155]
polarplot(g_3(x_3)*2*pi,i_3,'g','linewidth',1.2)
[0156]
hold on
[0157]
g_6＝y4500(:,17)/(2*pi)；
[0158]
x_6＝find(g_6》14.5&g_6《15.5)；
[0159]
g_6(x_6)；
[0160]
u＝x_6(1:1)；
[0161]
w＝x_6(end:end)；
[0162]
i_6＝(r_m-r)*(w_i4-y4500(u:w,18))-r*(y4500(u:w,20))；％转速4500r/min下滚动体的滑移速度
[0163]
polarplot(g_6(x_6)*2*pi,i_6,'k','linewidth',1.2)
[0164]
pax＝gca；
[0165]
pax.thetaaxisunits＝'radians'；
[0166]
pax.thetazerolocation＝'top'；
[0167]
pax.fontname＝'times new roman'；
[0168]
pax.fontsize＝12；
[0169]
f＝legend('1500r/min','2000r/min','3000r/min','4500r/min')；
[0170]
set(f,'fontname','times new roman','fontsize',12)
[0171]
legend('boxoff')
[0172]
rlim([-0.1 0.3])
[0173]
程序4
[0174]
％不同径向载荷下(1000n，2000n，3000n，4000n)，极坐标系下，滚动体滑移速度曲线随滚动体位置角变化的绘制程序：
[0175]
figure('color',[1 1 1])
[0176]
g_4＝y1000(:,17)/(2*pi)；
[0177]
x_4＝find(g_4》14&g_4《15)；
[0178]
g_4(x_4)；
[0179]
m＝x_4(1:1)；
[0180]
n＝x_4(end:end)；
[0181]
i_4＝(r_m-r)*(w_i-y1000(m:n,18))-r*(y1000(m:n,20))；％载荷1000n下滚动体的滑移速度
[0182]
polarplot(g_4(x_4)*2*pi,i_4,'b','linewidth',1)
[0183]
hold on
[0184]
g_5＝y2000(:,17)/(2*pi)；
[0185]
x_5＝find(g_5》14&g_5《15)；
[0186]
g_5(x_5)；
[0187]
o＝x_5(1:1)；
[0188]
p＝x_5(end:end)；
[0189]
i_5＝(r_m-r)*(w_i-y2000(o:p,18))-r*(y2000(o:p,20))；％载荷2000n下滚动体的滑移速度
[0190]
polarplot(g_5(x_5)*2*pi,i_5,'r','linewidth',1)
[0191]
hold on
[0192]
g_6＝y3000(:,17)/(2*pi)；
[0193]
x_6＝find(g_6》14&g_6《15)；
[0194]
g_6(x_6)；
[0195]
u＝x_6(1:1)；
[0196]
w＝x_6(end:end)；
[0197]
i_6＝(r_m-r)*(w_i-y3000(u:w,18))-r*(y3000(u:w,20))；％载荷3000n下滚动体的滑移速度
[0198]
polarplot(g_6(x_6)*2*pi,i_6,'g','linewidth',1)
[0199]
hold on
[0200]
g_7＝y4000(:,17)/(2*pi)；
[0201]
x_7＝find(g_7》14&g_7《15)；
[0202]
g_7(x_7)；
[0203]
u＝x_7(1:1)；
[0204]
w＝x_7(end:end)；
[0205]
i_7＝(r_m-r)*(w_i-y4000(u:w,18))-r*(y4000(u:w,20))；％载荷4000n下滚动体的滑移速度
[0206]
polarplot(g_7(x_7)*2*pi,i_7,'k','linewidth',1)
[0207]
pax＝gca；
[0208]
pax.thetaaxisunits＝'radians'；
[0209]
％pax.thetadir＝'clockwise'；
[0210]
pax.thetazerolocation＝'top'；
[0211]
pax.fontname＝'times new roman'；
[0212]
pax.fontsize＝12；
[0213]
f＝legend('1000n','2000n','3000n','4000n')；
[0214]
set(f,'fontname','times new roman','fontsize',12)
[0215]
legend('boxoff')
[0216]
rlim([-0.02 0.1])
[0217]
得到不同轴承转速、载荷下滚动体与轴承内外圈的滑移速度，如图5和图6 所示，当滑移速度不等于0时即表示滚动体出现打滑现象，滚动体处于(2π/3,3 π/2)的角位置时，滚动体相对于轴承内、外圈的滑移速度为0，滚动体处于纯滚动状态即没有发生打滑现象；当滚动体处于(0,2π/3)(3π/2,0)的角位置时，滚动体相对于轴承内、外圈的滑移速度不恒为0，此时滚动体产生打滑现象。滑移速度小于0时表示滚动体相对于内外圈，往运动的相反方向滑移，滑移速度大于0时表示滚动体相对于轴承内、外圈，往运动的相同方向滑移(由滚动体重力或滚动体碰撞作用产生)。
[0218]
以上述依据本发明的理想实施例为启示，通过上述的说明内容，相关工作人员完全可以在不偏离本项发明技术思想的范围内，进行多样的变更以及修改。本项发明的技术性范围并不局限于说明书上的内容，必须要根据权利要求范围来确定其技术性范围。