王 燕， 刘建新， 李 淼， 李奕璠， 蔡久凤(. 西南交通大学 机械工程学院， 成都 600; . 西南交通大学 牵引动力国家重点实验室， 成都 600;. 中车株洲电力机车研究所有限公司， 湖南 株洲 400)

齿轮啮合过程中，由于内部啮合刚度和传递误差的动态变化以及齿侧间隙的存在，导致系统表现出强烈的非线性动力学特性。该特性得到了国内外学者的高度关注,并开展了大量的研究[1-6]。而齿轮作为机车驱动系统的关键部件之一，其动力学特性不仅影响着机车驱动系统的振动性能，也影响着机车的动力传输和走行部的运行安全[7]。但是，目前针对计及齿轮啮合特性的机车驱动系统的振动研究工作仍较为少见。此外，机车服役过程中，由于轮轨黏着的独特性，使得机车驱动系统的外部激励具有一定的特殊性。

基于此，本文建立了轮轨黏着下计及齿轮啮合静传递误差、时变啮合刚度和齿侧间隙在内的机车驱动系统的动力学模型及方程，采用多尺度法进行求解，获得了系统在内部齿轮啮合动态激励和外部轮轨黏着力共同作用下的系统主共振的频率响应方程。并开展了实例研究，分析了系统参数变化对主共振频率响应曲线的影响，仿真了不同蠕滑率下轮轨黏着力的动态变化对机车驱动系统主共振响应的影响。可为轮轨黏着下计及齿轮啮合特性的机车驱动系统主共振问题提供参考依据。

1 机车驱动系统的动力学分析模型

机车驱动系统的动力学分析模型如图1所示：其中，符号Ii、Ri、θi、Ti分别为主动(i=1)、被动(i=2)齿轮端上的等效转动惯量、基圆半径、转动角位移和施加的力矩；km(t)、cm分别为齿轮啮合刚度与阻尼；b为齿侧间隙的一半；e(t)为齿轮静传递误差。

1.1 动力学方程

根据图1所示的动力学分析模型，建立相应的动力学方程

km(t)R1f(Δx)=T1，

km(t)R2f(Δx)=-T2

(1)

式中：Δx=R1θ1-R2θ2-e(t)；f(Δx)为齿侧间隙函数。

整理式(1)后得到

(2)

式中：me=I1I2/(I1R22+I2R12)；fm=I1R2R/(I1R22+I2R12)；F0=I2R1T1/(I1R22+I2R12)；R为车轮滚动圆半径；Fu为轮轨黏着力。

式(2)包含了系统在力矩平衡状态下的相对位移Δx0。因此，设静平衡位置处的轮轨黏着力为Fu0，并令x=Δx-Δx0，整理式(2)后，得到

(3)

式中：静平衡位置处的静态载荷F1=F0+fmFu0。

1.2 齿轮啮合的内部激励

1.2.1 齿轮时变啮合刚度

本文采用势能原理求解齿轮时变啮合刚度，因此，先将齿轮时变啮合刚度表示为关于主动齿轮转动角度θ1(单位°)的傅里叶级数的形式

(4)

式中:k0为平均啮合刚度；kj为各阶谐波项的幅值；φj为相位角；wn=2π/T，T=360/Z1，Z1为主动齿轮齿数。

进而通过变换得到km(t)关于时间的傅里叶级数形式

(5)

式中:we=2πZ1n1/60，n1为主动齿轮转速,r/min。

1.2.2 齿轮静传递误差

齿轮静传递误差e(t)采用余弦函数表示为[8]

e(t)=e0+ercos(wet+φ)

(6)

式中:e0、er分别为齿轮静传递误差的常数项和谐波项的幅值；φ为相位角。

1.2.3 齿侧间隙

选定齿侧间隙的一半b为特征尺寸，将齿侧间隙函数f(Δx)采用3次多项式表示为

f(Δx)/b=a1(Δx/b)+a3(Δx/b)3

(7)

式中:a1=0.167;a3=0.064。

1.3 轮轨黏着力

采用文献[9]中简化的轮轨黏着特性曲线，如图2所示。

图2中，um为最大黏着系数，sm为最大黏着系数对应的蠕滑率；ku为黏着特性曲线在AB段的斜率。OA段时，黏着系数处于黏着段；AB段时，黏着系数处于滑动段。

图2 简化的轮轨黏着特性曲线Fig.2 Simplified characteristic curve of wheel-rail adhesion

轮轨黏着系数u可表示为

(8)

式中:蠕滑率s是与轮对纵向速度v有关的函数，用图1中的符号可表示为

(9)

参考文献[10]中的处理方法，将蠕滑率进行泰勒级数展开，取其一次项后表示为

(10)

式中：w为平衡状态下被动齿轮的转动角速度；静平衡位置处的蠕滑率s0=1-v/(wR)。

轮轨黏着力Fu=uQ，其中Q为轴重(单位：N)。轮轨黏着力的动态变化值Fu-Fu0可表示为

1) 静平衡位置处的蠕滑率s0在OA段时

(11)

2) 静平衡位置处的蠕滑率s0在AB段时

(12)

2 方程的无量纲化

将式(11)、式(12)代入式(3)，整理好后得到

(13)

式中：c1，F2，F3的值如表1所示。其中，c1为轮轨黏着力变化引起的等效阻尼；F2为静态载荷的波动值；F3为误差速度项激励的系数。

表1 c1、F2、F3的值Tab.1 The values of c1、F2、F3

令x=by，Δx0=bΔy0(b为特征尺寸)，将式(5)～式(7)代入式(13)，整理后得到无量纲化后的方程为

a3(y+Δy0)3]=f1+f2+f3ϖsin(ϖτ+φ)+

f4ϖ2cos(ϖτ+φ)

(14)

3 系统主共振的频响方程

cos(ϖτ+φ-θ)

(15)

令ti=εiτ(i=0,1)，按照多尺度法[11-12]，设式(15)的近似解为

y=y0(t0,t1)+εy1(t0,t1)

(16)

将式(16)代入式(15)，令等式两端ε0和ε1的系数分别相等，得到

(17)

ε1次

(18)

式(17)的解为

(19)

将式(19)代入式(18)，整理可得：

(20)

式中：cc为前面各项的共轭复数。

当ϖ≈w1时，系统产生主共振。引入一个解谐参数σ=O(1)来表示激励频率ϖ与系统固有频率w1的接近程度。令ϖ=w1+εσ，消除式(20)中的久期项后得到

(21)

设A=a(t1)eiβ(t1)/2并代入式(21)，当ke2与ke1相比较小时，分离实部与虚部，整理后得到

φ-θ)

(22a)

(22b)

式中：γ=σt1-β。

当a′=γ′=0时，系统存在稳态运动。由式(22)可得到稳态响应的频率响应方程为

(23)

4 仿真分析

在某型机车驱动系统中，齿轮副的主要计算参数为：主动齿轮齿数Z1=32，被动齿轮齿数Z2=75，模数m=12，压力角a=22.5°，静传递误差的常数项、谐波项分别为e0=0、er=0.108 mm。轮轨黏着特性曲线中，选取um=0.25，sm=0.02，ku=-0.25。

基于势能原理[13-15]，综合考虑轮齿弯曲变形、剪切变形、轴向压缩变形、赫兹接触变形、以及轮体变形对应的啮合线上的等效刚度。进而获得齿轮综合时变啮合刚度，并将其展开为傅里叶级数形式，结果如图3所示。其中幅值kj、相位角φj(j=0～4)的值如表2所示。

表2 啮合刚度各阶谐波幅值与相位角Tab.2 Harmonic amplitude and phase angle of mesh stiffness

当φ1-φ=π时，忽略θ的微小变化，仿真得到不同f2、f3、f4、ke1、μ下，系统主共振的频率响应曲线如图4-图8所示。

由图4～图8可知，随着无量纲参数f2的增大，频响曲线呈现向右平移的趋势，且幅值有所变化。由此可见，轮轨黏着力的动态变化值影响着系统的主共振发生的频率和振幅。同时，系统的非线性特性使得主共振的频响曲线产生多值情况，且无量纲参数f3、f4、ke1的减小或μ的增大有利于减小多值区域的产生。此外，f3、f4、ke1的减小对振幅有一定的抑制作用。所以，误差速度项激励的系数F3、误差谐波项的幅值er、一次谐波刚度的比值ke1的减小不仅有利于减小多值区域的产生，也对振幅有一定的抑制作用。阻尼的增大也有利于减小系统多值区域的产生。

图4 不同f2时的频响曲线Fig.4 Frequency response curve with different f2

图5 不同f3时的频响曲线Fig.5 Frequency response curve with different f3

图6 不同f4时的频响曲线Fig.6 Frequency response curve with different f4

图7 不同ke1时的频响曲线Fig.7 Frequency response curve with different ke1

图8 不同μ时的频响曲线Fig.8 Frequency response curve with different μ

设系统的初始条件为：y(0)=dy(0)/dτ=0。当ϖ=w1，不计齿轮啮合阻尼cm，静平衡位置处的蠕滑率s0分别为0.012、0.035，得到考虑轮轨黏着力的动态变化时机车驱动系统主共振时的时域响应、频谱图、相平面图和轮轨黏着系数的变化曲线分别如图9、图11所示。当静平衡位置处的蠕滑率s0分别为0.012、0.035，轮轨黏着力为恒值时，得到机车驱动系统主共振的时域响应、频谱图和相平面图分别如图10、图12所示。

(9a)

(9b)

(9c)

(9d)图9 s0=0.012，轮轨黏着力变化时Fig.9 The results when s0=0.012 considerating the change of wheel-rail adhesion force

对比图9和图10可知，当s0=0.012时，静平衡位置处的黏着系数处于黏着段。在无量纲时间[2 000～2 500]范围内,考虑轮轨黏着力的动态变化时主共振位移的最大值较轮轨黏着力为恒值时减小了约35.35%。且考虑轮轨黏着力的动态变化与否对系统振动的频谱图和相轨迹具有巨大影响。

对比图11和图12可知，当s0=0.035时，静平衡位置处的黏着系数处于滑动段。在无量纲时间[2 000～2 500]范围内,考虑轮轨黏着力的动态变化时主共振位移的最大值较轮轨黏着力为恒值时增大了约115.55%。且考虑轮轨黏着力的动态变化与否时，两者的频谱图和相轨迹存在较大差异。

(10a)

(10b)

(10c)图10 s0=0.012，轮轨黏着力为恒值Fig.10 The results when s0=0.012 and wheel-rail adhesion force is a constant

(11a)

(11b)

(11c)

(11d)图11 s0=0.035，轮轨黏着力变化时Fig.11 The results when s0=0.035 considerating the change of wheel-rail adhesion force

(12a)

(12b)

(12c)图12 s0=0.035，轮轨黏着力为恒值Fig.12 The results when s0=0.035 and wheel-rail adhesion force is a constant

5 结 论

本文考虑了机车驱动系统内部齿轮啮合的静传递误差、时变啮合刚度和齿侧间隙，建立了内部齿轮啮合动态激励和外部轮轨黏着力激励共同作用下系统的动力学模型及方程，采用解析法和数值仿真相结合的方法，研究了机车驱动系统主共振时的振动特性，通过分析得出以下结论：

(1) 采用解析法—多尺度法进行求解，获得了轮轨黏着下计及齿轮啮合特性的机车驱动系统主共振时的频率响应方程，可为系统主共振的频响特性提供理论依据。

(2) 通过数值仿真系统参数变化对系统主共振频响曲线的影响得知：轮轨黏着力的动态变化值会影响系统主共振产生的频率和幅值；系统的非线性特性使得主共振的频响曲线产生多值情况，且误差速度项激励的系数F3、误差谐波项的幅值er、一次谐波刚度的比值ke1的减小或阻尼的增大有利于减小多值区域的产生，同时，F3、er、ke1的减小对振幅有一定的抑制作用。

(3) 通过数值仿真不同黏着系数s0下系统主共振的振动响应(不计齿轮啮合阻尼时)得知：在无量纲时间[2 000～2 500]范围内，当s0=0.012，考虑轮轨黏着力的动态变化时位移的最大值较轮轨黏着力为恒值时减小了约35.35%；当s0=0.035，考虑轮轨黏着力的动态变化时位移的最大值较轮轨黏着力为恒值时增大了约115.55%；且两种情况下，考虑轮轨黏着力的动态变化与否对系统振动的频谱图和相轨迹具有较大影响。

参 考 文 献

[1] AL-SHYYAB A, KAHRAMAN A.Non-linear dynamic analysis of a multi-mesh gear train using multi-term harmonic balance method: sub-harmonic motions[J].Journal of Sound and Vibration,2005, 279(1):417-451.

[2] SHEN Yongjun,YANG Shaopu,LIU Xiandong.Nonlinear dynamics of a spur gear pair with time-varying stiffness and backlash based on incremental harmonic balance method[J]. International Journal of Mechanical Sciences,2006,48(11): 1256-1263.

[3] 魏静,孙伟,褚衍顺,等.斜齿轮系统分岔与混沌特性及其参数影响研究[J].哈尔滨工程大学学报,2013,34(10):1301-1309.

WEI Jing,SUN Wei,CHU Yanshun,et al.Bifurcation and chaotic characteristics of helical gear system and parameter influences[J]. Journal of Harbin Engineering University,2013,34(10):1301-1309.

[4] 李应刚,陈天宁,王小鹏,等.外部动态激励作用下齿轮系统非线性动力学特性[J].西安交通大学学报,2014,48(1):101-105.

LI Yinggang,CHEN Tianning,WANG Xiaopeng,et al.Non-linear dynamics of spur gear pair under external periodic excitation[J]. Journal of Xi’an Jiaotong University,2014,48(1):101-105.

[5] FARSHIDIANFAR A,SAGHAFI A.Global bifurcation and chaos analysis in nonlinear vibration of spur gear systems[J]. Nonlinear Dynamics,2014,75(4):783-806.

[6] GOU Xiangfeng,ZHU Lingyun,CHEN Dailin.Bifurcation and chaos analysis of spur gear pair in two-parameter plane[J]. Nonlinear Dynamics,2015,79(3):2225-2235.

[7] 黄冠华,张卫华,宋纾崎,等.高速列车齿轮传动系统谐振分析[J].交通运输工程学报,2014,14(6):51-58.

HUANG Guanhua,ZHANG Weihua,SONG Shuqi,et al. Harmonic resonance analysis of gear transmission system for high speed train[J].Journal of Traffic and Transportation Engineering, 2014, 14(6): 51-58.

[8] 王晓笋,巫世晶,周旭辉,等.含侧隙非线性齿轮传动系统的分岔与混沌分析[J].振动与冲击,2008,27(1):53-56.

WANG Xiaosun,WU Shijing,ZHOU Xuhui,et al.Bifurcation and chaos in a nonlinear dynamic model of spur gear with backlash[J]. Journal of Vibration and Shock,2008,27(1):53-56.

[9] YAO Yuan,ZHANG Hongjun,LI Yeming,et al.The dynamic study of locomotives under saturated adhesion[J]. Vehicle System Dynamics,2011,49(8):1321-1338.

[10] 姚远,张红军,罗赟,等.机车传动系统扭转与轮对纵向耦合振动稳定性[J].交通运输工程学报,2009,9(1):17-20.

YAO Yuan,ZHANG Hongjun,LUO Yun,et al.Torsional-longitudinal coupling vibration stability of drive system for locomotive[J]. Journal of Traffic and Transportation Engineering, 2009,9(1):17-20.

[11] 张晨旭,杨晓东,张伟.含间隙齿轮传动系统的非线性动力学特性的研究[J].动力学与控制学报,2016,14(2):115-121.

ZHANH Chenxu,YANG Xiaodong,ZHANG Wei.Study on non-linear dynamics of gear transmission system with clearance[J]. Journal of Dynamics and Control,2016,14(2):115-121.

[12] THEODOSSIADES S,NATSIAVAS S.Non-linear dynamics of gear-pair systems with periodic stiffness and backlash[J].Journal of Sound and Vibration,2000,229(2):287-310.

[13] CHEN Zaigang,SHAO Yimin.Dynamic features of a planetary gear system with tooth crack under different sizes and inclination angles[J].Journal of Vibration and Acoustics,2013, 135(3):031004.

[14] LIANG X H, ZUO M J, PANDEY M.Analytically evaluating the influence of crack on the mesh stiffness of a planetary gear set[J].Mechanism and Machine Theory,2014,76:20-38.

[15] WAN Zhiguo,CAO Hongrui, ZI Yanyang,et al.An improved time-varying mesh stiffness algorithm and dynamic modeling of gear-rotor system with tooth root crack[J].Engineering Failure Analysis,2014,42:157-177.