行星滚柱丝杠小转角运动的动力学分析

王 乐，李海泓，张 玲

(中国空空导弹研究院，河南 洛阳 471009)

0 引 言

滚柱丝杠按其用途可分为行星滚柱丝杠(RSM)、差动式滚柱丝杠及循环式滚柱丝杠，其中行星滚柱丝杠具有很高的承载能力，能适应极高的速度加速度，且具有很高的稳定性［1］，因此在传动进给系统中得到广泛的应用。行星滚柱丝杠主要由五部分组成:丝杠、滚柱、螺母、持环及齿轮，如图1 所示。

行星滚柱丝杠是把转动运动形式转化为轴向运动形式的机构，其与行星滚珠丝杠(BSM)结构相似，区别在于行星滚柱丝杠的载荷传递元件为螺纹滚柱，而非滚珠，相比于行星滚珠丝杠，行星滚柱丝杠的基础表面具有更大的曲率半径，而且接触表面的数量也较多，这决定了行星滚柱丝杠拥有高寿命、大载荷及高速度运行能力。

图1 行星滚柱丝杠结构示意图

行星滚柱丝杠一般应用在关键的、高精度机械里，如医疗器械［2－4］、航空航天工业［5－6］、光学装置［7］、机器人及高精密度机床［8－9］。国内外对行星滚柱丝杠的研究日趋成熟，能够得到广泛的工业应用是建立在对行星滚柱丝杠效率、失效形式以及动态载荷试验的深入研究的基础之上;Velinsky［10］对行星滚柱丝杠的运动进行了研究，并给出了其效率随螺纹螺旋角及接触角变化规律。现有的研究主要集中在滚柱丝杠大转角运动特性［11－12］，对小转角运动中的动力学及受力特性研究还是一片空白。本文首先通过理论分析建立行星滚柱丝杠的数学模型，从弹性变形及传动关系研究行星滚柱丝杠运行过程;其次建立行星滚柱丝杠的有限元模型，通过有限元显式算法对行星滚柱丝杠进行动力学仿真，研究行星滚柱丝杠动力学响应。最终获得行星滚柱丝杠的运动特性和动力学特性，为其工业应用提供了依据。

1 行星滚柱丝杠数学模型

行星滚柱丝杠在运行时，尤其是在小信号的丝杠转动过程中，其轴向位移可以分为啮合接触弹性变形量和传动轴向位移两部分。

1.1 啮合接触弹性变形量

1.1.1 行星滚柱丝杠的赫兹理论模型

1896年，Hertz 在经典赫兹弹性接触理论基础上，进一步提出了关于两弹性体点接触的局部应力和变形的经典解［13］，经典赫兹理论建立在如下假设条件上:

a. 两个接触表面光滑，只能够产生弹性变形，且服从胡克定律。对于行星滚柱丝杠，可以假设各部件的加工表面是光滑的，在额定载荷工作条件下，接触点处产生的塑性变形量不超过各转动体当量半径的千分之一，这样可以认为行星滚柱丝杠工作在弹性与塑性的临界点处。

b. 等效接触面的尺寸远小于接触的转动体的表面曲率半径。行星滚柱丝杠的啮合运动中，丝杠与滚柱的接触、滚柱与螺母的接触都属于点接触，充分接触之后相当于小面接触，接触面相比接触体曲率半径相当小，满足假设。

c. 考虑接触表面间的法向压力载荷，不考虑接触表面之间的切向摩擦力。行星滚柱丝杠属高精度传动装置，可以忽略各传动部件之间的摩擦力，把啮合传动过程简化为正压力接触传动。

用赫兹接触理论方法求解点接触，可以采用查表法获得各参数数值，也可以根据接触模型计算出的主曲率F(ρ)的值，经过迭代或插值积分得到系数ma，mb，K(e)，L(e)的值，再代入应力应变表达式中求出接触问题的解。本文为了获得较精确的赫兹接触理论解，将结合数值积分、数值迭代两种方法对赫兹理论的偏心方程进行求解，得到等效接触椭圆的偏心率，进而得到接触椭圆长短半轴，再代入到应力应变公式中得到赫兹理论解。计算时，输入的参数有三个:行星滚柱丝杠基本结构参数、材料属性及点接触法相压力载荷。

对于行星滚柱丝杠，简化的赫兹点接触模型如图2 所示。在载荷Q 的作用下，螺纹啮合接触点将扩展成为一个接触面。丝杠－滚柱啮合点接触及滚柱－螺母啮合点接触部分可以简化为一个椭圆，椭圆长轴为2a，短轴为2b。

图2 行星滚柱丝杠赫兹点接触模型

长半轴a 及短半轴b 的表达式如下:

接触椭圆的长短半轴系数分别为

式中:k 为椭圆率，k=b/a;K(e)为第一类完全椭圆积分;L(e)为第二类完全椭圆积分。k(e)，L(e)可分别由下式求出:

根据赫兹理论，得到

式中:Q 为法向压力;δ 为接触弹性趋近量;E'为当量弹性模量;E1，E2为两接触物体的弹性模量;μ1，μ2为两接触物体的泊松比;∑ρ 为两物体在接触点处的主曲率的和;椭圆偏心率e 与椭圆率k的关系为

主曲率函数F(ρ)由下式计算得到:

也可以表示为

若已知接触物体在接触点处的各个主曲率，则可以由式(7)或式(8)得到主曲率函数，然后进行数值迭代得到偏心率e，再得到长短半轴系数ma，mb，最后得到弹性趋近量δ。

1.1.2 啮合接触面曲率计算

行星滚柱丝杠副中丝杠与滚柱接触时，其第一、二主曲率分别为

行星滚柱丝杠副中滚柱与螺母接触时，其第一、二主曲率分别为

行星滚柱丝杠结构基本参数如表1 所示。

表1 行星滚柱丝杠基本结构参数

由公式(9)得丝杠－滚柱接触主曲率ρ11=0.357 1，ρ12=0.246 6，ρ21=0，ρ22=0.108 9，进而可得丝杠－滚柱的主曲率函数如下:

由公式(10)得滚柱－螺母接触主曲率ρ11=0.357 1，ρ12=0.246 6，ρ21=0，ρ22=0.042 9，进而可得滚柱－螺母的主曲率函数如下:

1.1.3 弹性变形量计算

行星滚柱丝杠材料选用42CrMo4，其材料性能参数:93 ℃时，弹性模量210 GPa，密度7 800 kg/m3，泊松比0.29。

对第一、第二类完全椭圆积分K(e)，L(e)，采用龙贝格(Romberg)积分法进行数值积分，积分表达式如下:

式中:a=0;b=π/2;Tm，k为2m－2 阶牛顿－科特斯(Newton－Cotes)公式计算结果;对于第一类完全椭圆积分对于第二类完全椭圆积分

式中F(ρ)由式(11)～(12)计算得到，式(8)为关于椭圆偏心率e 的方程，编制数值积分及数值迭代的Matlab 程序求解得到e 值，将得到的e 值代入到式(5)中，可得最大弹性变形量δ。

a. 丝杠－滚柱点接触理论求解

不同载荷轴向Q 作用下，丝杠－滚柱接触椭圆面的长、短轴随载荷的变化曲线如图3 所示，接触点弹性趋近量随载荷的变化曲线如图4 所示。

b. 滚柱－螺母点接触理论求解

不同载荷轴向Q 作用下，滚柱－螺母接触椭圆面的长、短轴随载荷的变化曲线如图5 所示，接触点弹性趋近量随载荷的变化曲线如图6 所示。

图3 长、短轴变化曲线

图4 弹性变形量变化曲线

图5 长、短轴变化曲线

图6 弹性变形量变化曲线

1.2 传动轴向位移

行星滚柱丝杠三个转动部件的节圆半径影响相对转速的大小，各部件的螺纹导程及螺纹旋向决定螺母的轴向进给量。丝杠与滚柱之间存在轴向位移差，但由于螺母与滚柱螺纹节距匹配及螺母持环的约束，滚柱相对于螺母无相对轴向位移，只做平面转动。

1.2.1 行星滚柱丝杠无滑移传动

行星滚柱丝杠副无滑移传动形式如图7 所示，丝杠与滚柱运动过程中相切于节圆切点B，滚柱与螺母相切于节圆切点A。

图7 行星滚柱丝杠无滑移传动

由于螺母无周向速度，故A 点周向速度为零;由于丝杠有自转角速度ωS，故滚柱B 点的轴向速度如下式:

由运动关系可得，滚柱的圆心速度:

假设丝杠转过θS角度，其节圆上B 点转过的弧长为sB=rSθS，由圆心速度方程可知，滚柱中心O 转过弧长为

约束点B 相应的转角如下式:

式中:a 为丝杠与滚柱的中心距。

故行星滚柱丝杠的传动关系为

式中:l 为轴向传动位移;pS为丝杠导程;pR为滚柱导程;当丝杠与滚柱螺纹旋向相同时取正号，当丝杠与滚柱螺纹旋向相反时取负号。

行星滚柱丝杠的结构参数代入式(18)，且本文采用的行星滚柱丝杠副丝杠与滚柱螺纹旋向相反，取负号得

1.2.2 行星滚柱丝杠有滑移传动

当考虑行星滚柱丝杠副运动过程中丝杠与滚柱接触啮合滑动时，丝杠相对传动螺母转角为θ'SC，如图8 所示。

图8 行星滚柱丝杠有滑移传动

丝杠转动2π 角度时，丝杠、滚柱、螺母转角的关系满足下式:

式中:θ'R为滚柱转角;θ'SC为丝杠相对于螺母转角;θ'S为丝杠转角。

丝杠实际转过角度与丝杠滑动转角之间的关系如下式:

结合式(20)～(21)有:

行星滚柱丝杠运动过程中，考虑啮合滑动，滚柱转角应小于没有啮合滑动条件下的转角，丝杠滑动转角小于没有啮合滑动条件下的转角，用下面条件表示:

代入到上式有:0 ≤θ'SC≤θSC，查有关行星滚柱丝杠技术资料可知，在额定载荷下，丝杠滑动量与转动量之间满足关系式θslide= 0.005θS，有:

由式(24)可得行星滚柱丝杠副传动关系如下式:

当丝杠与滚柱螺纹旋向相同时取正号，当丝杠与滚柱螺纹旋向相反时取负号。结合行星滚柱丝杠副结构参数得

2 行星滚柱丝杠有限元分析

2.1 行星滚柱丝杠有限元模型

行星滚柱丝杠的结构复杂，螺纹数量多，考虑到有限元计算效率，本文对行星滚柱丝杠的模型进行了简化，如图9 所示，取等长度的丝杠、滚柱及螺母接触部分进行建模，模型中略去了滚柱中心圆柱部分及丝杠中心圆柱部分。

图9 行星滚柱丝杠有限元模型

考虑行星滚柱丝杠螺纹结构的复杂性，螺纹接触部分采用四面体网格划分，设置C3D4 网格单元。每条螺纹面上三层网格，尽量减少网格数量，总网格数量为67 944。

定义行星滚柱丝杠材料密度为7.8E－9 ton/mm3，弹性模量为210 000 MPa，泊松比为0.29。螺纹啮合接触采用“surface－to－surface”小滑移接触类型，切向库伦摩擦系数为0.3，法向设置为硬接触。

设置两个分析步，第一个分析步施加预紧拉力，如图10 所示，使螺纹啮合充分接触，第二个分析步施加丝杠扭矩，如图11 所示，在扭矩作用下使滚柱、螺母转动。第一个分析步中释放滚柱及螺母轴向方向自由度，约束三个部件其他所有自由度，对螺母施加轴向预紧力，滚柱和螺母产生轴向位移;第二个分析步在第一个分析步基础上，释放丝杠、滚柱及螺母转动方向自由度，在丝杠扭矩作用下三个部件产生角位移，滚柱和螺母有轴向位移。

图10 轴向预紧力曲线

图11 丝杠扭矩曲线

2.2 行星滚柱丝杠有限元仿真

对螺母施加如图10 所示的预紧力，对丝杠施加如图11 所示的扭矩载荷，总计算时间为0.001 s。计算完成后提取两个时间点行星滚柱丝杠应力分布云图，如图12 ～13 所示。

图12 0.000 313 6 s 时刻米塞斯应力分布

图13 0.000 704 s 时刻米塞斯应力分布

从行星滚柱丝杠应力分布云图可以看出，行星滚柱丝杠副运动过程中螺纹啮合处接触比较充分，各部件的应力分布比较均匀。从表2 中各部件最大米塞斯应力可知，丝杠与滚柱比螺母的最大应力大，这是因为丝杠与滚柱丝杠侧的接触螺纹的曲率较大，接触更趋近于点接触，0.000 313 6 s 时刻滚柱丝杠侧及0.000 704 s 时刻丝杠与滚柱最大应力超过了材料屈服应力，材料局部进入塑性屈服。

表2 行星滚柱丝杠最大应力分布

在0.000 2 s 之后，丝杠在扭矩作用下开始转动，图14(a)～(c)所示分别为丝杠角位移曲线、滚柱角位移曲线及螺母角位移曲线。滚柱的自转方向与丝杠转动方向相反，螺母与丝杠的转动方向相反。

图14 行星滚柱丝杠动力学响应曲线

图14(d)～(e)分别为滚柱轴向位移曲线、螺母轴向位移曲线。0.000 2 s 时刻，螺母参考点受到的轴向拉力为100 N，平均分布在每个螺纹接触点处的拉力为20 N，从图中数据可得0.000 2 s 时刻滚柱轴向位移为0.001 640 42 mm，与赫兹理论分析得到的弹性变形量0.001 6 mm 相比，相对误差为2.5%。0.000 2 s 时刻之后，滚柱及螺母最大正向轴向位移在0. 000 284 8 s 时刻达到，分别为0.001 767 59 mm，0.003 402 09 mm。在实际进给传动系统中，在预紧力作用下的螺母的正向轴向位移需作为补偿引入到运动模型中。

图14(f)为螺母轴向位移在0.000 2 s 之后随丝杠角位移变化曲线，即行星滚柱丝杠传动关系曲线。从图中得到曲线斜率为0.300 010，与理论分析得到的有滑移条件下的螺母轴向位移随丝杠角位移变化率0.300 882 4 符合很好，误差只有0.29%。

3 结 论

本文首先根据实际应用中的小转角运动建立了行星滚柱丝杠的数学模型，把行星滚柱丝杠的轴向位移分解为啮合接触弹性变形量和传动轴向位移，用Hertz 接触理论求解了不同载荷下啮合接触弹性变形量，并推导了行星滚柱丝杠传动关系;在此基础上建立了行星滚柱丝杠的有限元模型，通过有限元仿真获得了行星滚柱丝杠的动力学响应。主要研究结论如下:

a. 在一定预紧力下，行星滚柱丝杠螺纹啮合处发生弹性变形，在行星滚柱丝杠精确控制系统的实际应用中应对这部分弹性变形量进行补偿。

b. 通过有限元方法获得了行星滚柱丝杠螺纹啮合点接触的应力分布及各部件的动力学响应曲线，结合理论推导的传动比，为行星滚柱丝杠在进给系统的实际应用提供了依据。

［1］Lemor P C. The Roller Screw，an Efficient and Reliable Mechanical Component of Electro－Mechanical Actuators［C］// Energy Conversion Engineering Conference，IEEE，1996，1:215－220.

［2］Andrade A，Nicolosi D，Lucchi J，et al. Auxiliary Total Artificial Heart:A Compact Electromechanical Artificial Heart Working Simultaneously with the Natural Heart［J］.Artif Organs，1999，23(9):876 –880.

［3］Richenbacher W E，Pae W E Jr，Magovern J A，et al.Roller Screw Electric Motor Ventricular Assist Device［J］.Transactions on American Society for Artificial Internal Organs，1986，32(1):46 –48.

［4］Rosenberg G，Snyder A，Weiss W，et al. A Roller Screw Drive for Implantable Blood Pumps［J］.Transactions on American Society for Artificial Internal Organs，1982，28:123 –126.

［5］Gromov V V，Misckevich A V，Yudkin E N，et al. The Mobile Penetrometer，a“Mole”for Sub－Surface Soil Investigation［C］//Proceedings of the Seventh European Space Mechanisms and Tribology Symposium，Noordwijk，the Netherlands，1997:151 –156.

［6］Wander J，Byrd V，Parker J .Initial Disturbance Accommodating Control System Analysis for Prototype Electromechanical Space Shuttle Steering Actuator［C］//Proceedings of the American Control Conference，Seattle W A，1995，6:3961 –3964.

［7］Marks S，Cortopassi C，Devries J，et al. The Advanced Light Source Elliptically Polarizing Undulator［C］//Proceedings of the Particle Accelerator Conference，Canada，1997，3:3221 –3223.

［8］Brandenburg G，Bruckl S，Dormann J，et al. Comparative Investigation of Rotary and Linear Motor Feed Drive Systems for High Precision Machine Tools［C］//Proceedings of the Sixth International Workshop on Advanced Motion Control，Nagoya，Japan，2000:384 –389.

［9］Dupont P E. Friction Modeling in Dynamic Robot Simulation［C］//Proceedings of the IEEE International Conference on Robotics and Automation，Cincinnati，1990，2:1370 –1376.

［10］Velinsky S A，Chu B，Lasky T A. Kinematics and Efficiency Analysis of the Planetary Roller Screw Mechanism［J］. Journal of Mechanical Design，2009，131(1):1－8.

［11］陈芳. 两级行星滚柱丝杠传动设计与精确度研究［D］.南京:南京理工大学，2009.

［12］韦振兴. 行星滚柱丝杠副刚度及效率的分析与研究［D］.武汉:华中科技大学，2011.

［13］Harris T A，Kotzalas M N. Essential Concepts of Bearing Technology［M］. CRC Press，1978:106－109.