考虑修形的斜齿轮系统非线性激励与动力学特性研究
魏静 王刚强 秦大同
摘要： 斜齿轮的啮合刚度与轮齿误差的求解是三维空间问题，其修形后的啮合刚度计算方法不同于直齿轮，而传统解析方法在计算斜齿轮啮合刚度时没有考虑斜齿轮啮合线和啮合位置的三维空间位置，无法准确得到修形后的斜齿轮系统啮合刚度激励与误差激励。建立综合考虑齿廓修形和齿向修形的刚度与误差非线性耦合激励模型，研究不同齿廓修形参数与齿向修形参数对斜齿轮啮合刚度以及系统动力学特性的影响规律；以系统振动加速度幅值最小为优化目标，确定斜齿轮系统的最佳修形值，利用数值方法得到斜齿轮系统的振动加速度幅频响应曲线，研究结果发现：选取的最佳修形参数可有效降低斜齿轮齿数交替区啮合刚度的波动，大幅度降低共振点附近的振动加速度幅值；最后通过建立的齿轮传动系统实验平台进行系统动力学特性实验研究，验证了理论模型及分析结果的正确性。
关键词： 斜齿轮； 齿廓修形； 齿向修形； 非线性激励； 动力学特性
中图分类号： TH132.4； O322文献标志码： A文章编号： 1004-4523（2018）04-0561-12
DOI：10.16385/j.cnki.issn.1004-4523.2018.04.003
引言
齿轮传动系统由于在齿数交替區载荷波动较大，会引起系统的振动和噪声[1]；另一方面由于齿轮受载后产生弯曲变形和扭转变形，造成轮齿沿齿宽方向接触不均匀，出现偏载现象[2-3]。而齿廓修形可以补偿齿轮实际啮合过程中的基节偏差，减缓啮合刚度波动，减小啮合冲击[4-5]；齿向修形能有效改善载荷沿轮齿接触线的不均匀分布，避免边缘接触，从而提高齿轮承载能力。
刚度激励计算方法主要有材料力学方法、近似代替法、切片法、有限元方法等。Chaari等[6]利用材料力学方法，考虑直齿轮弯曲变形、剪切变形、径向压缩变形、赫兹接触变形和轮体变形计算直齿轮的啮合刚度。Liu等[7]采用近似代替法，利用斜齿轮啮合线时变性代替啮合刚度的变化。由于近似代替法存在较大误差，Ajmi等[8]采用切片法，将斜齿轮沿齿宽方向离散成一片片薄片斜齿轮，然后根据变形协调方程求得每片薄片斜齿轮弹性变形和啮合刚度，并与有限元方法所得结果对比验证。
齿轮系统动力学模型按自由度数可以分为：纯扭转模型[9-10]、弯扭耦合动力学模型[11]、弯扭轴摆耦合动力学模型[12]。其中弯扭轴摆耦合动力学模型比其他模型考虑的自由度更多，因此更能反映齿轮系统的真实情况。在早期的齿轮动力学建模过程中，将误差激励和刚度激励分别代入动力学模型中进行齿轮动力学分析[13]。随着研究的深入，Chen等[14]基于材料力学方法提出了新型的综合啮合刚度与传递误差耦合的非线性激励解析计算模型，建立了考虑轮齿误差影响的直齿轮啮合刚度模型。王奇斌等[15]在Chen的基础上，拓展了该模型，建立了考虑齿向修形的直齿轮啮合刚度解析模型，并与有限元方法作对比，结果表明在不同修形量下两种方法求解的啮合刚度相对误差不超过5%。
现有文献多关注于单一修形方式对齿轮动力学特性的影响，且研究对象一般为直齿轮。斜齿轮的啮合刚度与轮齿误差是一个三维空间问题，其修形后啮合刚度的计算方法不同于直齿轮，而传统解析计算方法在计算斜齿轮啮合刚度时没有考虑斜齿轮啮合线和啮合位置的三维空间位置，无法计算修形后的斜齿轮啮合刚度与轮齿误差。
齿廓和齿向综合修形示意图如图3所示。仅考虑齿廓修形时，斜齿轮修形后端面示意图如图3（a）所示，虚线表示斜齿轮理论齿廓，实线表示斜齿轮实际齿廓，齿廓修形曲线为直线。修形后齿廓曲线上任一点修形量为Cax=Ca（xLa）（6）式中Ca表示最大齿廓修形量，La表示齿廓修形长度，x为该点到齿廓修形起始点的距离。
随着修形量的增大，斜齿轮啮合刚度值逐渐减小，传递误差随之增大；在齿数交替区域，啮合刚度波动先逐渐变小到定值后又随着修形量增加刚度波动逐渐变大。
先逐渐变平缓到趋于直线后又随着齿廓修形长度增加刚度曲线斜率逐渐变大。结合图4（c）和（c）可知齿廓修形量为30 μm，修形长度为6.4 mm时，啮合刚度均方差最小，即啮合刚度波动最小，与未修形相比刚度均方差减小幅度为95.02%。
2.2仅考虑齿向修形参数对斜齿轮刚度和传递误差的影响图6给出了齿向修形长度Lc=10 mm，齿向修形量Cc分别为0，5，10，15，20，25 μm时， 斜齿轮综合啮合刚度，传递误差以及刚度均方差。对比图4和图6可知，与齿廓修形相同的是随着齿向修形量的增大，斜齿轮啮合刚度值逐渐减小，传递误差随之增大；啮合刚度波动先逐渐变小到定值后又随着修形量增加刚度波动逐渐变大。不同的是虽然合适的齿向修形量也减小了啮合刚度的波动，但是对齿数交替区域啮合刚度的影响没有齿廓修形对齿数交替区域啮合刚度的影响明显。
图7给出了齿向修形量Cc=15 μm，齿向修形长度Lc分别为0，5，10，15，20，25 mm时，斜齿轮综合啮合刚度、传递误差以及刚度均方差。 对比图6和图7可知，除幅值稍有不同之外，齿向修形长度对啮合刚度与误差的影响与齿向修形量对啮合刚度与误差的影响变化趋势一致，结合图6（c）和图7（c）可知，当修形量为15 μm，修形长度为10 mm时，刚度均方差最小，即啮合刚度波动最小，与未修形相比刚度均方差减小幅度为94.587%。
2.3综合考虑齿廓和齿向修形参数对斜齿轮刚度和传递误差的影响以齿廓修形量Ca=30 μm，齿廓修形长度La=6.4 mm作为最佳齿廓修形值，在此基础上，研究不同齿向修形量和不同齿向修形长度对斜齿轮刚度和传递误差的影响。
图8给出了齿向修形长度Lc=5 mm，齿向修形量Cc分别为0，2.5，5.0，7.5，10.0，12.5 μm时，斜齿轮综合啮合刚度、传递误差以及刚度均方差。从图8（a）和（b）知：不同齿向修形量对少齿啮合区刚度和传递误差几乎没有影响，但对齿数交替区和多齿区刚度和传递误差影响较大；随着齿向修形量增加，啮合刚度波动先减小到定值后增大，而传递误差变化趋势与啮合刚度相反。
图9 给出了齿向修形量Cc=5 μm，齿向修形长度Lc分别为0，2.5，5.0，7.5，10.0，12.5 mm时，斜齿轮综合啮合刚度、传递误差以及啮合刚度均方差。对比图8和9可以看出，在同一齿廓修形参数下，除幅值不同之外，不同齿向修形长度对啮合刚度与传递误差的影响与不同修形量对啮合刚度与传递误差的影响变化趋势一致。结合图8（c）和图9（c）可知，以齿廓修形量Ca=30 μm，齿廓修形长度La=6.4 mm作为最佳齿廓修形值，齿向修形量为5.0 μm，修形长度为5 mm时，啮合刚度均方差最小，即啮合刚度波动最小，与未修形相比刚度均方差减小幅度为96.69%。
3考虑修形的斜齿轮传动系统动力学响应3.1斜齿轮传动系统动力学模型用有限元法[17]建立如下式所示的斜齿轮系统动力学方程M+C+KX=F（10）式中M为质量矩阵；X为广义坐标；C为阻尼矩阵，包括转子陀螺力矩；K为刚度矩阵；F为外载荷向量。
根据动力学方程建立12自由度平行轴系斜齿轮转子系统动力学模型，如图10所示，由大小齿轮、一级平行轴系、轴承组成。可将系统分为轴段单元、斜齿轮啮合单元以及轴承单元，为考虑轴段单元的剪切变形，用改进的Euler-Bernoulli梁作为轴段单元的理论模型，详细单元模型和矩阵形式可参考文献[17]，其中轴1共由7个轴段组成，一共8个节点，轴2共由10个轴段组成，共11个节点。
将求得的齿轮啮合刚度均值K12=5.98×108 N/m代入系统动力学模型中，耦合系统的无阻尼固有频率如表2所示。
齿廓修形长度La=6.4 mm，修缘量Ca分别为0，10，20，30，40，50 μm时，特殊节点主齿轮节点的x和θz方向随转速变化的振动加速度幅频曲线如图11所示。θz方向振动加速度可由θz方向振动角加速度变换可得。主动轮x方向振动加速度幅频曲线在啮合频率fN等于系统第14，16和21阶固有频率时出现共振峰；当啮合频率等于f16/2时，出现二次谐波共振峰。主动轮θz方向振动加速度幅频曲线在啮合频率fN等于系统第13，14和16阶固有频率时出现共振峰；当啮合频率等于f14/3，f14/2和f16/2时，出现高次谐波共振峰。随着齿廓修形量的增加，主动轮振动加速度幅值整体呈减小的趋势，修形量为30 μm时，振动加速度幅值减小幅度最大，此时系统的共振峰也减少了。修形量继续增大时，振动加速度幅值反而上升。与未修形相比，修形后系统的共振峰出现了偏移，而且随着修形量的增加，偏移趋势越趋于明显。
齿廓修形量Ca=30 μm，齿廓修形长度La分别为0，1.6，3.2，4.8，6.4，8.0 mm时，特殊节点主齿轮节点的x和θz方向随转速变化的振动加速度幅频曲线如图12所示。对比图11和12可以看出，在相同位置、不同齿修形量和修形长度下，主动轮振动加速度幅频曲线整体变化趋势一致，只是幅值大小略有所不同。修形长度为6.4 mm时，振动加速度幅值降低幅度最大。
3.3仅考虑齿向修形的系统动力学响应
齿向修形长度Lc=10 mm，齿向修形量Cc分别为0，5，10，15，20，25 μm时，特殊节点主齿轮节点的x和θz方向随转速变化的振动加速度幅频曲线如图13所示。不同齿向修形量下主动轮振动加速度幅值变化曲线与不同齿廓修形参数下得到的振动加速度幅值曲线变化规律大致相同。齿向修形量为15 μm时，振动加速度幅值减小幅度最大。
3.4综合齿廓修形和齿向修形的系统动力学响应
取齿廓修形量Ca=30 μm，齿廓修形长度La=6.4 mm作为最佳齿廓修形值， 在此基础上，研究不同齿向修形量和不同齿向修形长度对斜齿轮系统振动响应的影响规律。
通过图11和12可知，经过合适的齿廓修形，斜齿轮系统的动力学性能已经得到很大的改善，因此齿向修形值不宜过大。取齿向修形长度Lc=5 mm，齿向修形量Cc分别为0，2.5，5.0，7.5，10.0，12.5 μm时，特殊节点主齿轮节点的x和θz方向随转速变化的振动加速度幅频曲线如图15所示。综合齿向修形和齿廓修形与只进行齿廓修形相比，主动轮振动加速度幅频曲线变化趋势一致；随着齿向修形量的增加，振动加速度幅值整体降低，当修形量为5.0 μm时，振动加速度幅值最小，修形量继续增大时，振动加速度幅值反而上升。
齿向修形量Cc=5 μm，齿向修形长度Lc分别为0，2.5，5.0，7.5，10.0，12.5 mm时，特殊节点主齿轮节点的x和θz方向随转速变化的振动加速度幅频曲线如图16所示。对比图15和16可知，不同齿向修形量和齿向修形长度下主动轮振动加速度幅频曲线变化较为一致，当修形长度为5.0 mm时，振动加速度幅值最小。由此确定齿轮综合最佳修形值为：齿廓修形量Ca=30 μm，齿廓修形长度La=6.4 mm，齿向修形量Cc=5 μm，齿向修形长度Lc=5 mm。
4斜齿轮传动系统动力学特性实验研究
4.1实验台搭建和实验设计根据高铁牵引齿轮传动系统的运行工况，设计不同的实验工况，并进行稳态实验，从而得到不同负载、不同转速下的斜齿轮传动系统的振动特性和规律。高铁牵引齿轮传动系统实验平台搭建和布置如图17所示。由于驱动电机的额定输出转速远低于实验所需的高铁牵引齿轮箱的输入转速，所以增加陪试齿轮箱，作为增速器使用。
为减少信号在传递过程中的衰减，应在靠近振動源的地方布置振动加速度测点，如齿轮箱轴承座处。布置振动加速度测点在输入轴承座处和两边的输出轴承座处，如图18所示。每个测点有3个方向（横向、垂向、轴向），测点编号为1x，1y，1z，2x，2y，2z，3x，3y，3z。其中测点1在输入轴承座处，测点2和测点3为输出轴承座处。
分别进行两组实验：（1） 空载实验：除去负载激励的影响，此时可以忽略电机转速波动影响，系统的振动响应信号可认为仅由齿轮传动系统的内部激励造成的；（2） 稳定工况实验：获得不同实验工况下的稳态振动响应信号，研究系统的动力学特性。稳定工况实验又可分为转速平稳工况和扭矩平稳工况，高铁牵引齿轮箱分为正转（从输入轴方向看，逆时针为正转）、反转两个方向，小齿轮轴额定输入转速4100 r/min，额定输出转速1688 r/min，额定输入扭矩为1300 N·m，额定输出扭矩为3160 N·m。
4.2實验结果与理论分析对比验证
由于理论模型没有对箱体进行有限元建模，而实验得到的振动加速度有效值是从箱体表面所测的，需做如下处理：实验结果取所有测点的总振动加速度平均值，理论计算结果取轴承和齿轮上具有代表性的节点总的振动加速度平均值，然后通过实验结果来验证理论模型及分析结果的正确性。
转速平稳实验工况下，测点1，2和3的x，y和z各方向的振动加速度平均值随负载的变化如图19所示。将齿廓修形量Ca=30 μm，齿廓修形长度La=6.4 mm，齿向修形量Cc=5 μm，齿向修形长度Lc=5 mm作为齿轮最佳修形量代入动力学模型，得到斜齿轮系统转速平稳工况下的振动加速度幅频曲线，理论结果与实验结果对比如图20所示。随着负载的增加，理论振动加速度幅频曲线和实验测得结果振动加速度曲线变化趋势一致，均呈下降趋势。不同的是空载到额定负载区间，实验值大于理论值；超过额定负载时，理论值大于实验值。理论值与实验值振动加速度误差如表3所示。在额定转速和额定负载情况下，理论值与实验值振动加速度误差约为5.4%，验证了理论模型的正确性。
扭矩平稳实验工况下，测点1，2和3的z方向的振动加速度有效值随转速的变化如图21所示。在转速为3500 r/min时，振动加速度幅值曲线出现了共振峰。现对比空载工况下和扭矩平稳工况下理论与实验结果，从图22（a）可知，空载工况下，理论得到的振动加速度幅值是随着转速的上升而近似线性上升，而实验结果是：当转速小于3000 r/min时，振动加速度幅值呈近似线性上升趋势，当转速从3000 r/min至3500 r/min时，振动加速度幅值上升幅度增加，这与实验测得转速为3500 r/min时箱体向加速幅值曲线出现共振峰结果是一致的。由图22（b）可知，扭矩平稳工况下振动加速度幅值随转速增加而增加，在转速为2500和3500 r/min时出现了峰值；理论计算得到振动加速度幅值曲线整体趋势与实验结果基本吻合，在转速为2500 r/min时，理论振动加速度幅值曲线出现了峰值，但却没有在3500 r/min时出现共振峰。主要原因由于建模时没有考虑箱体等结构件对系统振动影响导致结果出现一定的偏差。
5结论
（1） 建立斜齿轮刚度与误差非线性激励耦合模型，研究不同齿廓修形量和齿廓修形长度以及不同齿向修形量和齿向修形长度对斜齿轮啮合刚度和传递误差的影响。选取合适的修形值可有效减小斜齿轮齿数交替区啮合刚度波动。
（2） 考虑系统的内部激励和外部激励，研究不同齿廓修形量和齿廓修形长度以及不同齿向修形量和齿向修形长度时系统的振动响应规律；随着修形值的增加，系统振动加速度幅值减小，共振峰减少；当修形量继续增加至一定数值时，系统振动加速度幅值上升，说明系统存在最佳修形量。
（3） 通过建立的齿轮传动系统实验平台进行动力学特性实验研究，得到不同工况下振动加速度幅频曲线并进行对比分析，验证了理论模型及分析结果的正确性。
参考文献：
[1]Ni Gaoxiang， Zhu Caichao， Song Chaosheng， et al. Tooth contact analysis of crossed beveloid gear transmission with parabolic modification[J].Mechanism and Machine Theory，2017，113（07）：40—52.
[2]International Organization for Standardization. ISO 6336-1-2006，Calculation of load capacity of spur and helical gears-part1：basic principles， introduction and general influence factors[S].Geneva，Switzerland，2006.
[3]Alfonso Fuentes-Aznar， Ignacio Gonzalez-Perez. Mathematical definition and computerized modeling of spherical involute and octoidal bevel gears generated by crown gear [J]. Mechanism and Machine Theory， 2016， 106（12）： 94—114.
[4]Liu Hui， Zhang Chen， Xiang Chang L， et al. Tooth profile relief based on lateral-torsional-rocking coupled nonlinear dynamic model of gear system[J]. Mechanism and Machine Theory， 2016， 105（11）：606—619.
[5]Lin Tengjiao， He Zeyin. Analytical method for coupled transmission error of helical gear system with machining errors， assembly errors and tooth modifications [J]. Mechanical Systems and Signal Processing， 2017， 91（07）： 167—182.
[6]Chaari Fakher， Fakhfakh Tahar， Haddar Mohamed. Analytical modelling of spur gear tooth crack and influence on gearmesh stiffness [J]. European Journal of Mechanics - A/Solids， 2009， 28（3）：461—468.
[7]Liu Changzhao， Qin Datong， Lim C， et al. Dynamic characteristics of the herringbone planetary gear set during the variable speed process [J]. Journal of Sound and Vibration， 2014， 333（24）： 6498—6515.
[8]Ajmi M， Velex P. A model for simulating the quasi-static and dynamic behaviour of solid wide-faced spur and helical gears [J].Journal of Sound and Vibration， 2005， 40（2）：173—190.
[9]Hidaka Teruaki， Terauchi Yoshio， Fujii Makoto. Analysis of dynamic tooth load on planetary gear [J]. Bulletin of JSME， 1980， 23（176）： 315—323.
[10]Kahraman Ahmet. Free torsional vibration characteristics of compound planetary gear sets [J]. Mechanism and Machine Theory， 2001， 36（8）： 953—971.
[11]Kahraman Ahmet. Load sharing characteristics of planetary transmissions [J]. Mechanism and Machine Theory， 1994， 29（8）： 1151—1165.
[12]魏靜， 孙清朝， 孙伟， 等. 高速机车牵引齿轮传动系统动态特性及非线性因素影响研究[J]. 振动与冲击， 2012， 31（17）： 38—43.
Wei Jing， Sun Qing-chao， Sun Wei， et al. Dynamic analysis and effects of nonlinear factors of a gear transmission system for high speed locomotive[J].Journal of Vibration Shock， 2012， 31（17）：38—43.
[13]Jia S， Howard I， Wang J. The dynamic modeling of multiple pairs of spur gears in mesh， including friction and geometrical errors [J]. International Journal of Rotating Machinery， 2003， 9：437—442.
[14]Chen Zaigang， Shao Yimin. Mesh stiffness calculation of a spur gear pair with tooth profile relief and tooth root crack[J]. Mechanism and Machine Theory， 2013， 62：63—74.
[15]王奇斌，张义民.齿向修形直齿轮系统动力学特性分析[J].振动工程学报，2016，29（1）：61—68.
Wang Qibin， Zhang Yimin.Dynamic characteristics analysis of a spur gear system with the lead crown relif[J]. Journal of Vibration Engineering， 2016， 29（1）：61—68.
[16]Lin J， Parker R G. Planetary gear parametric instability caused by mesh stiffness variation [J].Journal of Sound and Vibration， 2002， 249（1）： 129—145.
[17]魏静，张爱强，秦大同，等. 考虑结构柔性的行星轮系耦合振动特性研究[J]. 机械工程学报， 2017， 53（1）：1—12.
Wei Jing， Zhang Aiqiang， Qin Datong， et al. Coupling vibration analysis for planetary gear system considering flexible structure[J].Journal of Mechanical Engineering， 2017，53（1）：1—12.
Abstract： Meshing stiffness and tooth error of helical gears are three-dimensional space problem， the calculation method of mesh stiffness after modification are different with spur gears. However， the traditional analytical method in the calculation of the helical gear meshing stiffness without considering the three-dimensional space position of helical gear meshing line and the meshing position， stiffness excitation and error excitation of a helical gear system after modification cannot be accurately obtained. A nonlinear coupling excitation model of stiffness and error considering profile relief and crown relief is established.Based on the model，the influences of profile relief parameters and crown relief parameters on the meshing stiffness of the helical gear and the dynamic characteristics of the helical gear transmission system are studied respectively. Taking the minimum amplitude of the vibration acceleration as the optimization objective， the optimum modification values of the system are determined. The amplitude and frequency response curves of vibration acceleration of helical gear system are obtained by numerical method. The research results show that selecting the optimum modification parameters can effectively slow down the meshing stiffness fluctuation of the alternating area of teeth， and greatly reduce the amplitude of system vibration acceleration near the resonance point. Finally， the dynamic characteristics of the helical gear are researched with the experimental platform of high speed rail traction gearbox transmission system， and the correctness of the theoretical model is verified.
Key words： helical gear； profile relief； crown relief； nonlinear excitation； dynamic characteristics
摘要： 斜齿轮的啮合刚度与轮齿误差的求解是三维空间问题，其修形后的啮合刚度计算方法不同于直齿轮，而传统解析方法在计算斜齿轮啮合刚度时没有考虑斜齿轮啮合线和啮合位置的三维空间位置，无法准确得到修形后的斜齿轮系统啮合刚度激励与误差激励。建立综合考虑齿廓修形和齿向修形的刚度与误差非线性耦合激励模型，研究不同齿廓修形参数与齿向修形参数对斜齿轮啮合刚度以及系统动力学特性的影响规律；以系统振动加速度幅值最小为优化目标，确定斜齿轮系统的最佳修形值，利用数值方法得到斜齿轮系统的振动加速度幅频响应曲线，研究结果发现：选取的最佳修形参数可有效降低斜齿轮齿数交替区啮合刚度的波动，大幅度降低共振点附近的振动加速度幅值；最后通过建立的齿轮传动系统实验平台进行系统动力学特性实验研究，验证了理论模型及分析结果的正确性。
关键词： 斜齿轮； 齿廓修形； 齿向修形； 非线性激励； 动力学特性
中图分类号： TH132.4； O322文献标志码： A文章编号： 1004-4523（2018）04-0561-12
DOI：10.16385/j.cnki.issn.1004-4523.2018.04.003
引言
齿轮传动系统由于在齿数交替區载荷波动较大，会引起系统的振动和噪声[1]；另一方面由于齿轮受载后产生弯曲变形和扭转变形，造成轮齿沿齿宽方向接触不均匀，出现偏载现象[2-3]。而齿廓修形可以补偿齿轮实际啮合过程中的基节偏差，减缓啮合刚度波动，减小啮合冲击[4-5]；齿向修形能有效改善载荷沿轮齿接触线的不均匀分布，避免边缘接触，从而提高齿轮承载能力。
刚度激励计算方法主要有材料力学方法、近似代替法、切片法、有限元方法等。Chaari等[6]利用材料力学方法，考虑直齿轮弯曲变形、剪切变形、径向压缩变形、赫兹接触变形和轮体变形计算直齿轮的啮合刚度。Liu等[7]采用近似代替法，利用斜齿轮啮合线时变性代替啮合刚度的变化。由于近似代替法存在较大误差，Ajmi等[8]采用切片法，将斜齿轮沿齿宽方向离散成一片片薄片斜齿轮，然后根据变形协调方程求得每片薄片斜齿轮弹性变形和啮合刚度，并与有限元方法所得结果对比验证。
齿轮系统动力学模型按自由度数可以分为：纯扭转模型[9-10]、弯扭耦合动力学模型[11]、弯扭轴摆耦合动力学模型[12]。其中弯扭轴摆耦合动力学模型比其他模型考虑的自由度更多，因此更能反映齿轮系统的真实情况。在早期的齿轮动力学建模过程中，将误差激励和刚度激励分别代入动力学模型中进行齿轮动力学分析[13]。随着研究的深入，Chen等[14]基于材料力学方法提出了新型的综合啮合刚度与传递误差耦合的非线性激励解析计算模型，建立了考虑轮齿误差影响的直齿轮啮合刚度模型。王奇斌等[15]在Chen的基础上，拓展了该模型，建立了考虑齿向修形的直齿轮啮合刚度解析模型，并与有限元方法作对比，结果表明在不同修形量下两种方法求解的啮合刚度相对误差不超过5%。
现有文献多关注于单一修形方式对齿轮动力学特性的影响，且研究对象一般为直齿轮。斜齿轮的啮合刚度与轮齿误差是一个三维空间问题，其修形后啮合刚度的计算方法不同于直齿轮，而传统解析计算方法在计算斜齿轮啮合刚度时没有考虑斜齿轮啮合线和啮合位置的三维空间位置，无法计算修形后的斜齿轮啮合刚度与轮齿误差。
齿廓和齿向综合修形示意图如图3所示。仅考虑齿廓修形时，斜齿轮修形后端面示意图如图3（a）所示，虚线表示斜齿轮理论齿廓，实线表示斜齿轮实际齿廓，齿廓修形曲线为直线。修形后齿廓曲线上任一点修形量为Cax=Ca（xLa）（6）式中Ca表示最大齿廓修形量，La表示齿廓修形长度，x为该点到齿廓修形起始点的距离。
随着修形量的增大，斜齿轮啮合刚度值逐渐减小，传递误差随之增大；在齿数交替区域，啮合刚度波动先逐渐变小到定值后又随着修形量增加刚度波动逐渐变大。
先逐渐变平缓到趋于直线后又随着齿廓修形长度增加刚度曲线斜率逐渐变大。结合图4（c）和（c）可知齿廓修形量为30 μm，修形长度为6.4 mm时，啮合刚度均方差最小，即啮合刚度波动最小，与未修形相比刚度均方差减小幅度为95.02%。
2.2仅考虑齿向修形参数对斜齿轮刚度和传递误差的影响图6给出了齿向修形长度Lc=10 mm，齿向修形量Cc分别为0，5，10，15，20，25 μm时， 斜齿轮综合啮合刚度，传递误差以及刚度均方差。对比图4和图6可知，与齿廓修形相同的是随着齿向修形量的增大，斜齿轮啮合刚度值逐渐减小，传递误差随之增大；啮合刚度波动先逐渐变小到定值后又随着修形量增加刚度波动逐渐变大。不同的是虽然合适的齿向修形量也减小了啮合刚度的波动，但是对齿数交替区域啮合刚度的影响没有齿廓修形对齿数交替区域啮合刚度的影响明显。
图7给出了齿向修形量Cc=15 μm，齿向修形长度Lc分别为0，5，10，15，20，25 mm时，斜齿轮综合啮合刚度、传递误差以及刚度均方差。 对比图6和图7可知，除幅值稍有不同之外，齿向修形长度对啮合刚度与误差的影响与齿向修形量对啮合刚度与误差的影响变化趋势一致，结合图6（c）和图7（c）可知，当修形量为15 μm，修形长度为10 mm时，刚度均方差最小，即啮合刚度波动最小，与未修形相比刚度均方差减小幅度为94.587%。
2.3综合考虑齿廓和齿向修形参数对斜齿轮刚度和传递误差的影响以齿廓修形量Ca=30 μm，齿廓修形长度La=6.4 mm作为最佳齿廓修形值，在此基础上，研究不同齿向修形量和不同齿向修形长度对斜齿轮刚度和传递误差的影响。
图8给出了齿向修形长度Lc=5 mm，齿向修形量Cc分别为0，2.5，5.0，7.5，10.0，12.5 μm时，斜齿轮综合啮合刚度、传递误差以及刚度均方差。从图8（a）和（b）知：不同齿向修形量对少齿啮合区刚度和传递误差几乎没有影响，但对齿数交替区和多齿区刚度和传递误差影响较大；随着齿向修形量增加，啮合刚度波动先减小到定值后增大，而传递误差变化趋势与啮合刚度相反。
图9 给出了齿向修形量Cc=5 μm，齿向修形长度Lc分别为0，2.5，5.0，7.5，10.0，12.5 mm时，斜齿轮综合啮合刚度、传递误差以及啮合刚度均方差。对比图8和9可以看出，在同一齿廓修形参数下，除幅值不同之外，不同齿向修形长度对啮合刚度与传递误差的影响与不同修形量对啮合刚度与传递误差的影响变化趋势一致。结合图8（c）和图9（c）可知，以齿廓修形量Ca=30 μm，齿廓修形长度La=6.4 mm作为最佳齿廓修形值，齿向修形量为5.0 μm，修形长度为5 mm时，啮合刚度均方差最小，即啮合刚度波动最小，与未修形相比刚度均方差减小幅度为96.69%。
3考虑修形的斜齿轮传动系统动力学响应3.1斜齿轮传动系统动力学模型用有限元法[17]建立如下式所示的斜齿轮系统动力学方程M+C+KX=F（10）式中M为质量矩阵；X为广义坐标；C为阻尼矩阵，包括转子陀螺力矩；K为刚度矩阵；F为外载荷向量。
根据动力学方程建立12自由度平行轴系斜齿轮转子系统动力学模型，如图10所示，由大小齿轮、一级平行轴系、轴承组成。可将系统分为轴段单元、斜齿轮啮合单元以及轴承单元，为考虑轴段单元的剪切变形，用改进的Euler-Bernoulli梁作为轴段单元的理论模型，详细单元模型和矩阵形式可参考文献[17]，其中轴1共由7个轴段组成，一共8个节点，轴2共由10个轴段组成，共11个节点。
将求得的齿轮啮合刚度均值K12=5.98×108 N/m代入系统动力学模型中，耦合系统的无阻尼固有频率如表2所示。
齿廓修形长度La=6.4 mm，修缘量Ca分别为0，10，20，30，40，50 μm时，特殊节点主齿轮节点的x和θz方向随转速变化的振动加速度幅频曲线如图11所示。θz方向振动加速度可由θz方向振动角加速度变换可得。主动轮x方向振动加速度幅频曲线在啮合频率fN等于系统第14，16和21阶固有频率时出现共振峰；当啮合频率等于f16/2时，出现二次谐波共振峰。主动轮θz方向振动加速度幅频曲线在啮合频率fN等于系统第13，14和16阶固有频率时出现共振峰；当啮合频率等于f14/3，f14/2和f16/2时，出现高次谐波共振峰。随着齿廓修形量的增加，主动轮振动加速度幅值整体呈减小的趋势，修形量为30 μm时，振动加速度幅值减小幅度最大，此时系统的共振峰也减少了。修形量继续增大时，振动加速度幅值反而上升。与未修形相比，修形后系统的共振峰出现了偏移，而且随着修形量的增加，偏移趋势越趋于明显。
齿廓修形量Ca=30 μm，齿廓修形长度La分别为0，1.6，3.2，4.8，6.4，8.0 mm时，特殊节点主齿轮节点的x和θz方向随转速变化的振动加速度幅频曲线如图12所示。对比图11和12可以看出，在相同位置、不同齿修形量和修形长度下，主动轮振动加速度幅频曲线整体变化趋势一致，只是幅值大小略有所不同。修形长度为6.4 mm时，振动加速度幅值降低幅度最大。
3.3仅考虑齿向修形的系统动力学响应
齿向修形长度Lc=10 mm，齿向修形量Cc分别为0，5，10，15，20，25 μm时，特殊节点主齿轮节点的x和θz方向随转速变化的振动加速度幅频曲线如图13所示。不同齿向修形量下主动轮振动加速度幅值变化曲线与不同齿廓修形参数下得到的振动加速度幅值曲线变化规律大致相同。齿向修形量为15 μm时，振动加速度幅值减小幅度最大。
3.4综合齿廓修形和齿向修形的系统动力学响应
取齿廓修形量Ca=30 μm，齿廓修形长度La=6.4 mm作为最佳齿廓修形值， 在此基础上，研究不同齿向修形量和不同齿向修形长度对斜齿轮系统振动响应的影响规律。
通过图11和12可知，经过合适的齿廓修形，斜齿轮系统的动力学性能已经得到很大的改善，因此齿向修形值不宜过大。取齿向修形长度Lc=5 mm，齿向修形量Cc分别为0，2.5，5.0，7.5，10.0，12.5 μm时，特殊节点主齿轮节点的x和θz方向随转速变化的振动加速度幅频曲线如图15所示。综合齿向修形和齿廓修形与只进行齿廓修形相比，主动轮振动加速度幅频曲线变化趋势一致；随着齿向修形量的增加，振动加速度幅值整体降低，当修形量为5.0 μm时，振动加速度幅值最小，修形量继续增大时，振动加速度幅值反而上升。
齿向修形量Cc=5 μm，齿向修形长度Lc分别为0，2.5，5.0，7.5，10.0，12.5 mm时，特殊节点主齿轮节点的x和θz方向随转速变化的振动加速度幅频曲线如图16所示。对比图15和16可知，不同齿向修形量和齿向修形长度下主动轮振动加速度幅频曲线变化较为一致，当修形长度为5.0 mm时，振动加速度幅值最小。由此确定齿轮综合最佳修形值为：齿廓修形量Ca=30 μm，齿廓修形长度La=6.4 mm，齿向修形量Cc=5 μm，齿向修形长度Lc=5 mm。
4斜齿轮传动系统动力学特性实验研究
4.1实验台搭建和实验设计根据高铁牵引齿轮传动系统的运行工况，设计不同的实验工况，并进行稳态实验，从而得到不同负载、不同转速下的斜齿轮传动系统的振动特性和规律。高铁牵引齿轮传动系统实验平台搭建和布置如图17所示。由于驱动电机的额定输出转速远低于实验所需的高铁牵引齿轮箱的输入转速，所以增加陪试齿轮箱，作为增速器使用。
为减少信号在传递过程中的衰减，应在靠近振動源的地方布置振动加速度测点，如齿轮箱轴承座处。布置振动加速度测点在输入轴承座处和两边的输出轴承座处，如图18所示。每个测点有3个方向（横向、垂向、轴向），测点编号为1x，1y，1z，2x，2y，2z，3x，3y，3z。其中测点1在输入轴承座处，测点2和测点3为输出轴承座处。
分别进行两组实验：（1） 空载实验：除去负载激励的影响，此时可以忽略电机转速波动影响，系统的振动响应信号可认为仅由齿轮传动系统的内部激励造成的；（2） 稳定工况实验：获得不同实验工况下的稳态振动响应信号，研究系统的动力学特性。稳定工况实验又可分为转速平稳工况和扭矩平稳工况，高铁牵引齿轮箱分为正转（从输入轴方向看，逆时针为正转）、反转两个方向，小齿轮轴额定输入转速4100 r/min，额定输出转速1688 r/min，额定输入扭矩为1300 N·m，额定输出扭矩为3160 N·m。
4.2實验结果与理论分析对比验证
由于理论模型没有对箱体进行有限元建模，而实验得到的振动加速度有效值是从箱体表面所测的，需做如下处理：实验结果取所有测点的总振动加速度平均值，理论计算结果取轴承和齿轮上具有代表性的节点总的振动加速度平均值，然后通过实验结果来验证理论模型及分析结果的正确性。
转速平稳实验工况下，测点1，2和3的x，y和z各方向的振动加速度平均值随负载的变化如图19所示。将齿廓修形量Ca=30 μm，齿廓修形长度La=6.4 mm，齿向修形量Cc=5 μm，齿向修形长度Lc=5 mm作为齿轮最佳修形量代入动力学模型，得到斜齿轮系统转速平稳工况下的振动加速度幅频曲线，理论结果与实验结果对比如图20所示。随着负载的增加，理论振动加速度幅频曲线和实验测得结果振动加速度曲线变化趋势一致，均呈下降趋势。不同的是空载到额定负载区间，实验值大于理论值；超过额定负载时，理论值大于实验值。理论值与实验值振动加速度误差如表3所示。在额定转速和额定负载情况下，理论值与实验值振动加速度误差约为5.4%，验证了理论模型的正确性。
扭矩平稳实验工况下，测点1，2和3的z方向的振动加速度有效值随转速的变化如图21所示。在转速为3500 r/min时，振动加速度幅值曲线出现了共振峰。现对比空载工况下和扭矩平稳工况下理论与实验结果，从图22（a）可知，空载工况下，理论得到的振动加速度幅值是随着转速的上升而近似线性上升，而实验结果是：当转速小于3000 r/min时，振动加速度幅值呈近似线性上升趋势，当转速从3000 r/min至3500 r/min时，振动加速度幅值上升幅度增加，这与实验测得转速为3500 r/min时箱体向加速幅值曲线出现共振峰结果是一致的。由图22（b）可知，扭矩平稳工况下振动加速度幅值随转速增加而增加，在转速为2500和3500 r/min时出现了峰值；理论计算得到振动加速度幅值曲线整体趋势与实验结果基本吻合，在转速为2500 r/min时，理论振动加速度幅值曲线出现了峰值，但却没有在3500 r/min时出现共振峰。主要原因由于建模时没有考虑箱体等结构件对系统振动影响导致结果出现一定的偏差。
5结论
（1） 建立斜齿轮刚度与误差非线性激励耦合模型，研究不同齿廓修形量和齿廓修形长度以及不同齿向修形量和齿向修形长度对斜齿轮啮合刚度和传递误差的影响。选取合适的修形值可有效减小斜齿轮齿数交替区啮合刚度波动。
（2） 考虑系统的内部激励和外部激励，研究不同齿廓修形量和齿廓修形长度以及不同齿向修形量和齿向修形长度时系统的振动响应规律；随着修形值的增加，系统振动加速度幅值减小，共振峰减少；当修形量继续增加至一定数值时，系统振动加速度幅值上升，说明系统存在最佳修形量。
（3） 通过建立的齿轮传动系统实验平台进行动力学特性实验研究，得到不同工况下振动加速度幅频曲线并进行对比分析，验证了理论模型及分析结果的正确性。
参考文献：
[1]Ni Gaoxiang， Zhu Caichao， Song Chaosheng， et al. Tooth contact analysis of crossed beveloid gear transmission with parabolic modification[J].Mechanism and Machine Theory，2017，113（07）：40—52.
[2]International Organization for Standardization. ISO 6336-1-2006，Calculation of load capacity of spur and helical gears-part1：basic principles， introduction and general influence factors[S].Geneva，Switzerland，2006.
[3]Alfonso Fuentes-Aznar， Ignacio Gonzalez-Perez. Mathematical definition and computerized modeling of spherical involute and octoidal bevel gears generated by crown gear [J]. Mechanism and Machine Theory， 2016， 106（12）： 94—114.
[4]Liu Hui， Zhang Chen， Xiang Chang L， et al. Tooth profile relief based on lateral-torsional-rocking coupled nonlinear dynamic model of gear system[J]. Mechanism and Machine Theory， 2016， 105（11）：606—619.
[5]Lin Tengjiao， He Zeyin. Analytical method for coupled transmission error of helical gear system with machining errors， assembly errors and tooth modifications [J]. Mechanical Systems and Signal Processing， 2017， 91（07）： 167—182.
[6]Chaari Fakher， Fakhfakh Tahar， Haddar Mohamed. Analytical modelling of spur gear tooth crack and influence on gearmesh stiffness [J]. European Journal of Mechanics - A/Solids， 2009， 28（3）：461—468.
[7]Liu Changzhao， Qin Datong， Lim C， et al. Dynamic characteristics of the herringbone planetary gear set during the variable speed process [J]. Journal of Sound and Vibration， 2014， 333（24）： 6498—6515.
[8]Ajmi M， Velex P. A model for simulating the quasi-static and dynamic behaviour of solid wide-faced spur and helical gears [J].Journal of Sound and Vibration， 2005， 40（2）：173—190.
[9]Hidaka Teruaki， Terauchi Yoshio， Fujii Makoto. Analysis of dynamic tooth load on planetary gear [J]. Bulletin of JSME， 1980， 23（176）： 315—323.
[10]Kahraman Ahmet. Free torsional vibration characteristics of compound planetary gear sets [J]. Mechanism and Machine Theory， 2001， 36（8）： 953—971.
[11]Kahraman Ahmet. Load sharing characteristics of planetary transmissions [J]. Mechanism and Machine Theory， 1994， 29（8）： 1151—1165.
[12]魏靜， 孙清朝， 孙伟， 等. 高速机车牵引齿轮传动系统动态特性及非线性因素影响研究[J]. 振动与冲击， 2012， 31（17）： 38—43.
Wei Jing， Sun Qing-chao， Sun Wei， et al. Dynamic analysis and effects of nonlinear factors of a gear transmission system for high speed locomotive[J].Journal of Vibration Shock， 2012， 31（17）：38—43.
[13]Jia S， Howard I， Wang J. The dynamic modeling of multiple pairs of spur gears in mesh， including friction and geometrical errors [J]. International Journal of Rotating Machinery， 2003， 9：437—442.
[14]Chen Zaigang， Shao Yimin. Mesh stiffness calculation of a spur gear pair with tooth profile relief and tooth root crack[J]. Mechanism and Machine Theory， 2013， 62：63—74.
[15]王奇斌，张义民.齿向修形直齿轮系统动力学特性分析[J].振动工程学报，2016，29（1）：61—68.
Wang Qibin， Zhang Yimin.Dynamic characteristics analysis of a spur gear system with the lead crown relif[J]. Journal of Vibration Engineering， 2016， 29（1）：61—68.
[16]Lin J， Parker R G. Planetary gear parametric instability caused by mesh stiffness variation [J].Journal of Sound and Vibration， 2002， 249（1）： 129—145.
[17]魏静，张爱强，秦大同，等. 考虑结构柔性的行星轮系耦合振动特性研究[J]. 机械工程学报， 2017， 53（1）：1—12.
Wei Jing， Zhang Aiqiang， Qin Datong， et al. Coupling vibration analysis for planetary gear system considering flexible structure[J].Journal of Mechanical Engineering， 2017，53（1）：1—12.
Abstract： Meshing stiffness and tooth error of helical gears are three-dimensional space problem， the calculation method of mesh stiffness after modification are different with spur gears. However， the traditional analytical method in the calculation of the helical gear meshing stiffness without considering the three-dimensional space position of helical gear meshing line and the meshing position， stiffness excitation and error excitation of a helical gear system after modification cannot be accurately obtained. A nonlinear coupling excitation model of stiffness and error considering profile relief and crown relief is established.Based on the model，the influences of profile relief parameters and crown relief parameters on the meshing stiffness of the helical gear and the dynamic characteristics of the helical gear transmission system are studied respectively. Taking the minimum amplitude of the vibration acceleration as the optimization objective， the optimum modification values of the system are determined. The amplitude and frequency response curves of vibration acceleration of helical gear system are obtained by numerical method. The research results show that selecting the optimum modification parameters can effectively slow down the meshing stiffness fluctuation of the alternating area of teeth， and greatly reduce the amplitude of system vibration acceleration near the resonance point. Finally， the dynamic characteristics of the helical gear are researched with the experimental platform of high speed rail traction gearbox transmission system， and the correctness of the theoretical model is verified.
Key words： helical gear； profile relief； crown relief； nonlinear excitation； dynamic characteristics
