1 引言
液压马达是一种将流体压力能转化为旋转机械能的动力设备，广泛应用于兵器装备领域[1-2]。相对于柱塞式、叶片式、摆缸式等类型的液压马达，非圆行星齿轮液压马达在体积、质量、抗油液污染等方面具有显著优势，其核心是一个无系杆、变中心距的非圆行星齿轮机构[3]。非圆行星齿轮机构的特殊结构在设计过程中受到诸多约束条件的限制，齿廓的生成是设计难点之一[4]。
目前关于非圆齿轮齿廓设计的方法主要有2种，数值计算方法和范成仿真方法。数值计算方法是以共轭啮合原理或运动关系为基础计算齿廓坐标的方法，郑方焱等[5]利用节曲线的法向量，提出了一种非圆齿轮齿廓的通用算法；刘永平等[6]计算了椭圆齿轮齿廓法线长度从而得出齿廓点坐标数据；牛子孺等[7]构建变比齿廓数字求解模型，将数据点云拟合后获得变比齿条的齿廓；Li Botao等[8]提出了一种基于Jarvis-March的非圆齿轮齿廓数值计算方法。范成仿真法是从非圆齿轮的加工原理出发，模拟刀具与齿轮毛坯之间的范成运动，再通过布尔运算求得非圆齿轮齿廓实体的方法，俞高红等[9]针对共轭凹凸型非圆齿轮，在Pro/E软件中以主动轮做刀具、从动轮做毛坯生成了从动轮齿廓；丁国龙等[10]在SolidWorks中模拟了齿条刀具对齿扇毛坯的加工过程，实现了非圆齿扇的齿廓设计；Hu Jian等[11]在CATIA中开发宏程序，进行了切削过程的仿真，获得变速比齿轮副的齿廓。
数值计算方法的精度高，但是设计算法的难度较大，在工程实际中不易推广；范成仿真法需要对三维软件进行二次开发，程序编制复杂，且会对齿坯形成切除痕迹，导致齿廓的精确度不理想，但是范成仿真法生成齿廓的过程直观简捷，若能直接采用合理的方法避免程序编制，同时对齿廓曲线进行光顺，就能够突破范成仿真法的局限性。本研究在分析非圆齿轮范成加工的基础上，考虑磨削加工提高轮齿精度的原理，综合应用计算机图形学中包络图形绘制以及图像平滑技术生成齿廓，旨在能够高效精准地生成非圆齿轮行星轮系，为非圆齿轮液压马达的研究提供技术支撑。
2 非圆齿轮制造过程分析
2.1 范成加工分析
范成法加工齿轮的基本原理是保持刀具和齿轮毛坯之间按无侧隙啮合的运动关系来进行切齿的[12-13]，插齿切削是常见的范成加工方法，针对插齿切削的过程进行分析。
非圆齿轮的插齿切削过程是由切削主运动、径向进给运动、分齿展成运动合成得到的，如图1(a)所示，插齿刀在非圆齿轮毛坯的轴线方向上做直线往复运动，辅助让刀运动，形成切削主运动；在切削主运动过程中，为了达到全齿深度，在径向方向，插齿刀需做进给运动；非圆齿轮毛坯每自转一个角度，则插齿刀必须回转相应的角度，产生分齿展成运动。非圆齿轮的齿廓就是由插齿刀多次切削后所形成的，如图1(b)所示。插齿加工过程中可以简化为插齿刀节圆与非圆齿轮节曲线进行纯滚动，在纯滚动过程中，插齿刀在非圆齿轮毛坯上插制出待加工的齿廓形状，这个过程与计算机图形学中的包络极为相似。
图1 非圆齿轮插齿加工示意图
Fig.1 Schematic diagram of non-circular pinion machining
根据齐次坐标变换方法分析包络过程。非圆齿轮液压马达中的行星齿轮机构，是由太阳轮、内齿圈和行星轮构成，其中太阳轮、内齿圈为非圆齿轮。如图2(a)所示，以太阳轮毛坯中心为原点建立空间固定坐标系S-Oxy，以插齿刀中心为原点建立坐标系Sd-Odxdyd、Sp-Opxpyp，其中Sp-Opxpyp与插齿刀固连，在插齿刀自转的过程中，Sp-Opxpyp随之旋转，Sd-Odxdyd是随动坐标系，其xd轴、yd轴一直与S-Oxy的x轴、y轴平行。插齿刀沿太阳轮节曲线的逆时针方向滚动，在初始位置1处，Sd-Odxdyd、Sp-Opxpyp完全重合，当插齿刀由初始位置1滚动至位置2时，Sp-Opxpyp相对于Sd-Odxdyd转过了φ角。插齿刀齿廓是由建立在Sp-Opxpyp下的方程决定的，将Sp-Opxpyp下的坐标旋转φ角可得到Sd-Odxdyd下的坐标，表示为
(1)
将Sd-Odxdyd下的坐标平移可得到S-Oxy下的坐标，假设O1点在S-Oxy中的坐标为(xd,yd)，表示为
(2)
则插齿刀齿廓在Sp-Opxpyp下的坐标(xp,yp)经过2次变换矩阵可转换为S-Oxy下的坐标(x,y)，表示为
(3)
每前进一个步进角φ，计算一次插齿刀齿廓在S-Oxy下的坐标，将每一个计算步的图形叠加在一起形成非圆齿轮齿廓的包络图形，步进角φ越小，则包络的精度越高。
根据式(3)可知，获得包络图需要明晰Od点在S-Oxy中的坐标(xd,yd)以及插齿刀齿廓在Sp-Opxpyp下的坐标 (xp,yp)。动点Od形成了插齿刀中心轨迹，插齿刀中心轨迹是太阳轮节曲线的法向等距线，获得太阳轮的节曲线就可得到(xd,yd)，对于插齿刀齿廓坐标(xp,yp)而言，构建齿廓方程即可获得。对于内齿圈而言，包络图形的获取途径与太阳轮相同，如图2(b)所示。
图2 非圆齿轮齿廓包络线形成原理
Fig.2 Forming principle of tooth profile envelope of non-circular gear
2.2 磨削加工分析
插齿加工完成后，为了提高齿轮精度和齿面粗糙度，可以对齿形进行精加工，利用砂轮磨削轮齿是一种代表性的精加工方法。非圆齿轮的磨削可以采用多轴联动的方式进行，如图3所示，砂轮绕着自身的轴线高速旋转，并沿着齿宽方向往复移动，形成磨削的主运动；非圆齿轮在水平面上的两个方向平移，保持齿廓与砂轮的相切，完成展成运动；磨削完一个齿槽后，非圆齿轮旋转一定的角度，将砂轮切换到下一个齿槽中，实现分度运动。非圆齿轮的磨削过程复杂，若采用齐次坐标变换方法模拟磨削过程来优化齿廓的方法效率较低，可以更换思路，从磨削原理出发来探索齿廓优化的方法。
图3 非圆齿轮磨削示意图
Fig.3 Schematic diagram of grinding non-circular gear
磨削是砂轮磨粒对齿面挤压、刻划、切削、摩擦抛光的综合作用，从垂直于轴线的平面上来看，磨削的本质是将范成加工在齿廓上形成的大量微线段进行了平滑处理。在计算机图形学中，非均匀有理B样条(NURBS)是曲线造型的数学基础，其表达式为：
(4)
式中：Bi,k(u)为k次B样条基函数；Wi为权重系数；V1为控制顶点。
从包络图形中提取的非圆齿轮齿廓存在多段微小的NURBS曲线，造成了齿廓曲率的不连续，在给定误差的条件下，通过导入或创建的型值点来插值曲线，把分段曲线连接成光顺相切的曲线，从而提高齿廓精度，改善传动质量。
3 三维数字化建模过程
3.1 包络图形的绘制
3.1.1 节曲线的数值计算
太阳轮节曲线边数是n1，内齿圈的节曲线边数是n3，一般采用的组合有2-3型、2-4型、3-4型、3-5型、4-5型、4-6型、5-6型、5-7型等[14]，以4-6型为例进行计算，即太阳轮节曲线的边数是4、内齿圈的节曲线边数是6。太阳轮齿数z1、行星轮齿数z2、内齿圈齿数z3之间的关系为：
(5)
10≤z2<0.5(z3-z1)
(6)
计算得到的非圆行星齿轮机构基本参数如表1所示。
表1 非圆行星齿轮机构基本参数
Table 1 Basic parameters of non-circular planetary gear mechanism
参数太阳轮行星轮内齿圈节曲线边数n1=4n2=1n3=6齿数z1=44z2=10z3=66模数/mmm1=2m2=2m3=2齿宽/mmb1=10b2=10b3=10
非圆行星轮系的基本运动关系如图4所示，r1、r2、r3分别是太阳轮节曲线、行星轮节圆、内齿圈节曲线对应的向径，θ1、θ3分别是太阳轮、内齿圈所对应的极角，ω1、ω2、ω3分别是太阳轮、行星轮、内齿圈所对应的瞬时角速度，e、f分别是行星轮节圆与太阳轮节曲线、内齿圈节曲线的切点，μ1是m点处太阳轮节曲线切向的正方向与向径r1的夹角，μ2是m点处内齿圈节曲线切向的正方向与向径r3的夹角。在非圆行星轮系当中，行星轮与太阳轮啮合时的节圆，和行星轮与内齿圈啮合时的节圆可以是不同的，若如图4所示的2个节圆是同一个圆时，μ1、μ2相等。
图4 基本运动关系
Fig.4 Basic motion relationship
太阳轮的节曲线是高阶椭圆齿轮，其节曲线方程为：
(7)
式中：k为偏心率；θ1为极角；A为长轴半径。太阳轮与内齿圈的节曲线封闭条件为：
(8)
由微分几何可知：
(9)
太阳轮与内齿圈的轮齿均布条件为：
(10)
根据非圆行星轮系的基本运动关系，内齿圈的节曲线方程为：
(11)
非圆行星齿轮机构在运行过程中不产生运动干涉的条件为：
r3min≥r1max+2ha
(12)
联立方程(7)～式(10)计算太阳轮节曲线数据，再通过方程式(11)求出内齿圈节曲线数据，代入方程式(12)中进行校验，得到节曲线的数值解，将节曲线的数值导入UG NX中，采用三次样条曲线进行拟合，得到的太阳轮、内齿圈节曲线如图5所示。
图5 太阳轮、内齿圈节曲线
Fig.5 Pitch curve of a sun wheel and its annulargear
3.1.2 插齿刀齿廓的确定
渐开线齿形是非圆齿轮最常见的齿形，但是渐开线齿廓在齿根处采用过渡圆角，受载时弯曲应力较大，为了提高非圆齿轮的啮合传动性能，近年来出现了很多新型齿廓形式。正弦曲线的上半部分为凸曲线，下半部分为凹曲线，且具有一定程度的可分性，符合齿条轮廓的基本要求[15]，如图6(a)所示。以正弦曲线作为齿廓的齿条刀具生成的齿轮具有较高的弯曲强度，本研究中尝试采用正弦曲线齿廓的插齿刀形成非圆齿轮的实体模型，由于内齿圈无法采用齿条型刀具进行加工，因此将正弦曲线环形折弯后变成圆柱型插齿刀具，如图6(b)所示。
图6 插齿刀齿廓
Fig.6 Profile of a gear shaper
在坐标系S3-O3x3y3中绘制正弦线，其方程为
y3=hsin(bx3)
(12)
式中:h为正弦线的振幅，也即插齿刀的齿顶(根)高，mm，取值为2.5；b为常数，取值为1。
将横轴x3绕成一个圆，形成插齿刀节圆，则插齿刀的一个齿距是一个完整的正弦线周期，在坐标系S1-O1x1y1下，插齿刀的齿廓方程为：
(13)
式中：m0为插齿刀的模数，mm，取值为2；z0为插齿刀的齿数，取值为10；t为变量，变化范围为0～1。
3.1.3 包络过程的实现
依据计算机图形学中的包络原理，在UG NX软件平台中进行范成过程仿真，主要涉及建模、装配、运动仿真3个模块，以太阳轮为例进行说明。
在建模模块中，将太阳轮的节曲线偏移一个齿顶高的距离形成齿顶曲线并进行拉伸，形成太阳轮齿坯实体；对于插齿刀齿廓图形同样进行拉伸操作，形成插齿刀实体。
在装配模块中，调入太阳轮齿坯实体与插齿刀实体，选择合理的切削起始位置，将两者进行装配，特别需要注意的是插齿刀节圆需要与太阳轮节曲线相切。此外，还需要在装配模块中用草图方式创建一个辅助连杆，便于在范成过程中控制插齿刀的运动。
在运动仿真模块中，将装配好的模型按照定义运动体—建立运动副—创设约束—添加追踪—计算求解的步骤进行。建立运动副与创设约束是其中的关键步骤，不仅需要在辅助连杆上施加旋转副，而且需要在插齿刀节圆与太阳轮节曲线设置线在线上副，在插齿刀节圆圆心与辅助连杆之间设置点在线上副。在计算过程中，系统追踪每一个计算步得到的插齿刀位置，并进行备份，形成了内齿圈齿廓的包络图形。综合考虑齿廓精度与计算机的运算能力，设置计算步为300，内齿圈的范成仿真计算结果如图7所示。
图7 太阳轮的包络图
Fig.7 Envelope diagram of sun wheel
内齿圈的范成过程与太阳轮相同，仿真计算结果如图8所示。
图8 内齿圈的包络图
Fig.8 Envelope diagram of annulargear
3.2 齿廓曲线的光顺
依据计算机图形学中的图像平滑原理，进行齿廓曲线的光顺。利用提升体命令生成的内齿圈及太阳轮实体模型上存在范成过程仿真形成的切削痕迹，在UG NX软件平台上，无法直接对齿廓曲面进行优化，因此将太阳轮、内齿圈进行正投影，得到齿廓曲线，通过对齿廓曲线的优化来达到的目的，正投影效果如图9所示。
图9 正投影效果
Fig.9 Positive projection effect
采用光顺曲线串的命令来处理齿廓曲线，光顺曲线串是光顺样条+连结曲线+曲率圆角的组合体命令，通过对角度阈值、距离阈值的控制，使NURBS曲线二阶连续。对齿廓曲线光顺后，得到优化后的齿廓曲线，如图10所示，与光顺之前的齿廓曲线进行对比，可以明显看出，范成后形成的齿廓曲线上存在局部曲率的骤增或骤降现象，表面形貌上有明显的接茬，而优化后的齿廓曲线避免了曲率的急剧过渡和突然变化，将优化后的齿廓曲线拉伸后得到的曲面比较光滑，表面质量得到显著改善。
图10 齿廓光顺效果
Fig.10 Smoothing effect of gear tooth profile
4 仿真验证
行星轮的齿廓与插齿刀相同，只是齿全高存在变化，将插齿刀的齿顶去除一部分即可实现行星轮的三维模型创建。将内齿圈、太阳轮、行星轮进行装配，形成非圆齿轮行星轮系。
在UG NX软件平台中设置行星轮与太阳轮之间、行星轮与内齿圈之间的运动关系为3D接触，令内齿圈静止不动，太阳轮的角速度为120(°)/s，模拟非圆齿轮行星轮系的运行过程。行星轮与太阳轮、内齿圈组成密闭腔容积每变化一个周期，太阳轮的转动角为：
(14)
在当前的结构条件和仿真条件下，太阳轮的转动角为150°，密闭腔容积的变化周期为1.25 s。以行星轮1、2为观察对象，如图11所示，可以明显看出行星轮1、2与太阳轮、内齿圈组成封闭区域面积不断发生变化，对封闭区域面积进行测量，测量结果如图12所示，0 s时封闭区域面积为320.65 mm2，在运行过程中封闭区域面积从大变小，再从小到大，到达顶峰后又开始回落，1.25 s时封闭区域面积又回到了0 s时的初始状态，以此循环往复，这种现象与非圆齿轮液压马达的工作过程是相同的。
图11 非圆齿轮行星轮系运行过程
Fig.11 Operating process of planetary gear train with non-circular gears
图12 封闭区域面积
Fig.12 Area of enclosed area
太阳轮旋转1周，封闭区域面积的变化次数为：
(15)
为了分析整周期条件下的行星轮运动规律，令太阳轮旋转2.5周，即设置仿真时间为7.5 s，在仿真运行过程中，行星轮与太阳轮及内齿圈的轮齿逐一啮合，没有干涉现象的发生，输出角速度曲线进行分析，如图13所示，行星轮的角速度呈现出明显的周期性变化，根据基本运动关系，行星轮角速度ω2与太阳轮角速度ω1之间存在关系式：
图13 行星轮与太阳轮的角速度
Fig.13 Angular velocity of planetary and solar wheels
(16)
式中，角度α1、α′的几何意义在图4中标出。
将行星轮的角速度曲线与式(16)进行对比，其变化趋势与理论分析是一致的，但是行星轮角速度曲线存在一定程度的波动，造成这种现象的原因主要有2个，一是轮齿啮合的动态激励，二是非圆齿轮的变速运动特性形成的振动脉冲，存在角速度的波动属于齿轮传动过程中的正常现象。
5 结论
1) 以非圆齿轮的范成加工过程为理论依据，利用UG NX软件平台的运动仿真功能实现太阳轮、内齿圈的插齿包络图形绘制，在此基础上根据磨削加工原理，采用计算机图形学中的图像平滑方法光顺曲线串可以消除插齿切削痕迹，提高齿廓质量。
2) 4-6型非圆行星齿轮机构在太阳轮齿数为44、行星轮齿数为10、内齿圈齿数为66、齿轮模数为2 mm的结构条件下，令太阳轮的角速度为120(°)/s进行仿真计算，计算结果表明2个行星轮与太阳轮、内齿圈之间的封闭区域面积在0 s时为320.65 mm2，在运行过程中封闭区域面积从大变小，再从小到大，1.25 s时又回到了初始状态，以此循环往复，且行星轮角速度呈现出明显的周期性变化，与非圆齿轮液压马达的运行情况相符。
3) 使用运动仿真与光顺曲线串的功能生成齿廓的过程直观、结果准确，更换节曲线数据及插齿刀齿廓方程后，即可生成不同规格、不同类型的非圆齿轮行星轮系，可快速实现系列化、变异化设计，为非圆齿轮液压马达的设计提供一种新的思路与实施途径。
[1] 张伟杰，汪国胜，郭勇.机-电-液馈能悬挂系统性能分析[J].兵器装备工程学报,2021,42(12):103-111.
Zhang W J,Wang G S,Guo Y.Performance analysis of mechanical-electro-hydraulic regenerative suspension system[J].Journal of Ordnance Equipment Engineering,2021,42(12):103-111.
[2] Ding H.Application of non-circular planetary gear mechanism in the gear pump[J].Advanced Materials Research,2012,591-593:2139-2142.
[3] 郭守存，李殿起，尚志强，等.基于MATLAB的非圆行星齿轮结构设计[J].重型机械,2016(04):49-52.
Guo S C,Li D Q,Shang Z Q,et al.Structural design of non-circular planetary gear based on Matlab[J].Heavy Machinery,2016(04):49-52.
[4] 李渤涛，陈定方.非圆齿轮设计、制造、检测及应用[J].机械工程学报,2020,56(09):55-72.
Li B T,Chen D F.Design,manufacture,inspection and application of non-circular gears[J].Journal of Mechanical Engineering,2020,56(09):55-72.
[5] 郑方焱，陈定方，刘有源.一种非圆齿轮齿廓的通用算法[J].机械传动,2013,37(04):64-66.
Zheng F Y,Chen D F,Liu Y Y.A generic algorithm of non-circular gear tooth profile[J].Journal of Mechanical Transmission,2013,37(04):64-66.
[6] 刘永平，王鹏，李佳.解析法设计椭圆齿轮齿廓[J].机械设计与制造,2015(02):214-216.
Liu Y P,Wang P,Li J.Design the elliptical gear tooth profile via analytical method[J].Machinery Design &Manufacture,2015(02):214-216.
[7] 牛子孺，辛振波，赵新学，等.交错轴变比齿条变比齿廓数字建模[J].计算机集成制造系统,2018,24(10):2514-2524.
Niu Z R,Xin Z B,Zhao X X,et al.Digital modeling of rack variable ratio tooth profile of rack and pinion steering gear[J].Computer Integrated Manufacturing Systems,2018,24(10):2514-2524.
[8] Li B T,Hu J Q,Chen D F,et al.Numerical algorithm of non-circular gear’s tooth profile based on jarvis march[C]//International Conference on Human Centered Computing.Springer,Cham,2016.
[9] 俞高红，虞佳萍，童俊华，等.一种共轭凹凸型非圆齿轮机构的设计[J].中国机械工程,2016,27(16):2155-2159.
Yu G H,Yu J P,Tong J H,et al.Design of a conjugate concave-convex non-circular gear mechanism[J].China Mechanical Engineering,2016,27(16):2155-2159.
[10] 丁国龙，沈观成，周红浪，等.汽车转向器摇臂轴非圆齿扇三维参数化建模方法[J].机械科学与技术,2020,39(08):1242-1247.
Ding G L,Shen G C,Zhou H L,et al.Three-dimensional parametric modeling of non-circular gear sector of antomobile steering arm[J].Mechanical Science and Technology for Aerospace Engineering,2020,39(08):1242-1247.
[11] Jian H U,Dong K,Niu Z,et al.Digital design method of variable ratio gear pair for automobile mechanical variable ratio steering gear based on enveloping simulation[J].Journal of the Balkan Tribological Association,2016,22(01):184-199.
[12] Lin C,Xia X G,Li P L.Geometric design and kinematics analysis of coplanar double internal meshing non-circular planetary gear train[J].Advances in Mechanical Engineering,2018,10(12):1-12.
[13] Han J,Li D Z,Tian X Q,et al.Linkage model and interpolation analysis of helical non-circular gear hobbing[J].Journal of the Brazilian Society of Mechanical Sciences and Engineering,2020,42(11):1-13.
[14] 姚文席.非圆齿轮设计[M].北京：机械工业出版社,2012.
Yao W X.Noncircular gear design[M].Beijing：China Machine Press,2012.
[15] 崔高尚，范元勋.基于正弦基准齿形的齿轮共轭齿廓设计与分析[J].机械传动,2015,39(08):59-62.
Cui G S,Fan Y X.Design and analysis of gear conjugate profile based on the sine-curved benchmark tooth profile[J].Journal of Mechanical Transmission,2015,39(08):59-62.