CN115795703A - 一种断齿故障下平行轴齿轮系统动力学仿真模型建立方法 - Google Patents

Abstract

本发明提出一种断齿故障下平行轴齿轮系统动力学仿真模型建立方法，包括如下步骤：根据齿条加工齿轮的展成法，将齿条与齿轮的加工坐标进行转换，确定齿轮的齿廓坐标；确定齿轮啮合状态，通过齿轮啮合状态判断齿对中的齿轮断齿故障，分别计算齿对在齿面接触过程中的潜在接触点，得出含异常啮合状态的啮合点分布数据；计算齿轮副中每一齿对之间的初始间隙数据，初始间隙数据为初始状态下的啮合间隙数据；将初始间隙数据代入轮齿承载接触分析模型中，得到齿轮副的空载传递误差、载荷分配系数和时变啮合刚度；建立断齿故障下平行轴齿轮系统动力学仿真模型。本发明考虑了断齿故障导致的异常啮合比忽略异常啮合的仿真响应更接近实验结果。

Description

一种断齿故障下平行轴齿轮系统动力学仿真模型建立方法

技术领域

本发明涉及齿轮副传动系统啮合特性分析及动力学分析领域，具体涉及一种断齿故障下平行轴齿轮系统动力学仿真模型建立方法。

背景介绍

作为最常用的动力传动系统，直齿轮副的工作状态关乎机械系统的安全与运行稳定。受到制造工艺以及突变载荷的影响，直齿轮传动系统极易发生断齿故障，断掉的齿楔在正常的齿之间会急速的恶化齿轮啮合状态，导致其他齿对接连发生弯曲折断。进而产生机械设备失效，甚至产生威胁人生命安全的重大危害。因此，准确的判定直齿轮副传动系统断齿故障，及时采取挽救措施，是具有非常重要的实际意义。

对于齿轮时变啮合刚度的计算方法，专利(CN201911120877.3)中介绍了一种采用能量法计算齿轮副时变啮合刚度的方法，其中考虑了齿轮误差的影响。专利(CN201811449533.2)中同样采用能量法介绍了一种直齿圆柱齿轮内啮合齿轮副时变啮合刚度的计算方法。专利(CN201910303911.4)中介绍了考虑延长啮合和齿圈柔性的内啮合齿轮副啮合刚度计算方法。从以上可以看出，关于时变啮合刚度计算方法多数基于能量法，且鲜有关于断齿故障下时变啮合刚度的研究，此外目前关于异常啮合的仿真响应大多为忽略断齿故障的仿真响应，导致仿真结果不能很接近真实的实验结果。针对上述存在的问题，研究设计一种新型的断齿故障下平行轴齿轮系统动力学仿真模型建立方法，解决现有齿轮副传动系统啮合特性分析及动力学分析技术中所存在的问题是十分必要的。

发明内容

本发明为解决现有齿轮副传动系统啮合特性分析及动力学分析技术中存在的缺少断齿故障下时变啮合刚度的研究以及异常啮合的仿真响应中忽略断齿故障等问题，提出了一种断齿故障下平行轴齿轮系统动力学仿真模型建立方法，本方法考虑了啮合发生在理论啮合线外部的异常情况，精确计算齿对间的初始间隙，并基于轮齿承载接触分析模型得到了齿轮副在断齿故障情况下的空载传递误差和载荷分配系数，进一步求得齿轮副的时变啮合刚度，利用Timoshenko梁理论模拟齿轮副的支承轴，结合齿轮副的时变啮合刚度和空载传递误差，建立了断齿故障下平行轴齿轮系统动力学仿真模型。

本发明提供了一种直齿轮传动系统啮合间隙计算方法及故障分析模型，包括如下步骤：

S1.根据齿条加工齿轮的展成法，将齿条与齿轮的加工坐标进行转换，确定齿轮的齿廓坐标；

S2.根据步骤S1中得到的齿轮的齿廓坐标确定齿轮啮合状态，通过齿轮啮合状态判断齿对中的齿轮断齿故障，分别计算齿对在齿面接触过程中的潜在接触点，得出含异常啮合状态的啮合点分布数据；

S3.根据步骤S2中得到的所述含异常啮合状态的啮合点分布数据，计算齿轮副中每一齿对之间的初始间隙数据，所述初始间隙数据为初始状态下的啮合间隙数据；

S4.将步骤S3中获得的所述初始间隙数据代入轮齿承载接触分析模型中，得到齿轮副的空载传递误差和载荷分配系数，进一步求得齿轮副的时变啮合刚度；

S5.利用Timoshenko梁理论模拟齿轮副的支承轴，结合步骤S4中获得的齿轮副的时变啮合刚度和空载传递误差，建立断齿故障下平行轴齿轮系统动力学仿真模型。

优选的，所述步骤S1中，将齿条齿面与齿顶圆弧的切点C定义为加工起始点，沿齿条节线建立x轴，齿条齿槽中心线建立y轴，得到齿条的坐标系S_t，齿顶中心(b,a)在S_t坐标系中的坐标表式为：

式中，r_t为齿条齿顶圆弧半径，h_a为齿顶高系数，c为顶隙，m为模数，α_t为齿条的齿形角；

当齿条沿x_t方向向前移动时，齿轮节圆在齿条节线上沿切线滚动，齿轮中心O_g在坐标系S_t中的横坐标x_g0和纵坐标y_g0分别为：

式中，r_b为齿轮基圆半径，r为齿轮节圆半径，χ为齿轮的变位系数，α_c为C点的压力角，α_t1＝acos[r_b/(r+χ·m)]为节点压力角；

当齿条沿x_t方向的正方向移动Δx时，齿轮沿其中心旋转-Δx/y_g0，齿条的初始位置坐标S_t、齿条的当前坐标S_t1、齿轮的初始位置坐标S_g、齿轮当前位置坐标S_g1的转换关系为：

式中：

在坐标系S_t中，齿条齿廓r_tla和齿条齿顶圆弧r_tlb定义为：

对于齿轮齿廓任意一点i的啮合压力角

其对应的Δx由下式得到：

式中，

为点C的啮合压力角，

齿侧i点的压力角α_i与啮合压力角

的关系如下：

式中，θ_b为渐开线延长线与基圆相交点到轮齿中心线的圆心角的一半：

式中，z为齿数。

对于过渡曲线，θ是齿条齿顶圆弧的角位置，其变化范围为(α_t1,π/2)，θ与Δx的关系为：

Δx＝x_g0-b-a/tanθ (13)

齿条加工齿轮展成法处理后，齿轮的齿廓坐标为：

优选的，所述步骤S2中，在进行所述齿面接触之前先对齿轮状态进行假设：

(1)假设所有齿轮均完美加工，即齿面无误差；

(2)假设从动齿轮是完全固定的，只考虑主动齿轮的转动；

(3)假设齿轮副润滑良好，无齿对间的摩擦，接触力方向始终沿齿侧法向。

优选的，步骤S2所述判断齿对中齿轮断齿故障后，将从动齿轮存在断齿故障齿对的前一齿对定义为齿对1，将从动齿轮存在断齿故障的齿对定义为齿对2，将从动齿轮存在断齿故障齿对的后两对齿对分别定义为齿对3和齿对4；所述初始状态为施加外力之前的状态，所述初始间隙数据为初始状态下2πε/z_p周期内，齿轮副中每一齿对之间的的啮合间隙数据，其中，ε为重合度，z_p为主动齿轮的齿数，分别计算所述齿对1的潜在接触点A_1p和A_1g、所述齿对2潜在接触点A_2p和A_2g、所述齿对3潜在接触点A_3p和A_3g和所述齿对4的潜在接触点A_4p和A_4g；具体计算方式如下：

啮合过程中任意齿对的状态可分为三类:(a)提前啮入状态，(b)啮合状态，(c)延长啮出状态。首先，给主动齿轮一个旋转角Δθ_p，对应的从动齿轮旋转角Δθ_g＝Δθ_pz_p/z_g，其中，z_g为从动齿轮的齿数；

在初始状态Δθ_p＜2π/z_p，所述齿对1处于提前啮入状态，所述齿对1潜在接触点A_1p和A_1g的表达式为：

式中，O_g为从动齿轮中心坐标，O_p为主动齿轮中心坐标，α为压力角，α_ag为从动齿轮齿轮齿顶的压力角，r_ag为从动齿轮齿轮齿顶圆的半径，r_bp为主动齿轮齿轮基圆的半径，

为从动齿轮齿轮起始啮入点的啮合压力角，

为从动齿轮齿轮齿顶啮合压力角，

为所述齿对1的主动齿轮齿轮的潜在接触点，

为所述齿对1的从动齿轮齿轮的潜在接触点，α_A1p为潜在接触点

处的压力角；

经过一个啮合循环后，所述齿对1开始进入啮合状态，当2π/z_p≤Δθ_p≤2πε/z_p，所述齿对1潜在接触点A_1p和A_1g的表达式为：

所述齿对2的断齿高度为h，断齿齿顶啮合压力角

为：

式中，r_bg为从动齿轮齿轮基圆半径；

从动齿轮断齿齿顶的压力角α′_ag和断齿齿顶半径r′_ag的表示式为：

所述齿对2主动齿轮的旋转角Δθ_p＜tan[arctan(r_bg/r_ag)]–tan[arctan(r_bg/r′_ag)]，则所述齿对2处于提前啮入状态，所述齿对2的潜在接触点A_2p和A_2g的表达式为：

当所述齿对2处于啮合状态时，tan[arctan(r_bg/r_ag)]–tan[arctan(r_bg/r′_ag)]≤Δθ_p≤2πε/z_p，潜在接触点A_2p和A_2g的表达式为：

所述齿对3会经历啮合状态和延长啮出状态，当所述齿对3主动齿轮齿轮的旋转角

时，所述齿对3处于啮合状态时，接触点沿啮合线变化，潜在接触点A_3p和A_3g的表达式为：

当所述齿对3处于延长啮出状态时，

时，所述齿对3的主动齿轮齿顶与从动齿轮齿廓的啮合发生在接触线外，所述齿对3的潜在接触点A_3p和A_3g的表达式为：

式中，所述齿对3的Δθ₃的表达式为：

由于所述齿对4远离理论啮合线，因此只能处于延长啮出状态，假设主动齿轮转动Δθ₄的角度接触到从动齿轮，则所述齿对4的潜在接触点A_4p和A_4g的表达式为：

式中，所述齿对4的Δθ₄的表达式为：

优选的，所述步骤S2中，每对齿对的初始间隙为关于从动齿轮旋转中心的两个潜在接触点的弧长，表达式为：

式中，ε_i为每对齿对在对应轮齿啮合状态的初始间隙。

优选的，所述步骤S3中，

齿轮副的相容性方程为：

式中，λ_b为结构柔度，λ_c为接触柔度；I_1×Np为拥有1×Np维度的元素全为1的数组；I_Np×1为拥有Np×1维度的元素全为1的数组；f为施加在齿轮副上的总啮合力；LSR为载荷分配系数矩阵；STE为静态传递误差；ε为包含所有齿对初始间隙的向量，

齿轮副的时变啮合刚度k_pg的表达式为：

其中NLTE为空载传递误差，为目标轮齿啮合状态下所有潜在啮合齿对的初始间隙的最小值。

优选的，所述步骤S4中，主动齿轮和从动齿轮的相对投影位移p_pg的表达式为：

p_pg＝[sinα cosα 0 0 0 r_bp -sinα -cosα 0 0 0 r_bg]^T (29)

式中，α为齿轮副的压力角；r_bp和r_bg分别为主动齿轮和从动齿轮的齿轮基圆半径，主动齿轮和从动齿轮的啮合刚度矩阵K_pg为：

其中k_pg为齿轮副的时变啮合刚度。

优选的，所述步骤S4中，断齿故障引起的空载传动误差会使齿轮副产生一定的刚性运动。为了计算方便，根据达朗贝尔原理，这种刚体运动等价于一个虚力，主动齿轮和从动齿轮啮合虚力F_pg为：

上式导出的虚力作用于主动齿轮和从动齿轮的节点上，外力向量F的表示式为：

式中T_in为输入扭矩、T_out为输出扭矩。

优选的，所述步骤S5中，基于Timoshenko梁理论分别求解主动齿轮-输入轴结构的质量矩阵M_i和刚度矩阵K_i，从动齿轮-输出轴结构的质量矩阵M_o和刚度矩阵K_o。

优选的，所述步骤S5中，所述断齿故障分析模型的表达式为：

式中，q为系统位移矢量，M为质量矩阵，K为刚度矩阵，C为阻尼矩阵和F为外力向量。

本发明的有益效果为，本发明考虑了断齿故障导致的异常啮合，计算了齿轮副中每一齿对之间的初始间隙数据，初始间隙数据代入轮齿承载接触分析模型中，得到齿轮副的空载传递误差和载荷分配系数，进一步求得齿轮副的时变啮合刚度，随后利用Timoshenko梁理论模拟齿轮副的支承轴，建立断齿故障下平行轴齿轮系统动力学仿真模型，得到了时变啮合刚度和空载传递误差共同作用下的振动响应，考虑断齿故障引起的异常啮合的仿真响应比忽略异常啮合的仿真响应更接近实验结果，为断齿故障的诊断提供更准确的依据。

附图说明

图1是本发明方法的流程图；

图2是展成法由齿条生成的齿轮原理图；

图3是齿轮副的有限元模型；

图4是齿轮副的初始啮合状态；

图5是一对齿的啮合状态：(a)提前啮入状态；(b)延长啮出状态；

图6是发生断齿故障的齿轮；

图7是在一个求解周期内齿轮副的初始间隙和潜在接触点(1mm断齿高度为例)：(a)初始间隙；(b)潜在接触点；

图8是考虑不同从动齿轮断齿高度的直齿轮副的时变啮合刚度：(a)健康；(b)0.5mm；(c)1mm；(d)1.5mm；(e)2mm；(f)2.5mm；

图9是将本发明提出的方法与参考文献[Journal of Mechanical Design.141(2019)]中方法得到的时变啮合刚度对比；

图10是考虑到不同的断齿高度的空载传递误差：(a)1.5mm；(b)2mm；(c)2.5mm；

图11是考虑到不同的断齿高度的荷载分配系数变化：(a)1.5mm；(b)2mm；(c)2.5mm；

图12是在P1位置断齿高度1.5mm的齿轮副的接触状态；

图13是在P2位置断齿高度2mm的齿轮副的接触状态；

图14是在P3位置断齿高度2.5mm的齿轮副的接触状态；

图15是考虑不同输入力矩下1.5mm的断齿故障情况的空载传递误差和时变啮合刚度的分布；

图16是考虑不同输入力矩下2.5mm的断齿故障情况的空载传递误差和时变啮合刚度的分布；

图17是考虑到不同的断齿高度，齿顶接触百分比随扭矩变化情况；

图18是耦合齿轮副-轴耦合系统方案图；

图19是系统刚度矩阵组集示意图；

图20是健康状态下的时域波形：(a)仿真结果；(b)实验结果；

图21是健康状态下的频域振幅图：(a)仿真结果；(b)实验结果；

图22是1.5mm断齿故障状态下的时域波形：(a)考虑啮合线外啮合的仿真结果；(b)不考虑啮合线外啮合的仿真结果；(c)实验结果；

图23是1.5mm断齿故障状态下的频域幅值：(a)考虑啮合线外啮合的仿真结果；(b)不考虑啮合线外啮合的仿真结果；(c)实验结果。

具体实施方式

下面结合附图和实施例对本发明的实施方式作进一步详细描述。以下实施例用于说明本发明，但不能用来限制本发明的范围。

实施例1

步骤1，根据齿条加工齿轮的展成法，将齿条与齿轮的加工坐标进行转换，确定齿轮的齿廓坐标，

表1齿轮副的参数

沿齿条节线建立x轴，齿条齿槽中心线建立y轴，得到齿条的坐标系S_t，则齿顶中心(b,a)在坐标系S_t中的坐标可以表示为：

式中，r_t为齿条齿顶圆弧半径，h_a为齿顶高系数，c为顶隙，m为模数，α_t为齿条的齿形角。

当齿条沿x_t方向向前移动时，齿轮节圆在齿条节线上沿切线滚动，齿轮中心O_g在坐标系S_t中的横坐标x_g0和纵坐标y_g0分别为：

式中，r_b为齿轮基圆半径，r为齿轮节圆半径，χ为齿轮的变位系数，α_c为C点的压力角，α_t1＝acos[r_b/(r+χ·m)]为节点压力角；

当齿条沿x_t方向的正方向移动Δx时，齿轮沿其中心旋转-Δx/y_g0，齿条的初始位置坐标S_t、齿条的当前坐标S_t1、齿轮的初始位置坐标S_g、齿轮当前位置坐标S_g1的转换关系为：

S_t1→S_t:M_tt1＝TRSX(Δx)

S_g→S_g1:M_g1g＝ROTZ(-Δx/y_g0)

式中：

在坐标系S_t中，齿条齿廓r_tla和齿条齿顶圆弧r_tlb定义为：

对于齿轮齿廓任意一点i的啮合压力角

其对应的Δx由下式得到：

式中，

为点C的啮合压力角，

齿侧i点的压力角α_i与啮合压力角

的关系如下：

式中，θ_b为渐开线延长线与基圆相交点到轮齿中心线的圆心角的一半：

式中，z为齿数，

对于过渡曲线，θ是齿条齿顶圆弧的角位置，其变化范围为(α_t1,π/2)，θ与Δx的关系为：

Δx＝x_g0-b-a/tanθ

齿条加工齿轮展成法处理后，齿轮的齿廓坐标为：

至此，齿轮的齿廓坐标被求得，根据求得的齿面坐标可以建立齿轮模型。

基于有限元的离散齿轮啮合模型如图3所示，齿轮副的参数见表1。

步骤2，将步骤1中得到的齿轮的齿廓坐标确定齿轮啮合状态，通过齿轮啮合状态判断齿对中的齿轮断齿故障，分别计算齿对在齿面接触过程中的潜在接触点，得出含异常啮合状态的啮合点分布数据；

规定轮齿啮合的三种状态：提前啮入、正常啮合、延长啮出。根据每种轮齿啮合状态并考虑从动齿轮的断齿故障，分别计算齿轮过程中的潜在啮合点。本发明在仿真过程中考虑了异常啮合现象。

根据前述推导，获得了齿轮的齿廓坐标，下面以此为依据进行齿面接触分析。在进行分析之前先进行一些合理的假设：

(1)假设所有齿轮均完美加工，即齿面无误差

(2)假设从动齿轮是完全固定的，只考虑主动齿轮的转动

(3)由于齿轮副润滑良好，忽略齿对的摩擦，即接触力方向始终沿齿侧法向。

判断齿对中齿轮断齿故障后，将从动齿轮存在断齿故障齿对的前一齿对定义为齿对1，将从动齿轮存在断齿故障的齿对定义为齿对2，将从动齿轮存在断齿故障齿对的后两对齿对分别定义为齿对3和齿对4；初始状态为施加外力之前的状态，初始间隙数据为初始状态下2πε/z_p周期内，齿轮副中每一齿对之间的的啮合间隙数据，其中，ε为重合度，z_p为主动齿轮的齿数，分别计算齿对1的潜在接触点A_1p和A_1g、齿对2潜在接触点A_2p和A_2g、齿对3潜在接触点A_3p和A_3g和齿对4的潜在接触点A_4p和A_4g；具体计算方式如下：

任意齿对的状态可分为三类：(a)提前啮入状态；(b)啮合状态；(c)延长啮出状态。首先，给主动齿轮一个旋转角Δθ_p，对应的从动齿轮旋转角Δθ_g＝Δθ_pz_p/z_g，其中，z_p为主动齿轮的齿数，z_g为从动齿轮的齿数，

在初始状态Δθ_p＜2π/z_p，齿对1处于提前啮入状态，齿对1潜在接触点A_1p和A_1g的表达式为：

式中，O_g为从动齿轮中心坐标，O_p为主动齿轮中心坐标，α为压力角，α_ag为从动齿轮齿轮齿顶的压力角，r_ag为从动齿轮齿轮齿顶圆的半径，r_bp为主动齿轮齿轮基圆的半径，

为从动齿轮齿轮起始啮入点的啮合压力角，

为从动齿轮齿轮齿顶啮合压力角，A_1p为齿对1的主动齿轮齿轮的潜在接触点，A_1g为齿对1的从动齿轮齿轮的潜在接触点，α_A1p为潜在接触点A_1p处的压力角；

经过一个啮合循环后，齿对1开始进入啮合状态，当2π/z_p≤Δθ_p≤2πε/z_p，齿对1潜在接触点A_1p和A_1g的表达式为：

考虑到齿对2的从动齿存在断齿故障，潜在接触点位置的确定需要特别讨论。假设齿对2的从动齿轮发生断齿故障，断齿高度为h，断齿齿顶啮合压力角

为：

式中，r_bg为从动齿轮齿轮基圆半径；

从动齿轮断齿齿顶的压力角α′_ag和断齿齿顶半径r′_ag的表示式为：

在理想条件下，齿对2应该从初始位置开始进入啮合状态，但齿对2存在断齿故障，延迟了啮合状态的发生。齿对2主动齿轮的旋转角Δθ_p≤tan[arctan(r_bg/r_ag)]–tan[arctan(r_bg/r′_ag)]，则齿对2处于提前啮入状态，齿对2的潜在接触点A_2p和A_2g的表达式为：

当齿对2处于啮合状态时，tan[arctan(r_bg/r_ag)]–tan[arctan(r_bg/r′_ag)]≤Δθ_p≤2πε/z_p，潜在接触点A_2p和A_2g的表达式为：

在求解过程中，齿对3会经历啮合状态和延长啮出状态，当齿对3主动齿轮齿轮的旋转角

时，齿对3处于啮合状态时，接触点沿啮合线变化，潜在接触点A_3p和A_3g的表达式为：

当齿对3处于延长啮出状态时，

时，齿对3的主动齿轮齿顶与从动齿轮齿廓的啮合发生在接触线外，齿对3的潜在接触点A_3p和A_3g的表达式为：

式中，齿对3的Δθ₃的表达式为：

由于齿对4远离理论啮合线，因此只能处于延长啮出状态，假设主动齿轮转动Δθ₄的角度接触到从动齿轮，则齿对4的潜在接触点A_4p和A_4g的表达式为：

式中，齿对4的Δθ₄的表达式为：

每对齿对的初始间隙为关于从动齿轮旋转中心的两个潜在接触点的弧长，表达式为：

式中，ε_i为每对齿对在对应轮齿啮合状态的初始间隙。

步骤3，根据步骤2中得到的含异常啮合状态的啮合点分布数据，计算齿轮副中每一齿对之间的初始间隙数据，初始间隙数据为初始状态下的啮合间隙数据；

齿轮副的相容性方程为：

式中，λ_b为结构柔度，λ_c为接触柔度；I_1×Np为拥有1×Np维度的元素全为1的数组；I_Np×1为拥有Np×1维度的元素全为1的数组；f为施加在齿轮副上的总啮合力；LSR为载荷分配系数矩阵；STE为静态传递误差；ε为包含所有齿对初始间隙的向量，

齿轮副的时变啮合刚度k_pg的表达式为：

其中NLTE为空载传递误差，为目标轮齿啮合状态下所有潜在啮合齿对的初始间隙的最小值。

步骤4，将步骤3中获得的初始间隙数据代入轮齿承载接触分析模型中，得到齿轮副的空载传递误差和载荷分配系数，进一步求得齿轮副的时变啮合刚度。

考虑不同断齿高度的时变啮合刚度曲线如图8所示，并利用有限元方法对提出的轮齿承载接触分析方法进行验证。基于200Nm输入力矩的情况，可以看出提出方法得到的时变啮合刚度曲线与有限元方法得到的时变啮合刚度曲线非常接近。因此验证本方法的有效性。

本例中齿轮副的重合度为1.712，表明一个齿在一个求解周期内经历了双齿接触、单齿接触和第二次双齿接触的交替过程。通过分析刚度曲线和空载传递误差曲线，可以总结以下规律：

(1)若断齿高度仅破坏第一个双齿区(图8b、图8c)，则断齿故障引起的时变啮合刚度减小仅体现在第一个双齿区，单齿区和第二个双齿区保持不变。在第一个双齿区，由于断齿故障的存在，只有一对齿参与啮合，齿轮副的承载能力被削弱。

(2)随着断齿高度增加，单齿区被破坏(图8d)或第二个双齿区被破坏(图8e、8f)。齿轮啮合状态发生了质的变化，开始出现了空载传递误差。具有断齿故障的齿轮副如果严格按照传动比旋转，受到空载传递误差的影响，会在单齿区发生“脱啮”的异常啮合状态。断齿故障影响下的单齿区时变啮合刚度为0(见图9中“o”标记点的曲线)，但如果齿轮副的刚性运动闭合了初始间隙，齿轮副最终会啮合在一起，导致该位置的时变啮合刚度值为非零(见图9实线)。

(3)当单齿区或第二个双齿区被破坏时(见图8d-f)。齿轮副的重合度被破坏。断齿高度越大，齿轮副的实际重合度越小。当齿轮副的实际重合度小于1时，齿轮副之间的空载传动误差不再为0(见图10)。

如果空载传动误差不为0，齿轮副的啮合必然发生在理论啮合线之外，导致异常啮合状态的出现。同时，由于断齿故障而造成的隔齿啮合的异常啮合状态(如齿对1和齿对3同时啮合，而齿对2完全不啮合)也有可能发生。断齿故障发生后，轮齿间的载荷分配系数也发生了变化。图11展示了1.5mm、2mm和2.5mm断齿高度下的荷载分配系数变化情况。同时为了更清楚地说明断齿状态下齿轮副的接触状态，在不同断齿高度下选取P1、P2、P3的位置，并绘出其承载状态图，如图12-14所示。

对于1.5mm断齿高度工况，断齿故障破坏单齿区。基于ANSYS软件，图12给出了P1位置静力求解后齿轮副的接触状态。从图中可以看出，三齿对同时啮合，其中齿对2的接触力最大，齿对3的接触力次之，齿对1的接触力最小。这与图11a中P1位置的载荷分配系数分布一致。值得注意的是，由于断齿故障，三齿均处于边缘接触状态，这对齿轮副的啮合稳定性影响很大，很可能导致齿轮传动系统完全失效。

对于2mm断齿高度，在P2位置，齿对2刚刚进入啮合状态，与齿对1、3一起形成三齿啮合状态(如图13所示)。同样，在这个位置，所有三对齿都处于边缘接触状态。其荷载分配系数如图11b中P2位置所示。此时，齿对1承受的接触力最大。齿对2刚刚进入啮合状态，分担的接触力最小。

对于2.5mm断齿高度，P3处的接触力状态如图14所示。由于断齿高度较大，齿对2只在第二双齿区末端承受载荷(见图11c)。齿轮副在理论单齿区附近有隔齿啮合现象(图14中P3位置为例)。在P3位置，齿对1和齿对3啮合，均属于边缘接触状态。由于两齿对的啮合位置严重偏离理论啮合线，导致啮合力偏离基圆切线，产生冗余的非切向分力，非切向分力会破坏齿轮传动系统的稳定性，应避免隔齿啮合。

综上，可以看出当发生断齿故障时，齿轮副的啮合特性发生了变化。当由于断齿故障导致实际重合度小于1时，开始出现非零空载传动误差。非零的空载传动误差诱导齿轮齿的接触点不再沿理论啮合线变化。同时，由于齿廓的缺失，齿轮副的接触状态恶化，容易出现三对边缘接触齿对同时啮合的状态，甚至隔齿啮合的现象的发生，最终破坏齿轮传动系统的稳定性。

在恶劣工况下，齿轮传动系统经常发生断齿故障。因此，在不同的输入转矩条件下，由于断齿故障导致的啮合特性变化值得探讨。下面的分析是基于1.5mm断齿和2.5mm断齿的工况。

考虑输入扭矩变化范围T＝(20、60、100、140、180)Nm，图15和图16分别为1.5mm和2.5mm断齿故障下的载荷分配系数曲线和时变啮合刚度曲线。从图中可以看出，时变啮合刚度值随着输入扭矩的增大而增大。这是由于齿轮副的接触柔度是与啮合力相关的。同时，随着输入扭矩的增大，延长啮合现象越来越明显，载荷分配系数曲线越来越光滑。

在1.5mm断齿故障下的载荷分配系数曲线和时变啮合刚度曲线中(图15)，单齿区被局部破坏。以T＝20Nm的工况为例，在阴影区“a”之前，齿对3(TP-3)单独啮合；在“a”区，TP-3逐渐退出啮合，齿对2(TP-2)逐渐进入啮合，但此时两齿均属于异常啮合，因此时变啮合刚度曲线出现波动峰值。经过“a”区后，齿轮副进入正常工作状态。

在2.5mm断齿故障下的载荷分配系数曲线和时变啮合刚度曲线中(图15)，第二个双齿区部分受损。在图15阴影区c的位置，TP-3退出啮合状态，TP-1提前啮入接触状态。与1.5mm断齿高度工况相似，时变啮合刚度也出现了波动峰。这也是两个非正常啮合齿对联合作用的结果。在图15中区域“b”和“c”之间的区域，只有TP-1参与啮合，与健康状态相比，这个区域的时变啮合刚度降低。在“c”区域TP-2逐渐参与啮合，与TP-1形成双齿区。在c区域后，齿轮副离开受断齿故障影响的区域，开始恢复正常啮合状态。

断齿故障的一个显著特征是齿顶参与啮合的比例增加。在不同输入扭矩条件下，在一个重合度乘以啮合周期的时间范围内，齿顶接触的百分比如图17所示。

在健康状态下，随着扭矩的增加，齿顶接触百分比从7％增加到18％。当发生断齿故障时，齿顶接触百分比明显增加。当断齿高度为2.5mm时，T＝20Nm的低扭矩工况下，齿顶接触百分比为26％，比健康状态下高19％。在T＝180Nm的大扭矩工况下，齿顶接触百分比为39％，比正常工况高21％。如此大的齿顶接触百分比增量会导致断齿故障发生后异常啮合状态的时间增加，使工作状态恶化，并可能导致齿轮传动系统的不确定性增加。

步骤5，利用Timoshenko梁理论模拟齿轮副的支承轴，结合步骤4中获得的齿轮副的时变啮合刚度和空载传递误差，建立断齿故障下平行轴齿轮系统动力学仿真模型。

图18所示为一直齿轮副试验台。从动齿轮的一个齿被人为植入1.5mm高的断齿故障。传感器安装在靠近轴承B3的箱体上。基于试验台实际参数，下面建立仿真动力学模型。

轴与齿轮副耦合系统的动力学模型描述如图18所示，主动齿轮是输入轴的一部分，从动齿轮是输出轴的一部分。四个轴承用于限制输入轴和输出轴的运动。基于Timoshenko梁理论分别求解主动齿轮-输入轴结构的质量矩阵和刚度矩阵M_i,K_i和从动齿轮-输出轴结构的质量矩阵和刚度矩阵M_o,K_o。各轴的参数见表2。

表2轴的参数

轴承的支承刚度值列在表3中。

表3支承刚度值

主动齿轮和从动齿轮的相对投影位移p_pg的表达式为：

p_pg＝[sinα cosα 0 0 0 r_bp -sinα -cosα 0 0 0 r_bg]^T

式中，α为齿轮副的压力角；r_bp和r_bg分别为主动齿轮和从动齿轮的齿轮基圆半径，主动齿轮和从动齿轮的啮合刚度矩阵K_pg为：

其中k_pg为齿轮副的时变啮合刚度。

系统刚度组集方案如图19所示。

对于输入轴结构和输出轴结构，采用瑞利阻尼法模拟内部阻尼：

C_j＝α_jM_j+β_jK_j

其中j＝I,O分别表示输入轴和输出轴。α_j和β_j为阻尼系数。它们可以由下式求得：

其中，f_1j和f_2j为轴结构j的前两个非零固有频率，ξ_1j和ξ_2j为结构j的固有频率f_1j和f_2j对应的模态阻尼值。这里模态阻尼值都取0.03。

齿轮副啮合阻尼计算公式为：

其中m_p和m_g分别近似为主动齿轮和从动齿轮的质量。

断齿故障引起的空载传动误差会使齿轮副产生一定的刚性运动。为了计算方便，根据达朗贝尔原理，主动齿轮和从动齿轮啮合虚力F_pg为：

系统的输入扭矩作用于主动轴施力节点的旋转自由度上，输出力矩作用于从动轴驱动节点的旋转自由度上。由上式导出的虚力作用于主动齿轮和从动齿轮的节点上。F的表示式为：

式中T_in为输入扭矩、T_out为输出扭矩。

至此，得到了质量矩阵M、刚度矩阵K、阻尼矩阵C和外力向量F。断齿故障分析模型的表达式为：

式中，q为系统位移矢量，M为质量矩阵，K为刚度矩阵，C为阻尼矩阵和F为外力向量。

对比例1

根据步骤5中的断齿故障分析模型仿真结果，通过试验验证，提出了对于一对直齿轮副传动系统，当从动齿轮上发生断齿故障时，在啮合频率及其倍频周围，出现了以从动齿轮转频为间隔的边频簇。

在输入扭矩为20Nm、输入转速为1000r/min的工况下，健康状态下仿真结果与实验结果对比如图20所示。从图中可以看出，健康状态下的时域波形呈现周期性波动，波动周期为齿轮副啮合周期，波动范围在±20m/s²以内。对于幅值谱，可以看出仿真结果中的峰值幅值主要为啮合频率f_m及其倍频，幅值谱中不存在边频(见图21a)。然而，在实验结果的幅值谱中(见图21b)，出现了一些f_p＝f_m/z_p间隔的边带。这是因为在安装试验台的过程中，电机轴与输入轴之间存在一定的同轴度误差，这种现象产生的不对中会导致输入轴的旋转频率f_p。调制系统啮合频率。尽管如此，仍然可以看出，仿真的幅值谱与试验的幅值谱在5f_m,6f_m以及7f_m处频率幅值较高。

考虑1.5mm高度的从动轮断齿故障，考虑在啮合线外发生的异常啮合状态的仿真结果和实验结果的时域波形如图22所示。从图上可以看出，二者的断齿故障引起的时域波形峰值幅值非常接近。为了表明考虑异常啮合状态的必要性，图中又画出了不考虑异常啮合状态的时域波形(如图22b)。也就是说，当空载传递误差不为0的时候，齿轮被判定为不发生啮合。可以看出不考虑异常啮合状态的时域波形中，由于时变啮合刚度突然变为0，断齿故障引起的时域波形峰值幅值要比实验结果高达4倍。这意味着提出方法考虑异常啮合状态是必要且准确的。

同时，又对比了本发明提出模型、有限元模型和实验得出的幅值谱，如图23所示。值得说明的是，为了更加清晰的描述断齿故障造成的边频分布规律，我们在5f_m左右放大了幅值谱。大体上看，仿真结果和试验结果在啮合频率及其倍频附近的幅值趋势有着很好的一致性。在试验结果的放大图中(图23c)，可以看出在5f_m的周围出现了以主动齿轮转频f_p为间隔的边频组，如在健康幅值谱分析中解释的那样，这是由于输入轴与电机轴的不对中造成的，但是在仿真过程中并没有考虑这一现象，因此，在仿真结果中没有发现这一规律。从动轮断齿故障导致了边频的产生，在仿真结果中(图23a,b),在5f_m附近存在着以f_g＝f_m/z_g为间隔的边频组。这与实验结果的断齿故障引起的边频分布规律是吻合的。对比图23a,b可以发现，当不考虑异常啮合状态时，从动齿轮断齿故障特征边频的幅值要明显的高于试验结果和考虑非正常啮合现象的仿真结果。这表明，考虑异常啮合状态的仿真结果更加准确。

1、通过扭矩变化分析发现，随着输入扭矩的增加，齿轮的载荷分配系数曲线逐渐变平，时变啮合刚度值逐渐增大。在齿轮断齿故障的工况下，当输入扭矩足够大时，齿轮副可能会出现三齿接触现象，且三对齿均处于边缘接触状态。如果断齿高度继续增加，齿轮副将出现隔齿啮合现象(两个不相邻的齿对同时啮合)。这种异常啮合状态意味着啮合力不在沿着基圆切向方向，齿轮传动稳定性被严重破坏。

2、在一个求解周期(重合度×啮合周期)内，齿顶参与啮合百分比受输入扭矩和断齿高度的影响。可以发现，当转矩增大时，无论有无断齿故障，齿顶参与啮合百分比都增大。当发生断齿故障时，断齿高度越高，齿顶参与啮合百分比越高。过高的齿顶参与啮合百分比会使齿轮副的运行处于不稳定状态，容易造成更严重的问题。

3、为了验证提出方法的优越性，建立了轴-齿轮副耦合动力学模型。该模型考虑了断齿故障下的空载传递误差和时变啮合刚度。通过验证发现，无论是时域波形还是幅值谱，所提出的考虑啮合线外异常啮合状态的仿真结果都与实验结果具有较强的一致性。如果不考虑异常啮合状态，断齿故障引起的时域波形峰值将非常高，且边频带幅值显著高于实验结果。通过对仿真结果和实验结果的幅值谱的对比分析，可以发现，当从动齿轮发生断齿故障时，啮合频率及其倍频附近会出现以从动齿轮转频为间隔的边频组。

本发明的实施例是为了示例和描述起见而给出的，而并不是无遗漏的或者将本发明限于所公开的形式。很多修改和变化对于本领域的普通技术人员而言是显而易见的。选择和描述实施例是为了更好说明本发明的原理和实际应用，并且使本领域的普通技术人员能够理解本发明从而设计适于特定用途的带有各种修改的各种实施例。

Claims (10)

1.一种断齿故障下平行轴齿轮系统动力学仿真模型建立方法，其特征在于，包括如下步骤：

S1.根据齿条加工齿轮的展成法，将齿条与齿轮的加工坐标进行转换，确定齿轮的齿廓坐标；

S2.根据步骤S1中得到的齿轮的齿廓坐标确定齿轮啮合状态，通过齿轮啮合状态判断齿对中的齿轮断齿故障，分别计算齿对在齿面接触过程中的潜在接触点，得出含异常啮合状态的啮合点分布数据；

S3.根据步骤S2中得到的所述含异常啮合状态的啮合点分布数据，计算齿轮副中每一齿对之间的初始间隙数据，所述初始间隙数据为初始状态下的啮合间隙数据；

S4.将步骤S3中获得的所述初始间隙数据代入轮齿承载接触分析模型中，得到齿轮副的空载传递误差和载荷分配系数，进一步求得齿轮副的时变啮合刚度；

S5.利用Timoshenko梁理论模拟齿轮副的支承轴，结合步骤S4中获得的齿轮副的时变啮合刚度和空载传递误差，建立断齿故障下平行轴齿轮系统动力学仿真模型。

2.根据权利要求1所述的一种断齿故障下平行轴齿轮系统动力学仿真模型建立方法，其特征在于，所述步骤S1中，将齿条齿面与齿顶圆弧的切点C定义为加工起始点，沿齿条节线建立x轴，齿条齿槽中心线建立y轴，得到齿条的坐标系S_t，齿顶中心(b,a)在S_t坐标系中的坐标表式为：

式中，r_t为齿条齿顶圆弧半径，h_a为齿顶高系数，c为顶隙，m为模数，α_t为齿条的齿形角；

当齿条沿x_t方向向前移动时，齿轮节圆在齿条节线上沿切线滚动，齿轮中心O_g在坐标系S_t中的横坐标x_g0和纵坐标y_g0分别为：

式中，r_b为齿轮基圆半径，r为齿轮节圆半径，χ为齿轮的变位系数，α_c为C点的压力角，α_t1＝acos[r_b/(r+χ·m)]为节点压力角；

当齿条沿x_t方向的正方向移动Δx时，齿轮沿其中心旋转-Δx/y_g0，齿条的初始位置坐标S_t、齿条的当前坐标S_t1、齿轮的初始位置坐标S_g、齿轮当前位置坐标S_g1的转换关系为：

式中：

在坐标系S_t中，齿条齿廓r_tla和齿条齿顶圆弧r_tlb定义为：

对于齿轮齿廓任意一点i的啮合压力角

其对应的Δx由下式得到：

式中，

为点C的啮合压力角，

齿侧i点的压力角α_i与啮合压力角

的关系如下：

式中，θ_b为渐开线延长线与基圆相交点到轮齿中心线的圆心角的一半：

式中，z为齿数，

对于过渡曲线，θ是齿条齿顶圆弧的角位置，其变化范围为(α_t1,π/2)，θ与Δx的关系为：

Δx＝x_g0-b-a/tanθ (13)

齿条加工齿轮展成法处理后，齿轮的齿廓坐标为：

3.根据权利要求1所述的一种断齿故障下平行轴齿轮系统动力学仿真模型建立方法，其特征在于，所述步骤S2中，在进行所述齿面接触之前先对齿轮状态进行假设：

(1)假设所有齿轮均完美加工，即齿面无误差；

(2)假设从动齿轮是完全固定的，只考虑主动齿轮的转动；

(3)假设齿轮副润滑良好，无齿对间的摩擦，接触力方向始终沿齿侧法向。

4.根据权利要求3所述的一种断齿故障下平行轴齿轮系统动力学仿真模型建立方法，其特征在于，步骤S2所述判断齿对中齿轮断齿故障后，将从动齿轮存在断齿故障齿对的前一齿对定义为齿对1，将从动齿轮存在断齿故障的齿对定义为齿对2，将从动齿轮存在断齿故障齿对的后两对齿对分别定义为齿对3和齿对4；所述初始状态为施加外力之前的状态，所述初始间隙数据为初始状态下2πε/z_p周期内，齿轮副中每一齿对之间的的啮合间隙数据，其中，ε为重合度，z_p为主动齿轮的齿数；分别计算所述齿对1的潜在接触点A_1p和A_1g、所述齿对2潜在接触点A_2p和A_2g、所述齿对3潜在接触点A_3p和A_3g和所述齿对4的潜在接触点A_4p和A_4g；具体计算方式如下：

首先，给主动齿轮一个旋转角Δθ_p，对应的从动齿轮旋转角Δθ_g＝Δθ_pz_p/z_g，其中，z_g为从动齿轮的齿数；

在初始状态Δθ_p＜2π/z_p，所述齿对1处于提前啮入状态，所述齿对1潜在接触点A_1p和A_1g的表达式为：

式中，O_g为从动齿轮中心坐标，O_p为主动齿轮中心坐标，α为压力角，α_ag为从动齿轮齿轮齿顶的压力角，r_ag为从动齿轮齿轮齿顶圆的半径，r_bp为主动齿轮齿轮基圆的半径，

为从动齿轮齿轮起始啮入点的啮合压力角，

为从动齿轮齿轮齿顶啮合压力角，A_1p为所述齿对1的主动齿轮齿轮的潜在接触点，A_1g为所述齿对1的从动齿轮齿轮的潜在接触点，α_A1p为潜在接触点A_1p处的压力角；

经过一个啮合循环后，所述齿对1开始进入啮合状态，当2π/z_p≤Δθ_p≤2πε/z_p，所述齿对1潜在接触点A_1p和A_1g的表达式为：

所述齿对2的断齿高度为h，断齿齿顶啮合压力角

为：

式中，r_bg为从动齿轮齿轮基圆半径；

从动齿轮断齿齿顶的压力角α′_ag和断齿齿顶半径r′_ag的表示式为：

所述齿对2主动齿轮的旋转角Δθ_p＜tan[arctan(r_bg/r_ag)]–tan[arctan(r_bg/r′_ag)]，则所述齿对2处于提前啮入状态，所述齿对2的潜在接触点A_2p和A_2g的表达式为：

当所述齿对2处于啮合状态时，tan[arctan(r_bg/r_ag)]–tan[arctan(r_bg/r′_ag)]≤Δθ_p≤2πε/z_p，潜在接触点A_2p和A_2g的表达式为：

所述齿对3会经历啮合状态和延长啮出状态，当所述齿对3主动齿轮齿轮的旋转角

时，所述齿对3处于啮合状态时，接触点沿啮合线变化，潜在接触点A_3p和A_3g的表达式为：

当所述齿对3处于延长啮出状态时，

时，所述齿对3的主动齿轮齿顶与从动齿轮齿廓的啮合发生在接触线外，所述齿对3的潜在接触点A_3p和A_3g的表达式为：

式中，所述齿对3的Δθ₃的表达式为：

由于所述齿对4远离理论啮合线，因此只能处于延长啮出状态，假设主动齿轮转动Δθ₄的角度接触到从动齿轮，则所述齿对4的潜在接触点A_4p和A_4g的表达式为：

式中，所述齿对4的Δθ₄的表达式为：

5.根据权利要求4所述的一种断齿故障下平行轴齿轮系统动力学仿真模型建立方法，其特征在于，所述步骤S2中，每对齿对的初始间隙为关于从动齿轮旋转中心的两个潜在接触点的弧长，表达式为：

式中，ε_i为每对齿对在对应轮齿啮合状态的初始间隙。

6.根据权利要求5所述的一种断齿故障下平行轴齿轮系统动力学仿真模型建立方法，其特征在于，所述步骤S3中，

齿轮副的相容性方程为：

式中，λ_b为结构柔度，λ_c为接触柔度；I_1×Np为拥有1×Np维度的元素全为1的数组；I_Np×1为拥有Np×1维度的元素全为1的数组；f为施加在齿轮副上的总啮合力；LSR为载荷分配系数矩阵；STE为静态传递误差；ε为包含所有齿对初始间隙的向量，

齿轮副的时变啮合刚度k_pg的表达式为：

其中NLTE为空载传递误差，为目标轮齿啮合状态下所有潜在啮合齿对的初始间隙的最小值。

7.根据权利要求6所述的一种断齿故障下平行轴齿轮系统动力学仿真模型建立方法，其特征在于，所述步骤S4中，主动齿轮和从动齿轮的相对投影位移p_pg的表达式为：

p_pg＝[sinα cosα 0 0 0 r_bp -sinα -cosα 0 0 0 r_bg]^T (29)式中，α为齿轮副的压力角；r_bp和r_bg分别为主动齿轮和从动齿轮的齿轮基圆半径，主动齿轮和从动齿轮的啮合刚度矩阵K_pg为：

其中k_pg为齿轮副的时变啮合刚度。

8.根据权利要求7所述的一种断齿故障下平行轴齿轮系统动力学仿真模型建立方法，其特征在于，所述步骤S4中，空载传动误差使齿轮副产生了刚性运动，根据达朗贝尔原理，主动齿轮和从动齿轮啮合虚力F_pg为：

上式导出的虚力作用于主动齿轮和从动齿轮的节点上，外力向量F的表示式为：

式中，T_in为输入扭矩、T_out为输出扭矩。

9.根据权利要求8所述的一种断齿故障下平行轴齿轮系统动力学仿真模型建立方法，其特征在于，所述步骤S5中，基于Timoshenko梁理论分别求解主动齿轮-输入轴结构的质量矩阵M_i和刚度矩阵K_i，从动齿轮-输出轴结构的质量矩阵M_o和刚度矩阵K_o。

10.根据权利要求9所述的一种断齿故障下平行轴齿轮系统动力学仿真模型建立方法，其特征在于，所述步骤S5中，所述断齿故障分析模型的表达式为：

式中，q为系统位移矢量，M为质量矩阵，K为刚度矩阵，C为阻尼矩阵和F为外力向量。

Images