基于延拓EMD和GDE识别瞬时阻尼比的行星齿轮箱故障诊断

李康强，冯志鹏，孙宏图，神克常

(1.鲁东大学 蔚山船舶与海洋学院，山东 烟台 264025；2.北京科技大学 机械工程学院，北京 100083)

行星齿轮箱由于其独特的结构和运行方式，相比于其他的机械旋转机构，其振动信号具有明显的调幅调频特性。此外，实际工程测试中，行星齿轮箱自身的多自由度和环境中的多振源也造成了信号的复杂非平稳多分量特性，常规方法对于提取行星齿轮箱的故障特征存在困难[1-2]。

目前，基于振动信号分析的行星齿轮箱的故障诊断，大多数研究集中在振动信号的时域统计指标分析和频域故障特征提取。Skrimpas等[3]利用不同阶啮合频率之间的残差信号的统计特征值检测行星齿轮与相关轴承的故障。Yoon等[4]利用单压电应变传感器的应变信号分析行星齿轮箱并提取故障特征。Wang等[5]利用K均值聚类方法对于不同的故障类型进行分类并判断设备健康状态。Feng等[6]利用本征时间尺度分解结合幅频解调提取行星齿轮箱的故障特征。针对于现有研究，常规方法在故障诊断分析中对创新性提出了更高的挑战。

通常行星齿轮箱在运行中，太阳轮以轴向自转作为动力输入传递扭矩，行星轮与太阳轮和固定在箱体上的齿圈同时啮合，设备运行时多个行星轮既有公转又有自转，种种因素造成行星齿轮箱振动信号的复杂性，同时给故障诊断带来了困难。在行星齿轮箱的稳定运行过程中，随着齿轮的转动，齿轮对的啮合位置不断变化，齿轮的变形量的不断变化引起的刚度随之变化，换言之，即系统的瞬时阻尼比是不断变化的。假想状态下的健康齿轮箱运行时，前一齿轮对的啮合状态与后一齿轮对是完全相同的，因此整个持续变化状态具有周期性。假如啮合齿面发生损伤时，系统运行会产生额外的能量冲击，因此，故障状态引起的能量变化区别于正常状态的能量变化，而能量变化会反映在瞬时阻尼比的变化上。基于此，我们提出齿轮齿面故障会引起阻尼比的瞬态变化，换言之，阻尼比瞬态变化与故障特征频率相关。那么如何从振动信号中准确获取瞬时阻尼比成为行星齿轮箱故障诊断的关键问题。

对于阻尼比的研究，一般针对于模态阻尼比识别方法，如对数衰减法、ITD(ibrahim time domain)法和STD(spare time domain)法等[7]，对于瞬时阻尼比的研究非常有限。Feldman[8]提出利用Hilbert变换估计瞬时阻尼比的方法。但是，该方法必须综合计算幅值、频率和瞬时无阻尼固有频率，多种类型的误差因素势必影响到瞬时阻尼比的准确度。Loh[9]利用经验模式分解(empirical mode decomposition,EMD)与Hilbert谱相结合来估算瞬时阻尼比，为识别非线性非平稳信号提供了一种潜在方法。但是上述研究只是提出估算瞬时阻尼比的方法，并没有提出具体应用。在现今研究中，尚未发现利用瞬时阻尼比进行故障诊断的案例。

Zayezdny等[10]提出了基于生成微分方程(generating differential equation,GDE)的信号表示方法，利用信号、信号微分以及信号各种变换组合之间的关系描述信号的结构和性质。李康强等[11-12]通过深入研究发现，生成微分方程在解调信号方面有其独特的优势，并可以应用于旋转机械故障诊断。通过继续研究发现，利用信号的各阶微分函数的非线性组合运算，可以方便高效的计算瞬时阻尼比等动力学参数信息。但是生成微分方程的非线性组合运算过程需要信号满足单分量要求，而实际工程测试的振动信号由于振源多杂造成其成分复杂。为了准确提取故障特征，本文发挥EMD在复杂信号分解方面的自适应性和高计算性，将行星齿轮箱振动信号分解为多个本质模式函数并从中挑选出满足生成微分方程需要的分量作为敏感单分量，继而利用生成微分方程准确估计敏感分量的瞬时阻尼比，最后通过识别谱图中的故障特征频率来诊断故障。利用不同的仿真信号与试验实测信号验证方法的有效性与实用性。

1 故障特征提取方法

1.1 延拓经验模式分解

为满足生成微分方程单分量的需要和增强故障特征，首先需要将复杂多分量信号进行分解。经验模式分解可以通过“筛选”将包含多个振荡模式的信号自适应的分解为符合以下条件的本质模函数(intrinsic mode function,IMF)[13]：

(1)整个数据长度下，IMF的极值数量等于过零点的数量或者至多相差一个；

(2)任意数据点，局部极大值的包络和局部极小值的均值等于0。

经验模式分解的简化步骤为：

步骤1计算源信号x(t)的所有局部极大值与局部极小值，利用三次样条插值函数可以获得上下包络线，源信号减去上下包络线的均值m1可以得到第一个分量h1=x(t)-m1；

步骤2如果上述第一个分量满足本质模函数的条件，则视为分解出的第一个本质模函数IMF1，如果不满足，则令源信号x(t)重复步骤1直至计算出的分量满足条件；

步骤3源信号减去第一个本质模函数得到的余量r1=x(t)-IMF1作为新的源信号重复步骤1和步骤2进行第二次“筛选”过程，重复循环n次，当残余量rn成为不能再继续分解的单调函数时，分解终止，残余量rn即源信号的趋势项；

可以看出，经验模式分解满足自适应性和完备性，而且分解的过程相对简便高效。但是实际应用中，EMD由于获取端点极值的拟合包络时会出现较大跳动，而且随着分解的进行，这种跳动会逐渐往中间传递，造成分解结果产生端点效应的问题，基于此干扰，我们可以利用信号对称延拓来解决这个问题。在EMD分解之前，对源信号进行如下对称延拓

(1)

式中：rnN为源信号x(t)的长度；Ne为延拓值,Nei=1,2,…,Ne。延拓之后数据长度增加2Ne个点，所以分解完成后数据需要对称剪切2Ne个点以便恢复原长度。

1.2 瞬时阻尼比计算

生成微分方程可以利用重定义的状态函数对信号进行映射变换，将信号分解为一系列的微分方程的解的合集形式，从而估算信号的动力学参数。不失一般性，对于一般的调幅调频信号可表示为

x=Acos(ωt+φ)

(2)

式中：A(t)为幅值调制信号；φ(t)为频率(相位)调制信号。

kx-(2δA+δω)δx-[kA-δA(2δA+δω)-ω2]=0

(3)

不失一般性，调制部分A(t)和φ(t)相对于载波频率是缓变的，因此可以近似为定值。从而关于其状态函数可以近似δA=0，δω=0，kA=0。继而推导瞬时频率为

(4)

系统在振动过程中的能量是周期性持续衰减的，就振动过程中某一时刻点来说，其能量衰减响应可以表述为形如自由衰减响应形式

x(t)=Ae-ξωntsin(ωdt+φ)

(5)

式中：ξ为阻尼比；ωn为无阻尼固有频率；ωd为有阻尼固有频率；φ为初相位。

如上所述，结合微分与状态函数可以推出式(5)的生成微分方程为

(6)

(7)

对于离散信号而言，实际计算过程可以由对称差分替代微分，由中心差商替代导数

(8)

同时离散信号的基本状态函数可以表示为

(9)

代入式(7)可以得到离散信号的瞬时阻尼比序列，其中步长Δt=1。由于计算过程只有数值的基本运算，因此运算简便高效。

1.3 特征提取

基于EMD和GDE识别瞬时阻尼比的故障诊断方法的步骤归纳如下，其流程图如图1所示。

图1 流程图

步骤1对待分析信号进行对称延拓；

步骤2利用经验模式分解延拓信号，得到数个本质模函数；

步骤3基于GDE计算各个本质模函数的瞬时频率，并选择最先分解出的且瞬时频率在啮合频率或其谐波上下波动的本质模函数作为敏感分量进行下一步分析[14]；

步骤4基于GDE计算敏感分量的瞬时阻尼比，并对称剪切恢复至原信号长度；

步骤5对瞬时阻尼比进行Fourier变换，得到瞬时阻尼比频谱；

步骤6根据阻尼比频谱峰值频率匹配故障特征频率，定位故障。

2 仿真验证

2.1 脉冲信号仿真

为直观显示瞬时阻尼比变化特征和故障特征，我们首先仿真建模一系列等衰减等间隔脉冲信号。建模生成

u(t-nT)+n(t)

(10)

图2(a)为仿真信号时域波形，可以看出时域波形为一系列等间隔等衰减的脉冲信号，其中间隔为T=0.1 s，一个间隔周期内可以发现信号从最大幅值衰减为0并且持续保持一段时间，因此瞬时阻尼比应为一段时间内为恒定值然后下降到0值。图2(b)验证了瞬时阻尼比的波形结构为方波，方波顶部为预设定阻尼比ξ=0.06，方波底部为0值。图2(c)为瞬时阻尼比频谱，可以明显看到突出的峰值频率对应于设定的故障特征频率及其倍频kfg，k=1,2,…。仿真结果证明了通过瞬时阻尼比可以提取故障特征频率。

2.2 齿轮箱信号仿真

根据行星齿轮箱的振动机理，其信号可以用调幅-调频信号模型描述[15]。不失一般性，只考虑齿轮啮合频率和单一故障齿轮的特征频率，则振动信号仿真形式设为

x(t)=[1-cos(2πfct)][1+Acos(2πfgt)]×

cos[2πfmt+Bsin(2πfgt)+φ]+n(t)

(11)

式中：fc为行星架旋转频率；fg为太阳轮故障特征频率；fm为齿轮啮合频率。设置fc=1.5 Hz，fg=27 Hz，fm=200 Hz，调幅系数A和调频系数B为1，初相位φ=0，仿真信号采样频率为1 500 Hz，加入信噪比为n(t)=6 dB的Gauss白噪声。

信号分析如图3所示。同脉冲信号分析流程一致，从图3(b)细节放大验证了瞬时阻尼比的突出幅值间隔为故障周期T=1/fg≈0.037 s。从图3(c)的瞬时阻尼比频谱中可以很明显的得到突出频率幅值匹配故障特征频率fg及其倍频，在噪声干扰下频谱底部细节放大图可以匹配行星架旋转频率fc及其倍频。通过仿真验证方法可以有效的提取故障，且抗噪性也较为明显。

图3 齿轮箱仿真信号分析

3 试验分析

1.驱动电机；2.信号调节器；3.转矩和速度传感器；4.行星齿轮箱；5.信号采集器；6.CZ-5型磁粉制动器。

图5 齿轮模拟故障

3.1 健康信号分析

图6 正常信号分析结果

3.2 太阳轮故障信号分析

图7 太阳轮故障信号分析结果

3.3 齿圈故障信号分析

由于三个行星轮同时与齿圈啮合，因此当齿圈齿面出现故障的时候，每个运转周期会出现三次啮合(三次故障脉冲)，所以我们需要重点关注齿圈故障特征频率的1/3倍频的幅值情况，齿圈故障信号的分析图，如图8所示。按照分析流程我们依然选择第一个IMF作为敏感分量。从图8(d)阻尼比频谱中可以看出，与正常信号(深色虚线)相比，浅色实线代表的齿圈故障信号显示出较多的峰值。这些峰值都与齿圈的故障特征频率fr及其1/3倍频相关，例如1/3fr，2/3fr，4/3fr，2fr，3fr等较高的突出幅值。通过匹配齿圈故障特征频率，可以验证齿圈出现了故障，与此同时频谱中并没有出现与太阳轮故障特征频率和行星轮故障特征频率的突出幅值，验证了其他齿轮的健康状态。

图8 齿圈故障信号分析结果

4 结 论

瞬时阻尼比能够反映系统能量的变化状态，而且波动特征对应于系统的特征频率。生成微分方程可以利用各阶次微分函数的非线性组合计算系统的瞬时阻尼比。基于此，本文提出结合经验模式分解和生成微分方程的行星齿轮箱故障诊断方法。本方法能够简便高效的计算瞬时阻尼比，并通过瞬时阻尼比频谱可以匹配故障特征频率，准确定位行星齿轮箱的故障。通过仿真信号验证和试验实测信号对比分析，证明了该方法在提取故障特征上的有效性和实用性。