采用线接触的圆柱滚子轴承具有刚性好、摩擦力矩低和承载能力强等优点，因此被广泛应用于航空航天、高铁、风电、工业机器人等重要工业领域^[1-4]. 滚动轴承内部滚动体与内外滚道发生相对滑动会造成轴承打滑现象，不仅会增加摩擦力矩，加剧内部接触体与保持架的碰撞、磨损，还会迫使轴承运行温度升高并发生热变形，进而令轴承内部游隙减小和运行精度下降^[5-6]. 因此，准确解析圆柱滚子轴承的打滑特性，对于提升圆柱滚子轴承的稳定性和可靠性是极为重要的.

学者围绕滚动轴承打滑特性分析开展了大量的研究工作. Poplawski^[7]在Harris^[8-9] 拟静力学模型的基础上，忽略滚子歪斜特征，优化滚子轴承拟静力学模型，并分析接触区全油膜润滑下的保持架及滚子的打滑特性. Guo等^[10]进一步构建圆柱滚子轴承分析模型，详细讨论径向载荷、润滑剂黏度及转速等对保持架打滑行为的影响. 上述研究多基于拟静力学或拟动力学方法，针对轴承保持架及滚动体进行稳态分析. 然而，滚动轴承内部各元件间常常存在非连续接触及碰撞之类的作用，这使得上述模型无法实现对滚动轴承瞬态打滑行为的精准预测.

Laniado-Jácome等^[11]建立滚子轴承的有限元仿真模型，并分析滚子与内外套圈间的相对滑动问题. Kim^[12-13]利用简化的动力学模型分析径向游隙及几何偏差对轴承内部保持架打滑的影响. Gupta^[14-15]对Walters^[16]所建立的动力学模型进行改进，在考虑弹流润滑的条件下，建立圆柱滚子轴承的完全动力学模型，以实现对轴承内部滚子及保持架运动及受力状态的全面分析. Liu等^[17-18]建立弹流润滑条件下的滚动轴承打滑动力学模型，提出用于轴承打滑行为改善的结构参数优化设计准则. 贾永川等^[19-20]、李猛等^[21]分别建立启停阶段及非稳态工况下的圆柱滚子轴承动力学模型，分析径向载荷、工作温度对滚子及保持架打滑率的影响. 韩勤锴等^[22]研究径向力、弯矩之类的工况参数对轴承打滑特性的影响. Tu等^[23]建立圆柱滚子轴承的非线性动力学模型，研究滚子的局部打滑行为，分析径向载荷、主轴转速及保持架兜孔间隙对于滚子局部打滑情况的影响. 目前，在圆柱滚子轴承打滑的动力学分析方面，研究主要关注的是高速轻载之类工况下的轴承保持架宏观打滑现象，对于滚动体随其位置变化的局域打滑现象还有待深入研究^[24]. 保持架打滑对应于轴承整体打滑行为，是轴承内部各滚动体局域打滑行为的“平均化”表征，无法真实反映轴承内部滚动体的实际运行状态. 因此，须进一步针对不同载荷工况下的圆柱滚子轴承内部滚动体的局域打滑特性与机理开展研究.

本研究基于圆柱滚子轴承完全动力学模型，针对其内部滚子的局域打滑特性进行系统性研究，讨论载荷工况（径向载荷、转速）对圆柱滚子轴承内部滚子局域打滑行为的影响，并探究中低速中等载荷条件下圆柱滚子轴承内部滚子局部打滑的影响因素及变化规律，同时分析滚子打滑（局域打滑）与保持架打滑（整体打滑）特性间的差异与内在联系.

为了对轴承各组件的受力及运动进行描述，建立用于圆柱滚子轴承理论分析的多坐标系统，包括惯性坐标系及若干局部零件坐标系. 如图1所示，建立全局惯性坐标系$\left( {O,X,Y,Z} \right)$、滚子坐标系$\left( {{o_{\text{R}}},{x_{\text{R}}},{y_{\text{R}}},{z_{\text{R}}}} \right)$、套圈坐标系$\left( {{o_{\text{i}}},{x_{\text{i}}},{y_{\text{i}}},{z_{\text{i}}}} \right)$、保持架坐标系$\left( {{o_{\text{c}}},{x_{\text{c}}},{y_{\text{c}}},{z_{\text{c}}}} \right)$和接触区局域坐标系$\left( {{o_{\text{e}}},{x_{\text{e}}},{y_{\text{e}}},{z_{\text{e}}}} \right)$.

圆柱滚子轴承的惯性坐标系为固定坐标系，其坐标原点位于轴承质心，其空间位置不随时间及轴承零件的运动发生改变. 而轴承内各零件的局域坐标系则为运动坐标系，其空间位置随着零件的运动而持续发生变化，零件坐标系的坐标原点固定于各零件质心. 为了建立轴承的运动微分方程，在各零件坐标系下对各零件自身的运动和受力情况进行描述. 为了简化计算，在构建滚动轴承内部各组件的运动方程的过程中，首先在局部坐标系下对各零件自身的运动和受力情况进行描述，然后再通过坐标变换将上述关系统一转换到惯性系下进行求解计算. 通过三维坐标转换矩阵，可以将任意两坐标系下的位置坐标或空间向量进行坐标变换. 其中，惯性系与任意零件坐标系间的转换关系如下：

(1)$ \begin{split} {{\boldsymbol{T}}_{{\mathrm{ib}}}} &= {{{\boldsymbol{T}}_{{\mathrm{ib}}}}\left( {\eta ,\beta ,\gamma } \right)} = \\& \left[ {\begin{array}{*{20}{c}} {\cos \;\beta\; \cos\; \gamma }&\begin{gathered} \cos \;\eta\; \sin\; \gamma + \\ \sin \;\eta\; \sin\; \beta\; \cos \gamma \\ \end{gathered} &\begin{gathered} \sin \;\eta\; \sin\; \gamma - \\ \cos \;\eta\; \sin\; \beta\; \cos\; \gamma \\ \end{gathered} \\ { - \cos \;\beta\; \sin\; \gamma }&\begin{gathered} \cos\; \eta\; \cos\; \gamma - \\ \sin\; \eta\; \sin\; \beta\; \sin\; \gamma \\ \end{gathered} &\begin{gathered} \sin\; \eta\; \cos\; \gamma + \\ \cos\; \eta\; \sin\; \beta\; \sin\; \gamma \\ \end{gathered} \\ {\sin\; \beta }&{ - \sin\; \eta\; \cos\; \beta }&{\cos\; \eta\; \cos\; \beta } \end{array}} \right]. \end{split} $

式中：${{\boldsymbol{T}}_{{\text{ib}}}}$表示坐标转换矩阵，$\left( {\eta ,\beta ,\gamma } \right)$表示全局惯性系与局部零部件坐标系间的转换角.

轴承内各零件的受力情况由其相对位置及运动关系所决定，本研究模型主要考虑滚子与套圈、滚子端面与套圈挡边、滚子与保持架、套圈与保持架等零件间相互作用，确定各零件的受力状态，以实现对圆柱滚子轴承动力学模型的求解.

滚子所受到的作用力主要包括滚子与套圈间的法向接触力和切向拖动力. 首先须确定滚子与滚道间的接触变形量，进而再计算滚子与套圈间的接触载荷. 滚子与套圈间的接触变形量可以根据Lundberg^[25]提出的经验公式给出：

(2)$ \delta = \frac{{2Q}}{{{\text{π }}l}}\left[ {\frac{{1 - {\nu _1}^2}}{{{E_1}}}+\frac{{1 - {\nu _2}^2}}{{{E_2}}}} \right]\left[ {\ln \;\frac{l}{b}+1.193\;2} \right] ,$

(3)$ b = {\left\{ {\frac{{4Q}}{{{\text{π }}l\;\Sigma \rho }}\left| {\frac{{1 - {\nu _1}^2}}{{{E_1}}}+\frac{{1 - {\nu _2}^2}}{{{E_2}}}} \right|} \right\}^{{1}/{2}}}. $

式中：${v_1}、{v_2}$表示泊松比，${E_1}、{E_2}$表示弹性模量，$Q$表示接触载荷，$l$表示接触长度，$b$表示接触半宽，$\sum \rho $表示曲率半径之和.

如图2所示为滚子与套圈间的相互作用示意图，采用切片法，滚子第k个切片所受的接触力^[26]可以表示为

(4)$ {q_{{k}}} = {K_{\text{E}}}{\delta _{{k}}}^{{{10}}/{9}}+{C_{\text{E}}}{\delta '_{{k}}}. $

式中：${\delta _{{k}}}$表示滚子第k个切片的弹性变形量；${\delta '_{{k}}} $表示滚子的第k个切片的弹性变形率；${K_{\text{E}}}$表示接触总刚度，其为入口区油膜刚度、赫兹接触区接触刚度以及油膜刚度构成的综合刚度；${C_{\text{E}}}$表示接触总阻尼，由于赫兹接触区油膜黏性阻尼和材料阻尼都较低，几乎可以忽略不计，因此总的接触阻尼为入口区油膜阻尼.

(5)$ {K_{\text{E}}} = {k_{{\text{ef}}}}+\frac{{{k_{\text{h}}}{k_{\text{f}}}}}{{{k_{\text{h}}}+{k_{\text{f}}}}}, $

(6)$ {C_{\text{E}}} = {c_{{\text{ef}}}}. $

式中：${k_{{\text{ef}}}}$表示入口区油膜刚度，${k_{\text{h}}}$表示赫兹接触区接触刚度，${k_{\text{f}}}$表示油膜刚度，${c_{{\text{ef}}}}$表示入口区油膜阻尼.

滚子和套圈滚道间的法向作用力^[27]为

(7)$ Q_{{\mathrm{NOR}}} = {{\frac{1}{N}\sum\limits_{{{k = 1}}}^N {{q_{{k}}}} }}. $

式中：N表示切片的数量.

假设滚子与套圈接触变形量一致. 基于五参数流变模型^[28]，滚子与套圈间的润滑拖动系数可根据滚子局域滑动速度及润滑剂类型进行选取，进而可计算出滚子的拖动力. 由此，将作用于滚子的法向力与拖动力合成为

(8)$ {\boldsymbol{F}} = \left[ {0,} \; {\kappa {{Q}_{{\mathrm{NOR}}}},} \; {{{Q}_{{\mathrm{NOR}}}}} \right] ^{\mathrm{T}}.$

式中：$\kappa $表示拖动系数.

滚子及套圈的外力矩分别为

(9)$ {{\boldsymbol{G}}_{\text{R}}} = \left( {{{\boldsymbol{r}}_{\text{p}}}+{{\boldsymbol{r}}_{{\text{Rg}}}}} \right) {\boldsymbol{F}}, \;\; {{\boldsymbol{G}}_{\text{i}}} = - \left( {{{\boldsymbol{R}}_{\text{p}}}+{{\boldsymbol{r}}_{{\text{ig}}}}} \right) {\boldsymbol{F}}. $

式中：${{\boldsymbol{r}}_{\text{p}}}$表示接触区中心在滚子坐标系下的位置向量，${{\boldsymbol{r}}_{{\text{Rg}}}}$表示滚子质心在滚子坐标系下的位置向量，${{\boldsymbol{R}}_{\text{p}}}$表示接触区中心在套圈坐标系下的位置向量，${{\boldsymbol{r}}_{{\text{ig}}}}$表示套圈质心在套圈坐标系下的位置向量.

由于滚子端面与套圈挡边间的相互作用较为复杂，为了进行简化，针对两者间的接触关系提出以下2个假设：1）当滚子端面与套圈挡边不发生接触时，两者间不存在接触力，即忽略油膜拖动作用；2）当套圈挡边与滚子在其圆角处接触时，可用Hertz点接触模型近似计算两者间的接触载荷，并计算拖动力.

基于上述2条假设，根据滚子端面与挡边间的几何相互作用关系，其弹性变形量为

(10)$ \delta_{\mathrm{f}} = {\left( {{\boldsymbol{r}}_{{\text{ef3}}}^{\text{f}}} \right)_{\max }} - {R_{{\text{ce}}}};\; {\delta > 0} . $

式中：$ {\boldsymbol{r}}_{{\text{ef3}}}^{\text{f}} $表示滚子圆角曲率中心相对挡边参考点的位置向量，${\left( {{\boldsymbol{r}}_{{\text{ef3}}}^{\text{f}}} \right)_{\max }} $表示滚子圆角曲率中心相对挡边参考点的最大距离，$ {R_{{\text{ce}}}} $表示滚子圆角半径.

滚子和套圈挡边间的接触载荷可根据Hertz接触理论得出，而润滑剂的拖动力可以表示为

(11)$ {F_{{\text{drag}}}} = {{\frac{1}{N}\sum\limits_{k = 1}^N {{\mu _k}{q_k}} }}. $

式中：$ {\mu _k} $表示滚子第k个切片的油膜拖动系数. $ {\mu _k} $可通过如下公式来确定：

(12)$ \mu_k = \left\{ \begin{gathered} {\mu _{\text{d}}}; \\ {\mu _{{\text{bd}}}} = \left( { - 0.1+22.8s} \right){{\mathrm{e}}^{ - 181.46s}}+0.1; \\ {\mu _{\text{l}}} = {{\left( {{\mu _{{\text{bd}}}} - {\mu _{{\text{hd}}}}} \right){{\left( {{\varDelta _{\text{r}}} - {\varDelta _{{\text{hd}}}}} \right)}^6}} \mathord{\left/ {\vphantom {{\left( {{\mu _{{\text{bd}}}} - {\mu _{{\text{hd}}}}} \right){{\left( {{\varDelta _{\text{r}}} - {\varDelta _{{\text{hd}}}}} \right)}^6}} {{{\left( {{\varDelta _{{\text{bd}}}} - {\varDelta _{{\text{hd}}}}} \right)}^6}}}} \right. } {{{\left( {{\varDelta _{{\text{bd}}}} - {\varDelta _{{\text{hd}}}}} \right)}^6}}}+{\mu _{{\text{bd}}}}; \\ {\mu _{{\text{hd}}}} = \left( {A+Bs} \right){{\mathrm{e}}^{ - Cs}}+D. \\ \end{gathered} \right. $

式中：${\mu _{\text{d}}}$表示干摩擦系数，${\mu _{{\text{bd}}}}$表示边界润滑摩擦系数，${\mu _{\text{l}}}$表示混合润滑摩擦系数，${\mu _{{\text{hd}}}}$表示弹流润滑系数，$ {\varDelta _{\text{r}}} $表示润滑参数，s表示滚滑比，${\varDelta _{\text{bd}}} $表示混合润滑临界润滑参数，${\varDelta _{\text{hd}}} $表示弹流润滑临界润滑参数. $ {\varDelta _{\text{r}}} $可以表示为

(13)$ {\varDelta _{\mathrm{r}}} = {{{h_{\min }}}}{{ {({\sigma _1^2+\sigma _2^2} })^{-1/2}}}. $

式中：${h_{{\text{min}}}}$表示接触区最小油膜厚度，${\sigma _1}、{\sigma _2}$表示两接触表面的表面粗糙度.

因此，滚子端面所受作用力的合力为

(14)$ {{\boldsymbol{F}}_{\text{b}}} = \left[ {0,} \; {\kappa {Q_{\text{f}}},} \; { - {Q_{\text{f}}}} \right]^{\mathrm{T}} .$

式中：$ {Q_{\text{f}}} $表示滚子端面与套圈挡边的接触载荷.

在圆柱滚子轴承运行过程中，滚子与保持架、套圈与保持架间可能存在频繁的碰撞，且单次作用时间较短，使得两者间的接触作用十分复杂. 为了简化接触模型，提出以下3点假设：1）假设保持架为刚性，忽略其柔性变形；2）假设滚子-保持架、套圈-保持架间相互作用不存在流体动压效应造成的弹性变形；3）假设滚子-保持架、套圈-保持架间的局域接触变形为弹性变形.

滚子与保持架的位置关系示意图如图3所示. 图中，C_P表示滚子与保持架兜孔间的距离，$\omega_{\mathrm{c}} $为保持架自转角速度.

滚子与保持架横梁的接触变形^[24]与碰撞力分别如下：

(15)$ {\delta _{\text{c}}} = \left| {{r_{\text{c}}}} \right| - {X_{{k}}}\tan\;\; \beta - C_{\mathrm{P}}, $

(16)$ {Q_{\text{c}}} = {K_{\text{E}}}{\delta _{\text{c}}}^{{{10}}/{9}}+{C_{\text{E}}}{\delta '_{\text{c}}}. $

式中：${\delta _{\text{c}}}$表示滚子与保持架横梁的接触变形，${r_{\text{c}}}$表示滚子与保持架兜孔中心的距离，$\beta $表示滚子的歪斜角，X_k表示兜孔变形量.

滚子与保持架之间的摩擦力为

(17)$ {F_{\text{c}}} = {\mu _{\text{c}}}{Q_{\text{c}}}. $

式中：${\mu _{\mathrm{c}}}$表示滚子与保持架兜孔间的摩擦系数.

滚子与保持架之间由于碰撞力而产生的附加力矩^[26]为

(18)$ {M_{\mathrm{c}}} = \frac{{{Q_{\mathrm{c}}}}}{N}\left( {\frac{l}{2} - \omega k} \right). $

式中：$\omega$表示滚子角速度.

由于保持架不直接承受载荷，须通过其他部件引导来带动其运动. 保持架常见的引导方式有3种：滚子引导、内圈引导以及外圈引导. 以外圈引导为例，保持架与套圈的相互作用如图4所示.

在径向平面上，假定保持架质心坐标为$\left( {{y_{\text{c}}},{z_{\text{c}}}} \right)$，则保持架和引导套圈之间的最小引导间隙为

(19)$ {h_0} = \frac{1}{2}\left( {{D_{{\text{rc}}}} - {D_{{\text{co}}}}} \right) - \left({y_{\text{c}}^2+z_{\text{c}}^{\text{2}}}\right)^{1/2} . $

式中：$ {D_{{\text{rc}}}} $表示引导面直径，$ {D_{{\text{co}}}} $表示保持架外径.

用${\varDelta _{\text{r}}}$作为阈值用来判断接触状态，当${h_0} > {\varDelta _{\text{r}}}$时，此时保持架和引导套圈之间存在间隙，润滑剂和空气的混合物将会进入引导间隙，对保持架和引导套圈产生作用力，将这个力学过程视为短滑动问题来等效处理. 此时，外圈引导面对保持架引导面产生的流体压力的2个正交分量和，可以表示为

(20)$ \left. \begin{gathered} {F_{{\text{c}}{{y}}'}} = - \frac{{\eta U{B_{{\text{rc}}}}^3{\varepsilon ^2}}}{{C_1^2{{\left( {1 - {\varepsilon ^2}} \right)}^2}}}, \\ {F_{{\text{c}}{{z}}'}} = - \frac{{{\text{π }}\eta U{L^3}\varepsilon }}{4C_1^2{\left( {1 - {\varepsilon ^2}} \right)}^{3/2}}. \\ \end{gathered} \right\} $

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

(21)$ {M_{{\text{rc}}x'}} = \frac{{{\text{π }}\eta VD_{{\text{rc}}}^{\text{2}}{B_{{\text{rc}}}}}}{{2{C_1}\sqrt {1 - {\varepsilon ^2}} }}. $

式中：$\eta $表示润滑剂动力黏度，${B_{{\text{rc}}}}$表示引导面宽，$\varepsilon $表示保持架质心的偏移量，${C_1}$表示保持架引导间隙，$U$表示润滑剂混合物的拖动速度，$V$表示相对滑动速度.

当${h_0} \leqslant {\varDelta _{\text{r}}}$时，此时保持架和引导套圈之间存在接触，可采用Hertz线接触问题来处理. 接触力可以表示为

(22)$ \left. \begin{gathered} {F_{{\text{c}}{{y}}'}} = - 0.356{{E'}_{{\text{rc}}}}{B_{{\text{rc}}}}^{{8}/{9}}{\delta ^{{{10}}/{9}}}, \\ {F_{{\text{c}}{{z}}'}} = {\mu _{\text{c}}}{F_{{\text{c}}{{y}}'}}. \\ \end{gathered} \right\} $

式中：${\mu _{\text{c}}} $表示保持架与引导套圈的摩擦系数.

保持架和引导套圈间的摩擦力矩为

(23)$ {M_{{\text{rc}}x'}} = 0.5{D_{{\text{rc}}}}{F_{{{{\mathrm{c}}z'}}}}. $

在阻力估计方面，可以根据圆柱体在流体中平移时所受阻力的表达式^[29]来确定. 其中，阻力系数${C_{\text{D}}}$为雷诺数Re的函数，其与阻力的关系如下：

(24)$ {C_{\text{D}}} = \frac{2D}{{\rho {V_{\mathrm{r}}^2}A'}}, $

(25)$ A' = Ld. $

式中：$D$表示阻力，$\rho $表示润滑剂质量密度，$V_{\mathrm{r}}$表示滚子平移速度，$A'$表示滚子沿流向方向的投影面积，$L$表示滚子长度，$d$表示滚子直径.

在涡动力矩的计算方面，由于滚子在其柱面和端面上均会发生涡动，涡动力矩可以利用Rumbarge^[30]提出的圆柱体涡动模型进行计算. 其中，滚子在其柱面上的涡动力矩${M_{{\mathrm{cc}}}}$可表示为

(26)$ {M_{{\text{cc}}}} = \tau Ar_\tau. $

式中：$\tau $表示柱面上的剪切力，$A$表示柱面面积，$r_\tau$表示柱面距离转动中心的距离.

滚子端面上的涡动力矩可以表示为

(27)$ {M_{\text{e}}} = \rho \omega {r^5}{C_{\text{n}}}/2. $

式中：$r$表示滚子半径；$\omega $表示滚子自转的角速度； ${C_{\text{n}}}$表示端面涡动系数，其取值与雷诺数有关.

根据轴承各零件间的相互作用力和运动状态，分别建立圆柱滚子轴承内部各零件的动力学微分方程组.

滚子的动力学微分方程为

(28)$ \left. \begin{split} {m_{\text{b}}}{{\ddot y}_{{\text{b}}j}} =& - N_j^{\text{i}}\cos\; {\varphi _j}+N_j^{\text{o}}\cos\; {\varphi _j} - {Q_{{\text{c}}j}}\sin\; {\varphi _j} - \\ &{F_{{\text{r}}j}}\cos\; {\varphi _j} - T_j^{\text{i}}\sin\; {\varphi _j}+T_j^{\text{o}}\sin\; {\varphi _j} - {F_{{\text{c}}j}}\cos\; {\varphi _j}, \\ {m_{\text{b}}}{{\ddot z}_{{\text{b}}j}} =& N_j^{\text{i}}\sin\; {\varphi _j} - N_j^{\text{o}}\sin\; {\varphi _j} - {Q_{{\text{c}}j}}\sin\; {\varphi _j}+ \\ &{F_{{\text{r}}j}}\sin\; {\varphi _{{j}}} - T_j^{\text{i}}\cos\; {\varphi _j}+T_j^{\text{o}}\cos\; {\varphi _j}+{F_{{\text{c}}j}}\sin\; {\varphi _j}, \\ {J_{{\text{b}}x}}{{\dot \omega }_{{\text{b}}x}} = &- T_j^{\text{i}}\frac{{{D_{\text{w}}}}}{2} - T_j^{\text{o}}\frac{{{D_{\text{w}}}}}{2}+{F_{{\text{c}}j}}\frac{{{D_{\text{w}}}}}{2}, \\ {J_{{\text{b}}y}}{{\dot \omega }_{{\text{b}}y}} =& - M_{{\text{T}}\;j}^{\;\text{i}}\sin\; {\varphi _j} - M_{{\text{T}} j}^{\text{o}}\sin\; {\varphi _j}+ \\& M_{{\text{C}} j}^{\;\text{i}}\cos\; {\varphi _j}+M_{{\text{N}} j}^{\;\text{i}}\sin\; {\varphi _j}+M_{{\text{N}} j}^{\text{o}}\sin\; {\varphi _j}, \\ {J_{{\text{b}}z}}{{\dot \omega }_{{\text{b}}z}} =& - M_{{\text{T}}j}^{\;\text{i}}\cos\; {\varphi _j} - M_{{\text{N}} j}^{\text{o}}\cos\; {\varphi _j} - \\& M_{{\text{C}} j}^{\;\text{i}}\sin\; {\varphi _j}+M_{{\text{N}} j}^{\;\text{i}}\cos\; {\varphi _j}+M_{{\text{N}} j}^{\text{o}}\cos\; {\varphi _j}. \end{split} \right\} $

式中：${m_{\text{b}}}$表示滚子质量，${\ddot y_{{\text{b}}j}}、{\ddot z_{{\text{b}}j}}$为滚子的质心加速度，${J_{{\text{b}}x}}、{J_{{\text{b}}y}}、{J_{{\text{b}}z}}$表示滚子的转动惯量，${\dot \omega _{{\text{b}}x}}、{\dot \omega _{{\text{b}}y}}、 {\dot \omega _{{\text{b}}z}}$表示滚子的角加速度，$N_j^{\;\text{i}} $、$N_j^{\;\text{o}} $分别表示滚子与内、外滚道之间的法向作用力，$Q_{{\mathrm{c}} j} $表示保持架与滚子之间的法向接触力，$F_{{\mathrm{c}} j} $表示保持架与滚子之间的切向摩擦力，$F_{\mathrm{r}j} $表示滚子离心力，$T_j^{\text{i}} $、$T_j^{\text{o}} $分别表示滚子与内、外滚道之间的油膜拖动力，$M_{{\text{T}} j}^{\;\text{i}}、M_{{\text{T}} j}^{\text{o}} $分别表示滚子与内、外滚道油膜拖动力产生的附加力矩，$M_{{\text{N}} j}^{\;\text{i}}、M_{{\text{N}} j}^{\text{o}}$分别表示滚子内、外滚道接触力产生的附加力矩，$M_{\mathrm{c} j}^{\mathrm{i}} $表示保持架与滚子接触力产生的附加力矩，$D_{\mathrm{w}} $表示滚子的阻力.

内圈的动力学微分方程为

(29)$ \left. \begin{split} {m_{\text{i}}}{{\ddot y}_{\text{i}}} =& \mathop \sum \nolimits_{j = 1}^Z {Q_{{\mathrm{c}}j}}\cos\; {\varphi _j}+{F_{{\mathrm{c}}j}}\sin\; {\varphi _j} - {F_{\text{r}}}, \\ {m_{\text{i}}}{{\ddot z}_{\text{i}}} =& \mathop \sum \nolimits_{j = 1}^Z {{Q_{{\mathrm{c}}j}}\cos\; {\varphi _j}+{F_{{\mathrm{c}}j}}\sin\;{\varphi _j}} , \\ {J_{{\text{i}}x}}{{\dot \omega }_{{\text{i}}x}} =& \mathop \sum \nolimits_{j = 1}^Z {F_{{\mathrm{c}}j}}{{{D_{\text{w}}}}}/{2}, \\ {J_{{\text{i}}y}}{{\dot \omega }_{{\text{i}}y}} =& \mathop \sum \nolimits_{j = 1}^Z \left( {M_{{\text{T}}j}^{\;\text{i}}\sin\; {\varphi _j}} \right), \\ {J_{{\text{i}}z}}{{\dot \omega }_{{\text{i}}z}} =& \mathop \sum \nolimits_{j = 1}^Z \left( {M_{{\text{T}}j}^{\;\text{i}}\cos\; {\varphi _j}} \right). \end{split}\right\} $

式中：${m_{\text{i}}}$表示内圈的质量，${\ddot y_{\text{i}}}、{\ddot z_{\text{i}}}$表示内圈的质心加速度，$Z$表示滚子个数，${J_{{\text{i}}x}}、{J_{{\text{i}}y}}、{J_{{\text{i}}z}}$表示内圈的转动惯量，${\dot \omega _{{\text{i}}x}}、{\dot \omega _{{\text{i}}y}}、{\dot \omega _{{\text{i}}z}}$表示内圈的角加速度.

保持架的动力学微分方程为

(30)$ \left.\begin{split} m_{\mathrm{c}} \ddot{y}_{\mathrm{c}}=&\sum_{j=1}^Z\left(Q_{\mathrm{c}j} \sin\; \varphi_j+F_{\mathrm{c}j} \cos\; \varphi_j\right)+ \\&F_{\mathrm{c}y'} {\cos\; \psi_{\mathrm{c}}+F_{\mathrm{c}z'} \sin\; \psi_{\mathrm{c}}-G_{\mathrm{c}} ,} \\m_{\mathrm{c}} \ddot{{z}}_{\mathrm{c}}=&\sum_{j=1}^Z\left(Q_{\mathrm{c}j} \cos\; \varphi_j-F_{\mathrm{c}j} \cos\; \varphi_j\right)+ \\&F_{\mathrm{c}y'} \sin\; \psi_{\mathrm{c}}-F_{\mathrm{c}z'} \mathfrak{\operatorname { c o s }\; \psi _ { \mathrm { c } },} \\J_{\mathrm{c}x} \dot{\omega}_{\mathrm{c}x}=&\sum_{j=1}^Z\left(F_{\mathrm{c} j} \frac{D_{\mathrm{w}}}{2}\right)-M_{\mathrm{c}x'}-T_{\mathrm{CDO}}-T_{\mathrm{CDs}} .\end{split}\right\} $

式中：${m_{\text{c}}}$表示保持架质量，${\ddot y_{\text{c}}}、{\ddot z_{\text{c}}}$表示保持架的质心加速度，${J_{{\text{c}}x}}$表示保持架的转动惯量，${\dot \omega _{{\text{c}}x}}$表示保持架的角加速度，${T_{{\text{CDO}}}}$表示保持架表面的阻滞力矩，${T_{{\text{CDs}}}}$表示保持架侧面的阻滞力矩，$M_{\mathrm{c}x'} $表示保持架与外圈产生的力矩.

针对上述动力学微分方法，采用Hilber-Hughes-Taylor变步长积分算法进行计算.

在圆柱滚子轴承发生打滑时，其滚子及保持架的实际转速降低，低于其无滑移条件下的理论转速，因此可以利用基于实际转速计算得到的保持架打滑率及滚子打滑率对轴承的打滑程度进行描述. 保持架打滑率是目前圆柱滚子轴承打滑分析所采用的主流评价指标，可以反映轴承的宏观打滑现象. 保持架打滑率的表达式为

(31)$ {{{S}}_{\text{c}}}{\text{ = }}\frac{{{n_{\text{c}}} - {n'_{\text{c}}}}}{{{n_{\text{c}}}}} \times 100{\text{%}}. $

式中：$ {n'_{\text{c}}} $表示保持架实际转速，$ {n_{\text{c}}} $表示保持架无滑移条件下的理论转速.

保持架的理论转速为

(32)$ {n_{\text{m}}} = \frac{1}{2}\left[ {{n_{\text{i}}}\left( {1 - \gamma } \right)+{n_{\text{o}}}\left( {1+\gamma } \right)} \right], $

(33)$ \gamma = {d}/{{{d_{\text{m}}}}}. $

式中：${n_{\text{i}}}$表示轴承内圈转速，${n_{\text{o}}}$表示轴承外圈转速，${d_{\text{m}}}$表示轴承节圆直径.

当轴承外圈不发生转动时，保持架理论转速可以表示为

(34)$ {n_{\text{c}}} = \frac{{{n_{\text{i}}}}}{2}\left( {1 - \gamma } \right). $

滚子打滑率可以更加直观地描述轴承的局域打滑现象. 与保持架打滑率表达式类似，滚子打滑率可以表示为

(35)$ {S_{\text{R}}} = \frac{{{n_{\text{r}}} - {n'_{\text{r}}}}}{{{n_{\text{r}}}}} \times 100{\text{%}} , $

(36)$ {n_{\text{r}}} = \frac{{{d_{\text{m}}}{n_{\text{i}}}}}{{2D}}\left( {1 - {\gamma ^2}} \right). $

式中：$ {n'_{\text{r}}} $表示滚子实际自转速度，$ {n_{\text{r}}} $表示滚子在无滑移条件下的理论自转速度.

受圆柱滚子轴承其自身结构及受载特性影响，其滚子在圆周上受载情况呈周期性变化. 当滚子进入承载区后，滚子所受接触力逐渐增大；反之则逐渐减小. 由于这种周期性受载特点，滚子转速n_r呈周期性变化（见图5 (a)），这与保持架转速n_c在一定时间范围内的小幅波动差异性较大（见图5(b)）. 因此，可以通过预测滚子的瞬时转速，来研究滚子打滑，即通过研究周期性波动的滚子打滑率，来分析轴承的局部打滑现象；对保持架打滑的研究则是通过对其平均转速的预测.

在后续分析中，统一选用N216 型圆柱滚子轴承，其主要结构参数如表1所示，轴承保持架材料为尼龙. 轴承内圈与主轴间采用过盈配合，过盈量为10 μm；轴承外圈与轴承座间采用间隙配合，过盈量为5 μm，原始径向游隙为1 μm，轴承的径向工作游隙为−4 μm. 选择润滑剂为美孚DTE型润滑油，其密度为0.85 kg/L，40℃时的运动黏度为31.0 mm²/s，100 ℃时的运动黏度为5.5 mm²/s，黏度指数为102. 相比于文献中设置的径向载荷为300 N、转速为4 500 r/min的高速轻载工况，本研究主要针对中低速工况条件展开研究，此时轴承保持架的宏观打滑现象并不明显，在传统研究过程中往往被忽略. 基于文献^[16-19]的载荷工况，选择径向载荷为1 000 N、转速为1 000 r/min的工况条件，作为本研究的中低速工况条件，转速在4 500 r/min及以上时，视为高速工况，并开展此条件下的圆柱滚子轴承内部滚子局部打滑的影响因素及变化规律研究.

对主轴转速为1000 r/min、径向载荷为1000 N工况下的仿真结果进行分析，可以得到滚子打滑率S_R在一个转动周期内的变化情况，即滚子打滑率随其相位角的变化情况，如图6所示. 滚子与内外圈相对滑动速度v_r、接触载荷L与相位角$\theta $的关系如图7所示. 图中，相位角范围80°~280°为承载区，滚子在承载区、非承载区及由非承载区进入承载区的“咬入”过程中，打滑情况存在较大差异. 滚子在非承载区打滑率较小，而在滚子进入“咬入”过程中打滑率急剧增大，并在完全进入承载区后迅速降低至零，打滑率在承载区边缘形成峰值. 结合图7分析，滚子与内、外圈之间均在打滑率峰值对应的位置存在较大的相对滑动速度，与Shao等^[31]的研究结论一致. 滚子在其进入承载区的跨区运动过程中出现打滑突增现象，可能是由于其受载情况发生剧烈变化，此时滚子的运行稳定性下降，滚子与内外圈之间的滑动加剧，导致滚子与内外圈之间同时发生打滑. 这表明滚子打滑率是滚子与内外圈间相对滑动的一个综合结果，可以在一定程度上反映轴承的局部打滑情况.

如图7所示，本研究所绘图像的深色区域代表轴承处于承载区. 在中低速中等载荷工况下，滚子在承载区内几乎不发生打滑，但滚子在非承载区及承载区两端的打滑情况仍然存在. 通过与文献[32]得到的研究结果进行对比，发现这与传统的基于保持架打滑率的轴承打滑研究结论存在较大差异，即当圆柱滚子轴承承受的径向载荷高于临界载荷时，保持架几乎不发生打滑. 因此，仅仅依靠保持架打滑率这一单一指标对圆柱滚子轴承的打滑情况进行描述是不全面、不可靠的. 同时，圆柱滚子轴承常常被应用在中低速、中等载荷工况下，有必要针对滚子局部打滑率进行计算.

滚子在进入承载区时，所受的载荷和转速具有更加显著的时变特性，这使得打滑机理较为复杂. 为了对滚子及保持架的打滑机理进行更深入的研究，对滚子与保持架兜孔间接触力L_cr进行分析，如图8所示.

滚子在进入承载区时，打滑率急剧增大，滚子与保持架兜孔间发生频繁碰撞. 受到碰撞冲击载荷影响，滚子的质心角加速度a_r出现振荡，如图9(a)所示，滚子的自转角速度$\omega_{\mathrm{r}} $如图9(b)所示. 而在其他区域内，滚子与保持架间的碰撞频率相对较低，滚子运动更为平稳. 滚子在承载区、非承载区以及跨区“咬入”过程中，滚子与保持架间的碰撞作用有较大差异，从而进一步验证了采用滚子打滑率来描述轴承的局部打滑特性的合理性.

滚动轴承的承载区位置和范围与轴承所受径向载荷的方向和大小有关，因此，轴承内部打滑的位置及剧烈程度应与轴承所受的径向载荷大小有关. 设置分析的工况条件如下：主轴转速为1000、5000 r/min，径向载荷为100~5000 N.

滚子在圆周上的打滑率随径向载荷变化规律如图10所示，无论轴承处于高速工况还是低速工况，滚子在圆周上的打滑率均随径向载荷F_r增大而减小. 相比于高速工况下，滚子整周出现打滑现象（见图10(b)），低速时滚子打滑率在承载区的变化更加剧烈（见图10(a)）.

当轴承承受的径向载荷小于500 N时，滚子在轴承整个周期内均存在打滑现象，滚子打滑率在进入承载区后减小，在离开承载区前达到最小值，并在离开承载区后逐渐增大. 当径向载荷大于500 N（即临界载荷）时，滚子在非承载区存在打滑现象，而在完全进入承载区后几乎不发生打滑. 此时轴承不存在整体打滑，即保持架打滑率S_C几乎为0，如图11所示.

当径向载荷为500~3 000 N时，滚子在跨区过程中出现打滑波动，滚子打滑率在进入承载区前迅速增大，在承载区入口达到峰值，并在进入承载区后迅速减小，打滑率波动范围均随径向载荷的增大而减小；当轴承受到的径向载荷大于3 000 N时，滚子在非承载区打滑率缓慢增大，随着滚子进入承载区，滚子打滑率迅速降低至0，滚子打滑率的波动范围较小，滚子的局部打滑现象不显著. 以往工作并没有聚焦于低速工况下轴承打滑特性的研究，然而对于圆柱滚子轴承而言，其多服役于较低转速区间的工况条件下. 因此，研究重载低速条件下的轴承打滑特性具有重要意义.

如图12所示为轴承打滑情况随径向载荷加载方向的变化. 图中，$\theta_{\mathrm{L}} $表示加载位置角度. 当圆柱滚子轴承承受较大载荷或者处于较低的转速运行状态时，保持架平均打滑率出现了负值，即保持架实际转速大于纯滚动条件下的理论转速，这种情况被称为保持架的过度打滑^[10]. 这种现象的发生可能与径向载荷的加载方向及保持架自身重力的影响有关. 当加载方向与重力方向一致时，保持架打滑率及滚子打滑率最小，且保持架发生过度打滑现象. 当加载方向与重力方向相反时，保持架打滑率最大. 同时，当加载方向为270°时，滚子在进入承载区过程中出现滚子打滑率峰值，从而验证了上述猜想. 由于主轴旋转方向为顺时针，在滚子进入承载区过程中，滚子运动受其自身重力作用影响，其速度会降低. 同时，滚子打滑率会出现局部增大的现象.

对于外圈静止，且内圈与主轴过盈配合的圆柱滚子轴承而言，滚子及保持架打滑率的大小均与主轴转速大小相关. 由于受离心力作用，在高速轻载工况下，随着转速的增大，滚子在非承载区与外圈的相对滑动速度增大，在所有滚子的共同作用下，保持架的打滑率随主轴转速增大而增大. 而在中低速工况下，保持架几乎不存在整体打滑，鲜有研究涉及主轴转速对滚子局部打滑情况的研究. 为了探究滚子打滑率随主轴转速的变化规律，设置分析的工况条件如下：径向载荷为300 、1000 N，主轴转速为1000~8000 r/min.

圆柱滚子轴承内部滚子打滑率在其圆周方向上随主轴转速变化的规律如图13所示. 对于中等载荷（1000 N）条件下的圆柱滚子轴承（见图13(a)），滚子在进入承载区前，其打滑率逐渐增大，并在咬入过程中发生突增并迅速达到峰值，在滚子完全进入承载区后迅速减小，其滚子在承载区内打滑率几乎为0. 随着主轴转速的增大，在非承载区内滚子打滑率减小；在非承载区内滚子几乎不发生打滑，也几乎不受到主轴转速变化的影响，但滚子在跨区运动过程中的打滑率峰值则随之降低. 轻载（300 N）条件下的滚子打滑规律则有所不同，滚子在整个轴承圆周方向上均存在打滑现象. 随着主轴转速的增大，滚子在承载区内的打滑率增大，而在非承载区内的打滑率减小，滚子打滑率的波动范围减小.

分析滚子与内圈间的相对滑动速度n_i，如图14所示，当轴承处于轻载条件时，在整个承载区内滚子与内圈间均发生了较大的相对滑动；在中等载荷条件下时，滚子与内圈间仅在承载区的边缘位置具有较大的相对滑动速度.

保持架打滑率与主轴转速n_a关系如图15所示. 当轴承承受的径向载荷小于临界载荷时，保持架打滑率随着主轴转速增大而增大；当主轴转速较小或径向载荷大于临界载荷时，保持架几乎不发生打滑. 因此，为了能够综合反映出圆柱滚子轴承的打滑特性，须选用滚子打滑率和保持架打滑率同时作为评价指标. 但对于不同主轴转速下的圆柱滚子轴承，应该优先选用不同的评价指标对轴承打滑情况进行描述. 对于应用于中低速工况下的圆柱滚子轴承，由于其保持架不表现出宏观打滑现象，优先选用滚子打滑率作为轴承打滑评价指标；对于高速圆柱滚子轴承，其滚子在圆周上的打滑情况波动较小，则应优先选用保持架打滑率作为评价指标.

无论圆柱滚子轴承处于轻载条件还是中等载荷条件，滚子打滑率波动范围随主轴转速的增大而减小. 当转速大于5000 r/min时，滚子打滑率波动性几乎消失，在整个周期上均不发生打滑. 这种现象可能与离心力作用有关. 随着主轴转速的增大，滚子所受离心力增大，导致滚子与外圈间接触载荷L_O增大，尤其是在非承载区（如图16所示），而滚子与内圈间接触载荷大小几乎不随主轴转速升高而发生变化.

滚子所受载荷大小受到内外滚道接触载荷的综合作用，随着主轴转速的增大，滚子在承载区及非承载区间的载荷波动减小. 当转速大于5000 r/min时，滚子在承载区及非承载区内所承受的载荷整体波动较小，进而缓解了滚子打滑率的波动. 同时，由于外圈接触载荷的增大，滚子所受拖动力增大，滚子与外圈间的相对滑动减少，滚子实际转速更接近无滑移条件下的理论转速，而保持架的打滑情况也验证了该结论. 此外，本研究研究结果与文献[33]实验结果进行了对比，可初步验证本研究结果的合理性.

通过对圆柱滚子轴承进行动力学建模及计算，分析轴承内部滚子的局部打滑现象，讨论滚子打滑与保持架打滑的关系，并揭示滚子打滑的产生原因与变化规律，分析滚子打滑率作为圆柱滚子轴承打滑评价标准的适用范围及必要性，以及保持架运动对于滚子打滑的影响. 结论如下：

(1) 滚子由于其承受的接触载荷的变化产生周期性的波动，其在承载区内打滑率较小，离开承载区后打滑率增大，在重新进入承载区过程中迅速增大，完全进入承载区后又迅速降低，形成打滑率峰值现象.

(2) 滚子打滑率在整个圆周上均随径向载荷的增大而减小；滚子打滑率在其非承载区及跨区过程中随主轴转速的增大而减小，当轴承处于高速工况下时，滚子在承载区内打滑率均随主轴转速的增大而增大.

(3) 对于应用于中低速工况下的圆柱滚子轴承，由于其保持架不表现出宏观打滑现象，优先选用滚子打滑率作为轴承打滑评价指标；对于高速圆柱滚子轴承，其滚子在圆周上的打滑情况波动较小，应优先选用保持架打滑率作为评价指标.

（4）本研究暂未详细介绍实验验证工作，本研究对应的实验验证工作主要是基于电磁感应原理，通过在轴承上安装径向充磁磁铁，切割环形线圈，进而测量感应电压并计算滚子打滑率，后续将对这部分内容进一步展开.