面齿轮传动啮合刚度数值计算

　　0、 引言

　　作为机械装置中的一个重要零部件,齿轮传动被广泛应用于航空、风电、汽车等领域。随着工作转速的逐步提高,齿轮传动的动力学性能越来越受到设计、制造及使用者的重视。相比较其他类型的传动系统,齿轮传动系统的主要不同之处在于: 它不但会因为外部激励而产生动态响应,同时会因为传动过程中啮合齿对数的改变、轮齿的弹性变形及轮齿误差而导致啮合刚度发生变化,从而产生轮齿动态啮合力,且此种由于啮合综合刚度的时变性引起的动态激励是齿轮传动中最主要的动态激励形式之一。因此,确定齿轮传动的时变啮合刚度一直是齿轮动力学研究中的重要问题。

　　对于圆柱齿轮的啮合刚度计算问题,已有较多的学者进行过深入研究,得到了一些比较成熟的计算方法。对于面齿轮及螺旋锥齿轮等结构相对复杂的齿轮传动形式的啮合刚度计算问题,由于其齿面为复杂曲面,要准确计算时变啮合刚度存在较大难度,因此相关研究文献较少。

　　Gosselin 等基于有限条法给出了计算螺旋锥齿轮刚度的方法,但也只得到了沿齿高和齿长方向的位移曲线,并没有得到刚度曲线。Mennem等使用有限单元法计算了不同载荷下轮齿接触柔度,获得了时变刚度。

　　面齿轮传动是近二十年才真正发展起来的一种主要用于航空领域的齿轮传动形式,针对其啮合刚度的计算方法,目前尚未见文献报道。在面齿轮的动力学研究中,齿轮啮合刚度都采用经过傅里叶变换后的一次正余弦函数来近似,与实际啮合刚度存在比较大的差异。有学者研究发现,采用近似时变啮合刚度得到的系统动态因子比采用实际啮合刚度时大,而且不能有效地得到系统在低频阶段的动态响应。因此,采用近似时变啮合刚度来拟合实际啮合刚度研究齿轮传动系统动力学行为,并不能正确反映齿轮系统的动态特性。
 

　　本文提出了一种新的齿轮传动时变啮合刚度数值计算方法。首先以直齿圆柱齿轮为例,建立合理的有限元模型,得到直齿圆柱齿轮时变啮合刚度曲线,并与 ISO6336方法计算结果进行对比,说明该啮合刚度计算方法的正确性及有限元模型的精确性; 然后应用该啮合刚度计算方法,研究面齿轮传动时变啮合刚度变化规律,得到精确的面齿轮传动时变啮合刚度曲线。研究结果为面齿轮传动的动力学分析及设计提供参考。

　　1 齿轮传动啮合刚度计算方法

　　在齿轮运转过程中,轮齿之间的啮合作用可以简化为一个在啮合线方向上的时变弹簧,弹簧刚度等于该时刻的啮合刚度,如图 1 所示。图中,ω 为主动轮转速,Rbp、Rbg分别为主从动轮基圆半径,Km为啮合刚度,T为负载扭矩,LOA(line of action) 表示啮合线。

　　齿轮啮合刚度的一般表达式为Kn= Fn/ δn(1)其中,Fn为法向接触力;δn为法向综合弹性变形,它应包括由于局部赫兹接触产生的接触弹性变形量 δH、轮齿弯曲产生的接触位置的位移 δF、轮毂变形产生的接触位置的位移 δR、轴和支撑结构变形对接触点位置的影响量 δA。由于考虑轴、轴承及支撑结构变形会使问题过于复杂,故本文暂只考虑 δH、δF及 δR,因此可以得到

　　在齿轮传动过程中,由于轮齿受载变形及齿轮误差,从动轮转动位置会偏离其不考虑受载变形及齿轮误差时所应处的理想转动位置。齿轮传递误差eT在1958 年被Harris定义为“从动轮实际转动位移与理想转动位移之差”,如图 2 所示,其一般表达式为

　　式中,θp、θg分别为主从动轮转角; zp、zg分别为主从动轮齿数。

　　假设齿轮没有承受载荷,此时轮齿没有弹性变形,但是由于轮齿误差的存在,依旧会产生无负载传递误差 eNLT。在齿轮传动的加载接触分析中,由于存在几何模型误差、有限元网格划分产生的几何误差以及网格划分产生的尺侧间隙,故无负载传递误差是不可避免的。无负载传递误差导致齿轮啮合过程中轮齿产生刚性位移,该位移在计算轮齿受载变形时必须减除。同时,为了方便理解,将传递误差转换成啮合线方向上的位移,因此可以得到齿轮副在一定负载下的受载传递误差eLT一般表达式

　　由此可以得到轮齿啮合过程中综合弹性变形

　　2、 啮合刚度数值计算方法验证

　　由于结构的复杂性,面齿轮传动啮合刚度计算方法目前尚未见相关文献报道。为了说明上述啮合刚度计算方法的正确性及有限元模型的精确性,本文先以直齿圆柱齿轮为例进行计算,并将其与 ISO6336: 2006 方法的计算结果进行对比。

　　2. 1 有限元建模

　　基于唐进元等提出的齿轮虚拟加工制造方法,采用表 1 所示齿轮传动参数,得到高精度的直齿圆柱齿轮三维几何装配模型。

　　

　　为了节约计算成本,忽略轴及轴承等支撑结构的影响。鉴于静力学分析中质量分布不会影响分析结果,因此在模型中仅考虑可能产生啮合的5 对轮齿,其他部分简化为一个轮缘圆柱,并分别在大小轮内孔面与其对应旋转轴线上的参考点之间建立刚性耦合约束,在大轮耦合节点施加负载扭矩,在小轮耦合节点施加转动位移,以此来模拟小轮通过与大轮啮合抵消大轮反抗扭矩带动大轮运转的实际工况,边界条件如图 3 所示。对可能发生接触的 5 对轮齿进行局部网格细化,齿厚方向最小网格尺寸为 0. 037 mm,齿高方向最小网格尺寸为 0. 142 mm,总节点数为180 462,总单元数为158 592。网格模型如图 4 所示。

　　将大轮及小轮的工作齿面定义为接触对,暂不考虑摩擦对啮合刚度的影响。采用隐式静力学分析算法,输出接触力及大小轮的绕自身旋转轴线的旋转位移历程曲线。

　　2. 2 直齿圆柱齿轮啮合刚度计算

　　理论上,标准渐开线齿轮的静态传递误差为零,但考虑到有限元模型不可避免地存在网格划分带来的几何误差,且有限元分析时必须要有完整的边界约束条件,因此,本文近似以大轮负载为13. 2 N·m 时的传递误差作为无负载传递误差,该负载与额定负载 660 N·m 相比较小,产生的轮齿接触变形及弯曲变形可以忽略不计。基于上述齿轮传动啮合刚度计算方法,计算得到该齿轮副时变啮合刚度曲线,并与 ISO6336: 2006 给出的计算公式得到的结果进行对比,如图 5 所示。

　　由图 5 可知,所提出的方法得到的时变啮合刚度曲线与 ISO6336: 2006 方法计算结果相比,平均值相对误差为 5. 86%。因此可以证明所提出的齿轮传动时变啮合刚度计算方法的正确性,同时说明有限元模型已足够精确。

　　3、 面齿轮传动啮合刚度数值计算

　　3. 1 有限元模型

　　面齿轮传动几何参数如表 2 所示,其中面齿轮轮缘系数定义为面齿轮轮齿高度与面齿轮齿根平面到面齿轮齿底平面之间距离的比值。

　　

　　采用上述有限元建模方法,将可能发生接触的 5 对轮齿之外的部分简化为一个阶梯圆柱,忽略轴及轴承等支撑结构变形对面齿轮传动接触位置的影响,在小轮内孔面与位于其旋转轴线上的参考节点之间建立刚性耦合约束,在面齿轮与位于其旋转轴线上的参考节点之间建立刚性耦合约束,即将轴及轴承简化为刚性支撑,在小轮参考节点施加旋转位移,在面齿轮参考节点施加负载扭矩,以模拟小轮通过与大轮啮合作用带动大轮抵消负载扭矩而运转的实际工况,边界条件设置如图 6 所示。

　

　　对于未修形的直齿面齿轮传动而言,面齿轮由齿顶进入啮合,由齿根退出啮合,其理论接触轨迹在齿高方向大致呈一条竖直线。经齿轮接触分析得到,该对齿轮理论接触轨迹大致位于半径为164. 4 mm 处,因此对面齿轮靠近内端的一半轮齿及配对小轮对应啮合位置进行网格细化,其网格密度与上述直齿圆柱齿轮模型中细化部分网格密度基本一致,总节点数为242 312,总单元数为201 420,网格模型如图 7 所示。接触定义、材料参数及分析步设置与上述直齿圆柱齿轮有限元模型保持一致。此时对于单个面齿轮轮齿而言,约有 54 个计算点。

　　

　　3. 2 面齿轮传动啮合刚度数值计算

　　采用上述有限元模型,在大轮额定载荷为1536. 8 N·m 条件下,取负载为 38. 42 N·m 时的传递误差为无负载传递误差。额定载荷作用下各齿对法向接触力历程曲线如图 8 所示。由图 8 可知,该对齿轮副存在双齿接触区(DT) 及三齿接触区(TT) ,并且,由于在三齿接触区时前一对轮齿即将退出啮合,后一对轮齿刚刚进入啮合,齿间载荷分布非常不均匀,因此该对轮齿法向接触力最大值出现在该啮合区域。

　　

　　进一步根据式(6) 计算得到面齿轮传动时变啮合刚度曲线,如图 9 所示。由图 9 可以看出,面齿轮传动时变啮合刚度曲线类似于一条正弦曲线,呈现明显的周期特征,其最小值大致出现在三齿接触区内前后两对轮齿受载相同的位置(图 8中点 A 处) 。此时该两对轮齿约各承担总载荷的6. 44% ,中间一对轮齿约承担总载荷的 87. 15% ,齿间载荷的分布非常不均匀,并且此时前后两对轮齿分别在靠近齿根及齿顶位置接触,因此该啮合位置齿轮副弹性综合变形最大,啮合刚度最小。

　　啮合刚度最大值出现在双齿啮合区内参与啮合的两对轮齿受载相同的位置(图 8 中点 B 处) ,此时该两对轮齿各承担总载荷的 50%,载荷分布最为均匀,因此该啮合位置齿轮弹性综合变形最小,啮合刚度最大。

　　图 9 时变啮合刚度曲线(T =1536. 8 N·m)

　　4、 结论

　　(1) 提出了一种新的齿轮传动时变啮合刚度数值计算方法,并给出了相应的加载接触分析有限元建模方法。

　　(2) 以直齿圆柱齿轮为例,采用该方法得到该齿轮副时变啮合刚度曲线,与 ISO6336 方法计算结果相对比,验证了上述方法的正确性及所建有限元模型的精确性。

　　(3) 应用该啮合刚度数值计算方法得到了面齿轮传动载荷分布规律及时变啮合刚度曲线。