旋转机械故障特征提取的全矢二元经验模态分解方法研究

2019-05-27黄传金雷文平李凌均孟雅俊

振动与冲击 2019年9期

黄传金，雷文平，李凌均，孟雅俊，赵静

(1.郑州工程技术学院机电与车辆工程学院，郑州 450044；2.浙江大学机械工程学院，杭州 310007；3.郑州大学机械工程学院，郑州 450052)

基于振动信号监测的旋转机械故障诊断是研究热点[1-2]。由于振动信号的非线性，一些先进信号处理方法如小波变换[3]、局部均值分解(Local Mean Decomposition,LMD)[4]、局部特征尺度分解(Local Characteristic Decomposition,LCD)[5]、变分模态分解(Variational Mode Decomposition,VMD)[6]和经验模态分解(Empirical Mode Decomposition,EMD)及其改进方法[7-8]等被用于旋转机械故障诊断。

这些方法促进了旋转机械故障诊断技术的进步，但也存在不足之处。如小波变换需预置与故障波形相匹配的基函数，不同的旋转机械可能存在不同的故障波形特征，即使是同一机械，随其故障的发展相应的故障波形也可能会变化，从而导致小波变换缺少自适应。另外，小波变换受Heisenberg测不准原理的制约，难以同时获得最好的时频分析结果。LMD、LCD和EMD均根据信号的极值点拟合信号包络线，从而将信号分解成固有模态函数，它们均具有较强的自适应能力，但均缺少数学理论支撑，而且也存在模态混叠、端点效应等问题。VMD通过较强的数学推导得出，具有完备的数学体系支撑，但其分解尺度的选择还需人工经验。

另外，由于故障位置、旋转机械各向异性、多轴系标高不一致等原因，旋转机械不同方向信号强度可能存在差异，导致基于单通道的信号分析方法提取的故障特征可能不准确。因此，一些同源信息融合技术被提出。如Bently提出的全频谱[9]，孟建等首创的全息谱[10]和韩捷等提出的全矢谱[11]。这些融合技术基于傅立叶变换基础之上，分析静态信号取得较好的效果。为融合非平稳信号特征，一些和EMD、LMD等相联合的同源信息技术被提出[12-13]。但由于EMD或LMD各自分解两个方向的信号，可能存在尺度不统一的问题，给信息融合带来难题[14]。

为了直接处理双通道和多通道信号，近年来有学者将LMD、EMD向二元、三元获多元领域做了扩展。比如复局部均值分解(Complex Local Mean Decomposition,CLMD)[15]、二元经验模态分解(Bivariate EMD,BEMD)[16]、三元经验模态分解(TEMD)[17]、多元经验模态分解(Multivariate EMD,MEMD)[18]等。由于多通道信号提供了更多的旋转机械故障特征信息(如轴心轨迹、进动方向等)，基于多通道信号的旋转机械故障诊断方法也成为一个热点。比如文献[14]提出了基于BEMD的旋转机械电气和机械故障诊断，文献[19]提出了基于CLMD全矢包络谱的转碰摩故障特征提取方法，但它们均不能判断旋转机械进动方向。为获取转子的进动方向信息，文献[18]和文献[20]分别将MEMD和全频谱、全息谱结合进行故障诊断。在二元领域即可判断进动方向，多元方法将原始信号投影到球面，增加了计算的复杂性。另外，全频谱基于正圆半径振动强度评价难于准确表达振动强度，全息谱的图谱直接用椭圆表示，其分辨率低且难以进行能量分析。

试验分析发现基于数据驱动的分解方法会发生模态混叠或模态分裂，既存在两种情况：①一个IMF包含不同的频率信号；②同一频率成分的信号被分解到不同的IMF。对于情况①，虽有模态混叠，但同一个信号相位、强度信息被完整保留；当发生模态分裂时(情况②)信号的相位和强度信息将被破坏。

基于以上原因，本文提出了全矢BEMD的旋转机械故障特征提取新方法。将正交方向的振动信号作为一个整体分解成系列CIMFs(Complex Intrisic Mode Functions)，为防止模态分裂情况发生，提出了基于相关系数的CIMFs组和方法获取新的CIMFs，然后运用全矢谱技术融合两个方向上振动信号的特征信息，获取的特征更为全面、清晰，可便捷的判断进动方向，并解决了一元信号分析方法(如EMD、LMD)处理二元信号容易遇到的尺度不统一问题。

1 全矢谱技术和BEMD原理

1.1 全矢谱技术

稳态时各谐波的运动轨迹是一个椭圆，而在运动学上转子的椭圆运动轨迹可以看成两个频率相同而进动方向相反的正圆轨迹分运动的合成，如图1所示。全频谱分别取两个圆的半径，分别按正、负频率轴展现出来，再通过比较同一谐波下正、负频率下的幅度大小，来判定该椭圆的进动方向以及定性的振动强度。由于全频谱的图谱直接用正圆半径表示，图谱中难以准确表达各谐波下的振动强度。

图1 两圆合成轴心轨迹图Fig.1 Schematic diagram of orbit composed by two circular

(1)

式中

xp=xpc+j×xpsxr=xrc+j×xrs

(2)

假定{zn}构成的复序列为

{zn}={xn}+j×{yn}

(3)

对其做傅里叶变换，可以得到{zn}的离散傅里叶变换{Zk},利用傅里叶变换性质和数学推导，可得：

(4)

式中，k=0,1,2,…,N/2-1,φpk、φrk分别是谐波ωk的正进动圆和反进动圆的初始相位角；RLk、RSk分别表示谐波ωk椭圆的长半轴和短半轴；α表示RLk与x-轴的夹角，称为主振矢角；φαk是圆盘中心沿椭圆轨迹运动时的相位角。

当Xpk>Xrk时合成后得到的盘心运动轨迹将与Xpk为半径的圆的进动方向一致，称为正进动；当Xpk

1.2 BEMD原理

BEMD认为二元信号为快速旋转信号和慢速旋转信号的叠加，对复信号z(t)，基于BEMD的分解过程如下[21]

步骤1 将每个复数信号投影到N个方向，并将同一时刻t内同一方向所有复数信号的投影值相加，形成N个随时间变化的信号，第k个方向的投影值可用公式(5)表示

pφk(t)=Re(e-jφ(k)×z(t)); 1≤k≤N

(5)

其中，φ(k)=2πk/N，N为投影方向的个数。

步骤4 计算所有切线均值

(6)

步骤5 从原始信号中减去均值信号获得

h(t)=z(t)-m(t)

(7)

理想情况下，h(t)满足复CIMF条件，认为h(t)是第一个复旋转分量c1。通常h(t)不满足，对h(t)重复以上过程，假设经过L次迭代得到hL(t)满足CIMF条件，则hL(t)为第一个复旋转分量c1。

步骤6 用原始信号z(t)减去c1得到一个剩余信号r(t)

r1(t)=z(t)-c1(t)

(8)

看r1(t)是否是旋转分量，如不是，则对其重复以上步骤，直到剩余信号没有旋转分量为止。

经过BEMD算法分解，一个复数z(t)可用下式表达

(9)

ci(t)为复数形式的旋转分量;M为旋转分量的个数，rk为剩余信号。

2 复旋转分量的相关系数和全矢BEMD

2.1 复旋转分量的相关系数

试验发现运用BEMD分解复数信号时，同一频率信号容易分解到不同的CIMF中。本文定义了复旋转分量相关系数以便将同一成分信号分解到一个CIMF。

假设第i阶复旋转分量

ci(t)=cix(t)+j×ciy(t)

(10)

cix(t)、ciy(t)分别是ci(t)的实部和虚部。

相关系数是用以反映变量之间相关关系密切程度的统计指标，可以表征X和Y之间线性关系紧密程度的量。当相关系数较大时，X和Y相关程度较好；当相关系数较小时，X和Y相关程度较差。对旋转机械振动信号而言，因不同方向的信号可能存在差异，因此可定义ci(t)的相关系数

ρi,i+1=max{|ρix|,|ρiy|}

(11)

其中

(12)

式中:1≤i≤M-1,ρix为相邻CIMFs实部的相关系数；ρiy为相邻CIMFs虚部的相关系数；Real为求实部；Imag为求虚部；Cov为求协方差；Var为求方差。

当同一频率信号分解到相邻的CIMF中，则相应的相关系数较大；相关系数较小时，则两个相邻CIMFs不存在频率相同的信号。

2.2 全矢BEMD

为防止同一频率信号被分解到不同的CIMF中，完整提取相应椭圆长、短轴信息、角度信息和进动方向，提出了全矢BEMD方法。步骤如下：

步骤1 运用BEMD算法将复数z(t)分解成系列复旋转分量ci(t)与一个剩余函数之和，如式(9)表示；

步骤2 用式(11)计算相邻复旋转分量的相关系数ρi,i+1；

步骤3 设置一阈值λ，当ρi<λ，原旋转分量不变，当ρi≥λ时将相邻的复旋转分量相加，得到一个新的复旋转分量；

步骤4 重新对复旋转分量排序，得到新的复旋转分量序列，即：

(13)

式中，W为重排后的复旋转分量的个数，W≤M。

步骤5 运用式(4)对重排后的复旋转分量ci的实部和虚部进行信息融合，得到相应椭圆的长轴、短轴和角度信息，进而判断进动方向和提取故障特征信息。

运用全矢BEMD可以对信号进行多尺度的自适应地分解，无需预设基函数，可以考察信号的局部信息；在全矢BEMD中，根据相关系数得到的新的复旋转分量克服了模态分裂现象，完整地表述了信号的相位和强度信息。

3 案例分析

3.1 松动故障

在北京东方所柔性转子试验台INV1602上模拟松动故障(转速为2 400 r/min，基频X=40 Hz)；采用正交采样技术获取的水平方向和垂直方向的振动信号分别为x和y(采样频率为2 048 Hz，采样时长0.5 s)以及频谱如图2所示。从图2可知，水平方向与垂直方向的振动信号频谱存在差异，而且高频信号和分频信号不太清晰。

图2 柔性转子试验台模拟的松动故障信号及其频谱Fig.2 The looseness fault signals and its spectrums

运用G.Riling等[22]编写的EMD程序，条件为默认值，分别对水平和垂直方向的松动信号x和y进行分解，分解结果如图3所示。结果发现信号x被分解成10个IMF分量和一个剩余信号，y被分解成9个IMF分量和一个剩余信号，分解得到的IMF分量个数不同。而全矢谱根据式(3)按一一对应原则得到复数序列，因此IMF分量个数不同将会给后续的基于全矢谱的信息融合带来难题。

图3 松动信号基于BEMD的分解结果Fig.3 The result of looseness signal based on BEMD

将水平和垂直方向的松动信号x、y组成一个复数z，即z=x+j×y。运用BEMD(投影方向N=8)分解复数信号，得到5个CIMFs，如图4所示。CIMFs的实部和虚部频谱如图5所示。从图5可知，c1和c2频率有重叠；c4与c3、c5均有频率重叠。频率重叠意味着同一频率信号被分解到不同的CIMFs中，此时用全矢谱获取的振幅强度和角度信息以及进动方向不能表示真实情况。

图4 松动信号基于EMD的分解结果Fig.4 The result of looseness signal based on EMD

图5 前5阶复旋转分量的频谱Fig.5 The first 5 order spectrum of complex rotational components

前5阶CIMFs的复相关系数如表1所示。从表1可知,ρ1，2、ρ3，4和ρ4，5较大，而ρ2，3可较小；设阈值λ=0.05。对复数z(t)进行全矢BEMD分解，则获得CIMF分量有c11、c22(c11=c1+c2;c22=c3+c4+c5)。松动信号基于全矢BEMD的分析结果和全频谱、傅里叶谱的对比图如图6所示。对比图6中同一频率下的长半轴RL1、RL2和c11、c12实部(c11x、c22x)与虚部(c11y、c22y)的幅值发现长半轴幅值较大；观察长半轴RL1、RL2与相应全频谱可知同一频率下在长半轴的幅值远大于全频谱的，如在长半轴中2.6X、3.6X对应的幅值分别为10.73、11.8，而在全频谱中分别为5.86、6.62。全频谱中通过比较正负频率下的幅值判断某一频率的合成轨迹的运动方向，比如2.6X、3.6X对应的幅值分别为5.86、6.62，-2.6X、-3.6X对应的幅值分别为4.86、5.18，可判定2.6X、3.6X合成的轨迹运动方向为正进动。而在椭圆短半轴中2.6X、3.6X的幅值分别为0.997、1.43，它们均大于零，则它们合成的轨迹运动方向为正进动，与采用全频谱的判定结果一致。分析c11的长半轴可发现存在分频，且两个相邻成分相差X，存在较弱的调幅现象，有冲击故障特征，但没有二倍频；观察c22可知0.45X～0.55X的分频幅值较大，与基频的比值在40.4%～50.5%之间，可确定发生了松动故障。另外，全矢BEMD还提供了主振矢的角度信息。

表1 前5阶CIMFs的相关系数Tab.1 Correlation coefficient of the first 5 order CIMFs

图6 松动信号基于全矢BEMD的分析结果Fig.6 The result of loose signals based on full vector BEMD

3.2 碰摩故障

在北京东方所柔性转子试验台INV1602上模拟碰摩故障(转速、基频、采样频率同上)。水平和垂直方向的信号x、y及其频谱如图7所示。

图7 柔性转子试验台模拟的松动故障信号及其频谱Fig.7 The rub impact signals and its spectrums

从图7可知，信号x与y的频谱强度存在较大差异，为融合二者包含的故障特征信息，运用全矢BEMD对其分析，分析结果与全频谱、傅里叶谱的对比图如图8所示。通过全矢BEMD，得到三个分量c11、c22和c33及其幅值、运转方向和角度信息。图8中c11是高频带信号，c22是基频和2倍频信号，c33主要是0.5X的低频信号。分别与之对应得实部、虚部和全频谱比较，发现同一频率信号在长半轴中的幅值大于其它任何一个；通过短半轴判断的各谐波合成后的运动轨迹进动方向和基于全频谱的判断结果一致。全矢BEMD还可得到各主振矢角。