基于改进变分模态分解的滚动轴承弱故障特征提取

2021-07-25朱群伟朱丹宸张明悦

机电工程技术 2021年6期

朱群伟，朱丹宸，张明悦

（1.海军装备部广州局驻湛江地区军事代表室，广东湛江 524005；2.海军士官学校，安徽蚌埠 233012；3.92601部队，广东湛江 524009）

0 引言

滚动轴承作为旋转机械的核心零部件，常工作在高温、高压和重载的条件下，故障时有发生，严重影响着整个设备的安全稳定运行。因此，开展针对滚动轴承的状态监测和故障诊断方法研究就有着很强的现实意义。基于振动信号进行轴承故障诊断是常见的方法，但由于实际诊断过程中，为了不破坏设备的原有结构，振动信号测量所用传感器大多只能布置在设备外侧的壳体上，轴承故障冲击信号在向外传递过程中需经过复杂的传递路径，造成较大的能量损失，加上多振源导致的强背景噪声干扰，给准确有效的滚动轴承故障特征提取带来极大的难度。

由于经验模式分解（EMD）能够处理非平稳、非线性信号，因此在滚动轴承故障诊断中得到了广泛使用，它将信号分解成一系列从高频到低频的本征模函数（IMF），通过选取合适的IMF进行信号重构，可以去除高频噪声信号，提高原始信号的信噪比。EMD在滚动轴承故障诊断上得到了较为成功的运用［1－2］，然而，端点效应和模态混叠一直影响着EMD的处理效果。镜像延拓［3］、匹配延拓［4］等方法都在解决端点效应问题上取得了不错的进展。集合经验模式分解［5］（EEMD）的提出较好地解决了模态混叠对EMD的影响。EEMD虽然有着不错的抗混叠效果，但其结果受到总的集合次数以及所添加白噪声大小的影响，同时EEMD的计算效率也要低于EMD。为了提高EEMD的有效性，Wang等［6］引入AR模型去除信号中的低频干扰成分；刘永强等［7］在对原始信号进行EEMD分解的基础上，利用峭度与自相关函数峰态系数选取信号重构所需的IMF分量，突出了故障特征信息。在EEMD的基础上，为了减少误差提高算法的效率，又提出了如总体平均经验模式分解［8］以及完全噪声辅助集总经验模态分解［9］等方法。

变分模态分解（VMD）作为一种新的信号自适应分解方法，它可以将信号分解问题转化为约束优化问题，其最大特点就是摒弃了EMD和EEMD中的递归筛分思想［10］，具有精度高和收敛速度快的优点，近年来在降低背景噪声，提取轴承故障特征方面得到了广泛运用［11－12］。但VMD的信号分解结果受到分解层数K以及二次项惩罚系数α的影响较大，针对这个问题，Huang等［13］提出基于空间尺度的改进VMD方法，实现分解层数及其他参数的优化选取。Yan等［14］提出利用布谷鸟搜索算法对VMD中的参数K及α进行自动选取。Gu等［15］提出将最小平均包络熵作为目标函数，实现了VMD中参数的自动寻优，Gong T等［16］和Jiang X等［17］对VMD的分解策略进行改进，实现参数的较优选取。

基于滚动轴承故障信号在时域内呈现周期性冲击的特点，考虑到背景噪声的干扰，结合VMD算法的优势，本文提出一种改进VMD方法并将其用于滚动轴承故障特征提取。利用VMD处理原始信号，以最大加权频域相关峭度为目标函数，对参数K和α进行寻优，获取故障特性最为显著的分量；通过Teager能量算子进一步增强故障冲击特征，并结合快速傅里叶变换获得Teager能量谱，实现滚动轴承故障特征的有效提取。

1 基本原理

1.1 VMD 的基本原理

VMD是一种自适应信号分解方法，可以将原始信号分解为多个互不相关的稀疏子信号uk，也就是本征模态函数，VMD的实质为构造和求解变分问题的过程。

VMD求解的约束变分问题可以表示为：

式中：uk＝｛u1，u2，…，uk｝，为各模态函数；ωk＝｛ω1，ω2，…，ωk｝，为各模态函数的中心频率。

引入二次惩罚因子α和拉格朗日乘法算子λ（t）用于求解上述的约束变分问题，α可以用于确保重构精度，λ（t）用于保证约束的精确执行。由此得到的增广拉格朗日表达式为：

采用乘子交替方向法求解等式（3），代表式（1）的最优解。在确定分解模态个数的基础上，初始化子模态、对应的中心频率和拉格朗日算子λ1。子模态及其中心频率ωk分别按照以下公式进行更新：

拉格朗日算子也可根据下式进行更新：

式中：τ为噪声容忍度。

迭代终止的收敛条件为：

式中：ε为收敛误差。

通过分析VMD的分解过程可知，分解层数K、二次惩罚因子α、噪声容忍度τ和收敛误差ε必须预先进行确定，其中，K和α对分解结果有较大影响，而τ和ε对结果结果的影响相对小很多，因此，本文对前两个参数进行了优化选取，另两个参数采用默认值。

1.2 加权频域相关峭度

峭度（Kurtosis）和相关峭度（Correlated Kurtosis，CK）作为一种表征信号中冲击特性的指标，常被用于选取合适的分析频带，进而准确提取出滚动轴承故障特征，且相关峭度指标在探测周期性冲击成分时比峭度指标更为有效。然而，峭度指标容易受单一冲击噪声的干扰，相关峭度指标也容易受到低频谐波信号的影响，这些成分的存在均会导致指标值较大，影响最佳分析频带的准确和有效选取。因此，本文考虑到将两种指标相结合，并充分利用信号的时域和频域特性，提出了加权频域相关峭度指标，并将其用于选取VMD中的参数K和α。

1.3 Teager能量算子

Teager能量算子通过信号的瞬时值和其微分形式来估计信号源产生的动态信号的总能量，常被用于增强信号中的瞬态特征，且由于只需要知道信号的幅值和导数，计算复杂性较低，效率高，使用较为方便。

假设一个连续信号x（t），其Teager能量算子ψ可以表示为：

式中：（t）和（t）分别为信号x（t）对时间t的1阶和2阶导数。分析Teager能量算子的特性，假设线性无阻尼振动系统由质量为m和刚度为k的弹簧组成，其运动方程可以表示为：

式（12）的解为：

对其做能量算子计算可得：

分析式（14）可知，Teager能量算子的输出为信号瞬时振幅与瞬时频率的平方之积，与传统的将能量定义为幅值平方相比，增加瞬时频率的平方项，由于冲击成分的瞬时频率较高，因此，利用Teager能量算子能够增强瞬态冲击成分，适用于增强滚动轴承故障冲击成分。

对于离散时间信号x（n），其Teager能量算子为：

2 故障特征提取流程

利用前一节提出的加权频域相关峭度值，可对VMD的参数进行优化选取，并借助Teager能量算子对信号中的故障特征进行增强，本文所提算法具体流程如图1所示。

图1 本文算法的流程Fig.1 The flow chart of the proposed method

（1）给定VMD的分解层数K和二次惩罚因子α的选取范围及变化步长；

（2）利用所有K和α的组合进行VMD计算，计算每个组合下分解得到的模态函数的加权频域相关峭度值，选取其中最大的作为该组合下的结果，最终，能够输出最大加权频域相关峭度的就是最优参数；

（3）利用最优K和α对信号进行VMD分解，得到多个本征模态函数；

（4）将最大加权频域峭度对应的分量作为有效分量，并利用Teager能量算子进行分析；

（5）利用快速傅里叶变换分析步骤（4）的处理结果，提取出滚动轴承的故障特征。

3 仿真信号分析

为了验证本文提出的算法在滚动轴承故障特征提取上的有效性，本节以滚动轴承内圈故障为例，构造了仿真信号，表达式如下：

式中：转频fr＝30 Hz；Ai为以1／fr为周期的幅值调制；fn＝3 200，为系统的共振频率；r为系统的阻尼系数，r＝0.05；T为相邻两个冲击成分的间隔，T＝1／185 s，也即是故障特征频率为185 Hz；ti为第i个周期内，由滚动体滑移引起的延迟，且ti＝0.01T～0.02T；B（t）为谐波分量，用以模拟信号中的干扰成分，且fm＝100；n（t）为随机白噪声，本文构造仿真信号时，将信噪比设定为－10 dB，信号的采样频率为32 768 Hz。

仿真信号如图2所示。从图2（a）的仿真信号时域波形中能够看出，受到背景噪声和谐波干扰的影响，周期性冲击成分几乎难以识别，且从图2（b）所示包络谱中，也未能准确提取出故障特征频率及其倍频成分，说明仅通过包络分析难以从仿真信号中有效提取出故障信息。

图2 仿真信号Fig.2 The simulation signal

进一步采用本文提出的算法分析仿真信号，为了避免信号分解不完全且尽量减小计算量，根据仿真信号特点，将K的取值范围确定为2～8，变化步长为1，将α的取值范围确定为100～6 000，且步长为100。通过计算各参数组合下的结果，确定的最优参数组合为K＝7，α＝2 200，此时加权频域相关峭度值可取到最大1.62×10－9，该参数组合下原信号的VMD分析结果如图3所示。

图3 VMD分析仿真信号的结果Fig.3 The decompotion results of the simulation signal using VMD

各分量的加权频域相关峭度值如图4（a）所示，从图中可以看出，第3阶模态分量u3的值明显大于其他各分量，因此选取该分量进行后续分析，该分量的时域波形如图4（b）所示，通过观察可以看出，经过VMD分析，原始信号中的噪声干扰得到了一定的去除。

图4 析仿真信号的结果Fig.4 Analyzed results of the simulation signal using the proposed method

采用Teager能量算子和快速傅里叶变换分析该分量得到的最优Teager能量谱如图4（c）所示，图中准确提取出了转频fr＝42 Hz及其倍频成分84 Hz，轴承内圈故障特征频率fi＝185 Hz及其2倍、3倍频也能够得到明显识别，由此可以判断出该轴承存在内圈故障。图4（d）为所选分量的包络谱，对比Teager能量谱可以发现，通过Teager能量算子计算，轴承故障特征得到了显著增强，各特征频率处的幅值都明显增大，且由于Teager能量算子计算方便，因此，采用Teager能量谱就有明显的优势。

为了体现VMD参数选取的重要性，随机选取VMD的参数，将分解层数K设定为4，将二次惩罚系数α设定为500，在分解得到的各分量中，选取加权频域相关峭度最大的分量进行Teager能量谱计算，结果如图5（a）所示，图中，转频及其倍频成分能够得到识别，轴承故障特征频率185 Hz得到了提取，但不明显，且仅有其二倍频得到了有效提取。相比于本文算法的分析结果，图5（a）中各特征频率处的幅值都显著降低，且提取出的故障信息不如图4（c）明显和丰富，说明了本文算法的有效性。

图5 仿真信号的对比分析结果Fig.5 Comparison results of the simulation signal

为了进一步体现出本文算法的有效性，下面利用EMD对仿真信号进行预处理，借助加权频域相关峭度选取故障信号最为丰富的模态分量，并计算其Teager能量谱，结果如图5（b）所示。图中，故障特征频率能够得到识别，但不够突出，其二倍频成分较为明显，但整个能量谱中噪声干扰较为明显，滚动轴承内圈故障信息不够丰富，说明利用EMD难以有效抑制噪声的干扰，突出了利用改进VMD进行故障特征提取的有效性。

4 实验数据验证

为了进一步验证本文算法在处理实际滚动轴承故障信号时的有效性，本文利用实验室滚动轴承故障模拟平台进行了实验测试，该实验台的结果如图6所示，考虑到实际结构中，滚动轴承一般位于设备内部，而传感器只能布置在设备外侧，信号在传递过程中将经过较为复杂的传递路径，影响弱故障冲击特征的有效传递。因此，本次实验时，将测点选取在远离故障轴承的轴承支撑结构外侧，振动测量的方向为径向，从而增大了信号的复杂程度，增加了故障特征提取的难度。

图6 实验台Fig.6 The test rig

实验过程中使用的轴承型号为6010，通过激光切割的方法在轴承内圈沿轴线方向切一个宽0.2 mm，深0.2 mm的窄缝，用于模拟滚动轴承内圈故障，该轴承的具体参数如表1所示。实验过程中，轴承的转速为3 000 r／min，信号的采样频率为32 768 Hz，轴承内圈故障特征频率可以通过以下公式进行计算：

表1 测试轴承参数Tab.1 Parameters of test bearing

式中：Z、d、D、a和fr分别为滚动体个数、滚动体直径、节径、接触角和转频。

通过计算可知，6010轴承在3 000 r／min的转速下，其内圈故障特征频率的理论值为370 Hz。

实测滚动轴承故障信号如图7所示，从图7（a）的时域波形中能够看出，受到噪声干扰的影响，滚动轴承故障导致的周期性冲击成分几乎难以识别，且图7（b）的包络谱中，仅有转频成分和故障特征频率成分处的谱线幅值较为突出，故障信息不够丰富，提取效果欠佳，说明仅利用包络分析难以取得很好的效果。

图7 实测信号Fig.7 The experimental signal

当利用本文提出的方法处理该实测信号时，将K的取值范围确定为2～10，变化步长为1，将α的取值范围确定为100～6 000，且步长为100。经过多次计算，最终确定的最优参数组合为K＝4，α＝500，此参数组合下的加权频域相关峭度值可以取到最大1.03×10－8，此时原始信号的VMD分析结果如图8所示。

图8 VMD分析实测信号的结果Fig.8 The decompotion results of the experimental signal using VMD

计算各分量的加权频域相关峭度值，结果如图9（a）所示，根据计算结果，选取第4阶模态函数u4进行后续分析，该分量的时域波形如图9（b）所示，相比于原始信号，该分量中的噪声成分得到了明显抑制。分别计算该分量的Teager能量谱和包络谱，结果如图9（c）和图9（d）所示。Teager能量谱准中能够明显识别出转频50 Hz及其倍频，滚动轴承内圈故障特征频率368 Hz（与计算得到的理论值相接近，加工和转速误差可能导致了偏差）及其2倍和3倍频均得到了准确提取，围绕特征频率及其倍频成分，间隔为转频的调制边频带也得到了有效提取，可以判断，该滚动轴承存在内圈故障。相比而言，虽然内圈故障特征频率及其倍频成分的谱线在包络谱中也较为明显，但幅值明显较低，部分调制边频带处的幅值不够突出，由此进一出体现出利用Teager能量谱进行故障特征提取的有效性。

图9 本文算法的分析实测信号的结果Fig.9 Analyzed results of the experimental signal using the proposed method

对比过程与仿真类似，在分析实测信号时，随机选取VMD的参数，将K和α分别设定为7和3 000，选取加权频域相关峭度最大的分量进行分析，其Teager能量谱如图10所示。对比图10（a）和图9（c）可知，图10（a）中仅能识别转频分量，轴承内圈故障特征频率及其2倍频，各频率分量处的幅值明显降低，部分不够突出，同时，故障特征频率的3倍频，转频的倍频成分未能得到有效识别，提取出的故障信息不够丰富，说明了对VMD算法的参数进行优化选取的重要性。利用EMD分析实测信号得到的Teager能量谱如图10（b）所示，图中，转频、故障特征频率及其倍频成分能够得到识别，但幅值并没有十分突出，干扰频率成分的幅值较大，通过对比进一步突显了本文所提算法从强背景噪声中提取弱故障冲击时的有效性。