裂纹故障对斜齿轮时变啮合刚度及振动响应的影响分析

2019-08-31林腾蛟郭松龄赵子瑞

振动与冲击 2019年16期

林腾蛟，郭松龄，赵子瑞，魏静

(重庆大学机械传动国家重点实验室，重庆 400044)

齿轮传动系统作为机械设备中重要的运动和动力传递部件，广泛应用于航空航天、交通运输、船舶海洋等领域。由于齿轮系统结构复杂、工作环境较为恶劣，使其在服役期间容易发生故障，影响设备的可靠运行，严重时还会引起生产安全事故[1]。齿轮系统动力学分析可为早期故障诊断提供良好的理论基础，而准确有效的时变啮合刚度计算方法又是进行动力学分析的关键。因此，研究含裂纹故障的齿轮副时变啮合刚度，分析齿根裂纹对齿轮系统动态特性的影响规律，对齿轮系统的早期故障诊断具有重要意义。

针对裂纹故障齿轮副的啮合刚度计算及动态特性分析，目前国内外学者已经做了诸多研究工作。Yang等[2]从能量角度提出了运用势能法来计算齿轮副的时变啮合刚度；Chaari等[3]运用解析法求得了两种裂纹参数下单对直齿轮副的时变啮合刚度，分析了裂纹参数对齿轮啮合刚度的影响规律；万志国等[4]在运用势能法求解直齿轮副啮合刚度时，考虑了基圆与齿根圆不重合的问题，提出了一种改进算法；Chen等[5-6]基于能量法计算了裂纹直齿轮副的时变啮合刚度，对比分析了不同裂纹参数下齿轮副的啮合刚度和振动响应；王旭等[7]对含齿根裂纹故障的齿轮系统进行动力学仿真，而后采用正交小波包和倒频谱相结合的方法对故障特征进行了提取与分析；马锐等[8]建立了单对齿轮副扭转振动参数化动力学模型，对裂纹故障的非线性机理进行了研究；Saxena等[9]采用势能法求得了直齿轮副的啮合刚度，研究了不同裂纹长度对齿轮转子系统模态特性及频率响应特性的影响；Pandya等[10]采用有限元法探究了不同齿轮参数下裂纹扩展路径对啮合刚度的影响；冯刚等[11]建立了含裂纹弧齿锥齿轮三维接触模型，采用有限元法研究了裂纹对齿轮啮合刚度的影响；Wan等[12]采用势能法求解了裂纹斜齿轮啮合刚度，分析了不同裂纹深度和长度对啮合刚度的影响，并求解了故障状态下齿轮副的振动响应；胡兴龙[13]针对风电增速箱中的斜齿轮副，采用能量法计算啮合刚度，并对带裂纹损伤的风电传动系统进行了动力学仿真。以上文献大多是针对裂纹直齿轮开展的研究工作，而很少有文献对裂纹斜齿轮的啮合刚度及动态特性进行研究。Wan等和胡兴龙在研究裂纹斜齿轮的啮合刚度时，将轮齿简化为基圆上的变截面悬臂梁，忽略了裂纹区及基圆至齿根圆段轮齿的变形能，导致计算结果存在一定误差。

本文基于上述研究成果，提出一种含裂纹故障斜齿轮副时变啮合刚度的改进算法，将轮齿简化为齿根圆上的变截面悬臂梁，在考虑基圆与齿根圆不重合因素的同时，计入裂纹区的变形能，可有效减小计算误差；然后建立“齿轮-轴-轴承”耦合动力学模型，模拟齿根裂纹故障的振动响应，并分析不同裂纹参数对传动系统动态特性的影响规律。

1 无故障斜齿轮时变啮合刚度计算

齿轮时变啮合刚度是由啮合齿数和轮齿接触位置的周期变化所引起的时变函数，是引起齿轮副内部动态激励的重要原因，正确有效的时变啮合刚度算法有助于更好地探明齿轮系统的动力学特性。

由于斜齿轮存在螺旋角，啮合过程由轮齿一端面进入啮合到另一端面退出啮合，其齿面接触线长度具有时变性。斜齿轮啮合刚度算法与直齿轮不同，基于积分方法，将其轮齿沿齿宽方向切分成若干薄片，每一薄片可近似视为直齿轮，通过计算各薄片啮合刚度并积分，即可得到斜齿轮啮合刚度。

综合考虑轮齿的接触、弯曲、剪切、轴向压缩及基体弹性刚度，采用能量法计算齿轮副啮合刚度，其综合啮合刚度可表示为[14]

(1)

式中：kh为接触刚度；kb为弯曲刚度；ks为剪切刚度；ka为轴向压缩刚度；kf为基体弹性刚度。

1.1 赫兹接触刚度

根据赫兹接触理论，啮合斜齿轮副的赫兹接触刚度kh可表示为

(2)

式中：E为杨氏模量；L为接触线长度；ν为泊松比。

1.2 弯曲、剪切及轴向压缩刚度

Wan等计算斜齿轮啮合刚度时，将轮齿简化为基圆上变截面悬臂梁模型，当齿根圆半径小于基圆半径时，没有计算基圆与齿根圆间轮齿部分的变形能；而当齿根圆半径大于基圆半径时，多计算了基圆与齿根圆间轮齿部分的变形能，这将导致啮合刚度计算产生一定的误差。为此本文对其作了改进，将轮齿假设为齿根圆上的悬臂梁，如图1所示。基圆以上齿廓为渐开线(C1M和D1N段)，为便于公式推导，采用直线CC1和DD1简化表示齿根处齿廓。

(a)Rb>Rr

(b)Rb

采用积分方法，厚度为dy的每一片轮齿薄片的弯曲变形能可表示为

(3)

式中：d(y)为啮合点和齿根在齿高方向的距离；h(y)为啮合点和齿轮中心线的距离；dIx为距离齿根x处截面的面积惯性距。

当Rb>Rr时，上述变量的表达式为

(4)

将式(4)代入式(3)，积分后得到啮合力F作用下轮齿的弯曲刚度kb为

(5)

其中，

由于式(5)的分母不可积分，采用求和方法替代积分方法来求解，于是式(5)可重新表示为

(6)

式中：Δy=l/N；l为接触线长度L在齿宽方向上的投影，表示为l=Lcosβb；N为斜齿轮沿齿宽方向切分的切片份数。

同理可导出剪切及轴向压缩刚度的计算式为

(7)

其中，

(8)

当Rb

(9)

(10)

(11)

其中，

(12)

(13)

(14)

式中：α和αf分别为分度圆和齿根圆压力角。

1.3 基体弹性刚度

将每一片轮齿薄片的基体刚度对dy进行积分，即可得到轮齿基体刚度，其表达式为

(15)

式中：L*,M*,P*,Q*为尺寸系数，其计算式参见文献[15]；uf和Sf如图1(a)所示，其表达式为

uf=Rb[(α1+α2)sinα1+cosα1]-Rr

Sf=2Rrα3

2 裂纹故障斜齿轮时变啮合刚度计算

当齿轮副存在齿根裂纹时，赫兹接触刚度、轴向压缩刚度和基体刚度不变，仅弯曲刚度和剪切刚度会受到齿根裂纹的影响，故需重新推导裂纹轮齿的弯曲刚度和剪切刚度，进而得出裂纹故障斜齿轮的啮合刚度计算式。考虑到Rb>Rf与RbRf时的推导过程。图2为斜齿轮的齿根裂纹模型，图2中,ha为裂纹终止点到轮齿中心线的距离，αa为当载荷F的作用点移动到K点时的齿轮啮合角。根据F的作用点是否位于裂纹影响区域内，分两种情况来讨论齿根裂纹。

2.1 F的作用点位于裂纹影响区域以外

由于齿根裂纹的存在，不仅轮齿有效截面的面积发生改变，而且轮齿悬臂梁的有效长度也从齿根处延伸到了裂纹终止点的位置。因此在推导裂纹轮齿啮合刚度时，从裂纹终止点开始积分，计入裂纹区轮齿的变形能，以提高计算的准确性，此时轮齿有效截面的面积和惯性矩计算式为

图2 斜齿轮裂纹轮齿模型Fig.2 Cracked tooth model of helical gear

(16)

(17)

其中，

ha=Rbsinα2-q1sinv

(18)

(19)

将式(18)、式(19)代入式(3)，积分后可得裂纹轮齿的弯曲刚度表达式为

(20)

其中，

同理可导出剪切刚度的计算式为

(21)

2.2 F的作用点位于裂纹影响区域以内

当载荷F的作用点位于裂纹影响区域内时，轮齿有效截面的面积和惯性矩计算式为

Ax=(ha+hx)L

(22)

(23)

同理可得裂纹轮齿的弯曲刚度、剪切刚度的表达式分别为

(24)

(25)

分别采用有限元方法、Wan等研究中的未改进算法和本文中的改进方法计算基圆半径大于和小于齿根圆半径情况下的裂纹斜齿轮的啮合刚度。表1给出了两个算例的齿轮副相关参数，主、从动轮参数相同，计算结果如图3所示。由图3可知，本文提出的计算裂纹斜齿轮啮合刚度方法与有限元方法计算结果较为接近，算例一中由于Wan等的研究中没有计算基圆与齿根圆之间的变形能，导致啮合刚度偏大；算例二中由于Wan等的研究中多计算了基圆与齿根圆之间的变形能，导致啮合刚度偏小。因此，采用本文的改进方法可以更准确的计算含裂纹斜齿轮的啮合刚度。

表1 裂纹齿轮副计算参数

(a)z=22

(b)z=58图3 裂纹斜齿轮啮合刚度对比曲线Fig.3 Comparison of mesh stiffness of cracked helical gear

3 裂纹参数对齿轮时变啮合刚度的影响

齿根裂纹可通过裂纹角v、前端面裂纹长度q、后端面裂纹长度q0及裂纹长度Lc等参数加以确定。本节针对表2所示的齿轮副，通过调整裂纹参数v,q,q0,Lc的值来研究不同齿根裂纹角度、深度、长度对斜齿轮副啮合刚度的影响，本算例中主动轮的输入转速为1 500 r/min。

表2 裂纹齿轮副基本参数

设定q为2 mm，q0为0，Lc为30 mm，即裂纹为贯穿型且裂纹深度在齿宽方向上线性变化，针对v为30°,45°及60°三种不同裂纹角度情况进行斜齿轮副时变啮合刚度计算，求解结果如图4所示。由图4可知，不同裂纹角度对齿轮副啮合刚度的影响较小，随着裂纹角度的增大，啮合刚度略微增大。

图4 不同裂纹角度下的齿轮副啮合刚度Fig.4 Mesh stiffness of gears with different crack angles

设定v为45°，Lc为30 mm，q0为0，针对q为1 mm,2 mm,3 mm三种不同裂纹深度情况进行齿轮副时变啮合刚度计算，求解结果如图5所示。由图5可知，不同裂纹深度对齿轮副啮合刚度的影响较大，随着裂纹深度的增大，啮合刚度随之减小，且减小幅度有所减小。

图5 不同裂纹深度下的齿轮副啮合刚度Fig.5 Mesh stiffness of gears with different crack depths

设定v为45°，q为3 mm，q0为0，针对Lc为10 mm,20 mm,30 mm三种不同裂纹长度情况进行齿轮副时变啮合刚度计算，其计算结果如图6所示。由图6可知，当裂纹长度Lc为10 mm时，齿根裂纹对综合啮合刚度的影响较小，随着Lc的增大，啮合刚度随之减小，且减小幅度明显增大。

4 传动系统耦合动力学模型

以单级裂纹故障斜齿轮副传动系统为研究对象，为了考虑传动轴的轴向、横向及扭转变形，采用Timoshenko梁单元将各轴离散化，轴段单元如图7所示。由图7可知，A，B两节点都有六个自由度(三个沿坐标轴的平移自由度以及三个绕坐标轴的转动自由度)，轴段单元i的位移向量可表示为

qi={xi,yi,zi,θxi,θyi,θzi,xi+1,yi+1,zi+1,θxi+1,θyi+1,θzi+1}

(26)

图7 轴系单元坐标系Fig.7 Coordinate system for shafting element

轴段单元的自由振动方程为

(27)

式中：Mi,Ki,Ci,Gi分别为第i个轴段单元的质量矩阵、刚度矩阵、阻尼矩阵、陀螺力矩阵；Mi,Ki,Gi的具体形式可参见文献[16]。阻尼矩阵Ci采用瑞利阻尼的形式，即

Ci=αMi+βKi

(28)

式中：α,β分别为质量比例系数和刚度比例系数。

将齿轮副视为一对通过弹簧和阻尼器连接的刚性圆盘，考虑齿轮的质量、转动惯量及陀螺效应，采用啮合单元来模拟主、从动轮的啮合关系，如图8所示。啮合单元的动力学方程为

(29)

图8 斜齿轮副动力学模型Fig.8 Dynamic model of helical gears

式中：Mpg,Cpg,Kpg,Gpg分别为啮合单元的质量矩阵、阻尼矩阵、刚度矩阵、陀螺力矩阵；F为载荷向量；qpg为齿轮副的位移向量，其表达式为

qpg={xp,yp,zp,θxp,θyp,θzp,xg,yg,zg,θxg,θyg,θzg}

(30)

式中：x，y为横向自由度；z为轴向自由度；θx,θy为摆动自由度；θz为扭转自由度。

考虑轴承支撑刚度时，由于文中未考虑箱体的柔性，因此直接将支撑刚度矩阵耦合到轴承所对应的轴节点上，轴承的支撑刚度矩阵为

KB=diag(kxx,kyy,kzz,kθxθx,kθyθy,0)

(31)

式中：kxx,kyy为径向支撑刚度；kzz为轴向支撑刚度；kθxθx,kθyθy分别为绕x轴与y轴的扭转刚度。

将轴段单元与齿轮副啮合单元的动力学方程耦合起来，并考虑轴承支撑刚度，即可得到整个系统的动力学方程为

(32)

式中：M,C,G,K分别为整个系统的质量矩阵、阻尼矩阵、陀螺力矩阵、刚度矩阵;q为系统的节点位移向量；F(t)为系统的载荷向量。

采用表2中的齿轮参数进行动力学仿真，传动系统模型如图9所示，输入轴与输出轴几何尺寸相同，输入功率为40 kW，输入转速为1 500 r/min，轴段参数如表3所示，轴承支撑刚度如表4所示。

图9 传动系统耦合动力学模型Fig.9 Coupling dynamic model of the transmission system

Tab.3 The parameters of the input and output shaftmm

表4 轴承支撑刚度

5 传动系统动态特性分析

针对上节给出的传动系统参数及模型，采用轴系单元法建立传动系统的“齿轮-轴-轴承”耦合动力学模型，将计算所得的齿轮副时变啮合刚度与静态传动误差合成为齿轮副内部动态激励，作为传动系统的输入激励，采用Newmark-β数值法求解系统的振动响应。

静态传动误差e(t)可表示为

e(t)=ersin(2πfmt+φ)

(33)

式中：er为轮齿误差幅值；fm为啮合频率；φ为相位角，取φ=0。

假设主动轮中仅有一个轮齿存在齿根裂纹，从动轮为正常齿轮，当q=1时的齿轮副时变啮合刚度,如图10所示。由图10可知，当裂纹轮齿进入啮合时，齿轮副的时变啮合刚度出现局部减小。

图10 齿轮副时变啮合刚度(q=1)Fig.10 Time-varying mesh stiffness of gears when q=1

图11和图12分别给出了无故障以及前端面裂纹长度q=1时的输出轴左端轴承(轴承3)节点处的x向振动加速度的时域曲线及频域曲线，其中fm表示啮合频率。由图可知，无裂纹故障时的加速度时域曲线存在周期性冲击响应，其频谱仅含啮合频率及其倍频成分；含裂纹故障的响应曲线除了存在正常齿轮的啮合冲击外，当裂纹轮齿参与啮合时，还会产生更加明显的冲击响应，相邻两个冲击的间隔为0.04 s，为故障齿轮的转动周期。频域响应谱中出现了以啮合频率为中心的调制边频带，边频间隔Δf为25 Hz，即故障齿轮的转频，该边频成分可用于诊断传动系统的裂纹故障。

(a) 无裂纹故障

(b)含裂纹故障图11 轴承3的x向加速度时域响应Fig.11 Time domain responseof x direction acceleration for bearing 3

(a)无裂纹故障

(b)含裂纹故障图12 轴承3的x向加速度频域响应Fig.12 Frequency domain response of x direction acceleration for bearing 3

6 裂纹故障对传动系统动态特性的影响

为了研究不同裂纹尺寸对传动系统动态特性的影响，针对“3”节中的裂纹参数，计算不同裂纹深度和不同裂纹长度下的传动系统动态响应。图13为不同裂纹深度下轴承3处x向加速度的时域及频域曲线。由图13可知，当裂纹深度为1 mm时，轴承3处的加速度时域图中出现了冲击现象，随着裂纹深度的增加，时域冲击的幅值增加，冲击现象愈加明显。从加速度频域图中可知，当存在齿根裂纹时，啮合频率及其倍频附近出现了频率间隔为25 Hz的边频带，其幅值随着裂纹深度的增加而增大。

图14为不同裂纹长度下轴承3处x向加速度的时域及频域曲线。由图14可知，当裂纹沿齿宽方向的长度为10 mm时，轴承3处的加速度响应与正常齿轮传动系统的差别不大，随着裂纹长度的增加，时域冲击的幅值增大，当裂纹长度扩展到整个齿宽(Lc=30 mm)时，冲击现象十分显著。从频域图中可知，随着裂纹长度的增加，峰值频率及其边频带处振动加速度幅值均有较大的增加。