基于改进连续变分模态分解的旋转机械故障诊断方法

摘要

基于改进连续变分模态分解的旋转机械故障诊断方法，包括：S1：采集振动信号；S2：设置平衡参数α的变化范围和增大步长；S3：以α增大的每一步值执行连续变分模态分解SVMD；S4：计算每一步所得模态分量的能量位置指标值EP；S5：找出每一步的最小EP指标值及其对应目标模态分量的序号；S6：画出α与每一步最小EP指标值及目标模态分量序号之间的关系曲线；S7：找出α的最优值；S8：找出α最优值对应的目标模态分量序号；S9：以α最优值执行SVMD，选取最优目标模态分量；S10：提取故障特征频率，实现旋转机械故障诊断。本发明克服了现有的SVMD中平衡参数难以确定、目标模态分量难以选取的问题，可以方便有效地实现旋转机械故障诊断。

说明书

技术领域

本发明属于机械故障诊断技术领域，尤其涉及一种基于改进连续变分模态分解的旋转机械故障诊断方法。

背景技术

旋转机械设备中通常包含有轴承、齿轮等旋转零部件，它们对旋转机械设备的正常运行有至关重要的影响。然而，在实际中，这些零部件经常要在动态载荷、交变载荷甚至过载的工况下运行，容易出现各种类型的故障损伤，如点蚀、裂纹、剥落和磨损等，严重影响旋转机械设备的使用性能和经济效益。因此，有必要采取行之有效的技术措施对旋转机械设备进行状态监测与故障诊断，以便及早发现损伤零部件，并对其及时维修或更换，防患于未然，旋转机械设备的维护维修成本。

旋转机械设备状态监测与故障诊断的一大常用方法是采集设备的振动信号，对振动信号进行处理与分析，识别出其中的故障特征，进而对设备的故障状态进行诊断。但由于旋转机械设备中的旋转零部件通常有很多个，不同旋转零部件激励的振动信号会相互耦合叠加，而损伤零部件激励的振动信号，在激励源到振动信号采集点的传递路径上又通常会存在较大程度的衰减，再加上背景噪声的干扰，要从传感器采集的振动信号中识别出微弱的故障特征，并不是一件容易的事情，这方面的研究亦是旋转机械设备故障诊断研究的一项重要内容。

目前，从旋转机械设备振动信号中提取故障特征一个常用方法是变分模态分解(VMD)。VMD算法由学者Dragomiretskiy和Zosso于2014年提出，其将振动信号分解成一系列具有不同中心频率的窄带本征模态分量，具有完备的数学理论支撑，抗噪性好，在旋转机械故障诊断领域得到了广泛的关注与研究。但VMD也存在一些不足之处，如模态分解层数和平衡参数难以优化设置、分解结果中目标模态分量难以选取，为此，学者Nazari和Sakhaei又在VMD的基础之上，于2020年提出了连续变分模态分解(SVMD)算法。SVMD在执行时不需要知道模态分解的层数，能够自适应地将信号分解成一系列模态分量。相比于VMD，SVMD少了模态分解层数这个参数，因而更容易执行，且具有更高的计算效率。但即便如此，SVMD仍然存在平衡参数难以设置、分解结果中目标模态分量难以选取的问题，这两个问题对于SVMD的性能有着至关重要的影响，尤其是在旋转机械振动信号处理与分析中，这两个问题的解决是提取故障特征、实现故障诊断的关键。

发明内容

本发明要克服现有技术中连续变分模态分解在旋转机械振动信号故障特征提取中存在的不足，提出一种基于改进连续变分模态分解的旋转机械故障诊断方法。

本发明的基于改进连续变分模态分解的旋转机械故障诊断方法，包括以下步骤：

S1：在旋转机械设备上，以f

S2：设置平衡参数α的变化范围为[α

α

其中，i＝1,2,3,…,N

S3：以α

S4：计算模态分量u

S5：从K

S6：根据最终得到的N

S7：利用α

S8：利用α

S9：以α

S10：计算最优目标模态分量u

优选地，步骤S4、步骤S5中，模态分量的能量位置指标值计算包括以下步骤：

S2-1：针对某一模态分量u(t)，计算其平方包络谱SES(f)，其中f表示频率；

S2-2：对平方包络谱SES(f)进行归一化处理，得到归一化平方包络谱NS(f)，具体方法为：

其中，max[SES(f)]表示求平方包络谱SES(f)中的最大幅值；

S2-3：利用归一化平方包络谱NS(f)计算故障特征频率成分的能量聚集性指标EC；

S2-4：利用归一化平方包络谱NS(f)计算故障特征频率的位置精度指标PA；

S2-5：计算模态分量u(t)对应的能量位置指标EP，公式为：

其中，p为调整系数，且0＜p≤1，β为平衡系数，β的计算公式为：

进一步优选地，步骤S10、步骤S2-1中，模态分量的平方包络谱计算包括以下步骤：

S3-1：针对某一模态分量u(t)，计算其平方包络信号u

其中，j为虚数单位，

S3-2：计算u

其中，

进一步优选地，步骤S2-3、步骤S2-5中，故障特征频率成分的能量聚集性指标EC的计算，包括以下步骤：

S4-1：对归一化平方包络谱NS(f)的幅值从大到小进行重新排列，得到重排的幅值序列SNS(n)(n＝1,2,…,N)，其中N表示频率点数；

S4-2：计算SNS(n)中第一个幅值与随后M个幅值之间差值的平均值，得到能量聚集性指标EC，公式为：

进一步优选地，步骤S2-4、步骤S2-5中，故障特征频率的位置精度指标PA的计算，包括以下步骤：

S5-1：找出归一化平方包络谱NS(f)中最大幅值所对应的频率值f

S5-2：计算位置精度指标PA，公式为：

其中，C为旋转机械中所包含旋转零部件的个数，f

本发明取得的积极效果在于：本发明利用连续变分模态分解方法，避免了现有的变分模态分解方法中模态分解层数难以确定的问题；利用所提出的能量位置指标克服了现有的连续变分模态分解方法中目标模态难以选取的问题；利用改进连续变分模态分解方法，可以得到包含充分而完整故障特征信息的最优目标模态分量，再结合平方包络谱分析，可以方便有效地提取出故障特征频率，实现旋转机械故障诊断。

附图说明

图1为本发明的实施流程图；

图2为本发明实施例中，滚动轴承振动信号的时域波形图；

图3为本发明实施例中，滚动轴承振动信号的平方包络谱。

图4为本发明实施例中，连续变分模态分解中，平衡参数α

图5为本发明实施例中，连续变分模态分解中，平衡参数α

图6为本发明实施例中，最优目标模态分量的时域波形图；

图7为本发明实施例中，最优目标模态分量的平方包络谱。

具体实施方式

下面结合附图和实施例对本发明作进一步说明。

参照图1，基于改进连续变分模态分解的旋转机械故障诊断方法，包括以下步骤：

S1：在旋转机械设备上，以f

S2：设置平衡参数α的变化范围为[α

α

其中，i＝1,2,3,…,N

S3：以α

S4：计算模态分量u

S5：从K

S6：根据最终得到的N

S7：利用α

S8：利用α

S9：以α

S10：计算最优目标模态分量u

步骤S4、步骤S5中，模态分量的能量位置指标值计算包括以下步骤：

S2-1：针对某一模态分量u(t)，计算其平方包络谱SES(f)，其中f表示频率；

S2-2：对平方包络谱SES(f)进行归一化处理，得到归一化平方包络谱NS(f)，具体方法为：

其中，max[SES(f)]表示求平方包络谱SES(f)中的最大幅值；

S2-3：利用归一化平方包络谱NS(f)计算故障特征频率成分的能量聚集性指标EC；

S2-4：利用归一化平方包络谱NS(f)计算故障特征频率的位置精度指标PA；

S2-5：计算模态分量u(t)对应的能量位置指标EP，公式为：

其中，p为调整系数，且0＜p≤1，β为平衡系数，β的计算公式为：

步骤S10、步骤S2-1中，模态分量的平方包络谱计算包括以下步骤：

S3-1：针对某一模态分量u(t)，计算其平方包络信号u

其中，j为虚数单位，

S3-2：计算u

其中，

步骤S2-3、步骤S2-5中，故障特征频率成分的能量聚集性指标EC的计算，包括以下步骤：

S4-1：对归一化平方包络谱NS(f)的幅值从大到小进行重新排列，得到重排的幅值序列SNS(n)(n＝1,2,…,N)，其中N表示频率点数；

S4-2：计算SNS(n)中第一个幅值与随后M个幅值之间差值的平均值，得到能量聚集性指标EC，公式为：

步骤S2-4、步骤S2-5中，故障特征频率的位置精度指标PA的计算，包括以下步骤：

S5-1：找出归一化平方包络谱NS(f)中最大幅值所对应的频率值f

S5-2：计算位置精度指标PA，公式为：

其中，C为旋转机械中所包含旋转零部件的个数，f

将上述发明应用于滚动轴承的振动信号处理中。轴承数据来源于西储大学轴承数据中心，所用数据为滚动轴承外圈故障损伤的振动数据，其具体故障为外圈6点钟位置处有直径为0.021英寸、深度为0.011英寸的单点损伤，加载电机施加的载荷为3HP，利用编码器测得驱动电机的输出转速为1721r/min，对应的转频为f

S1：从该滚动轴承对应的数据集中选取一段长度为L＝4000的振动信号x(t)，其时域波形如图2所示，其平方包络谱如图3所示，可见其中有较多的大幅值干扰分量，f＝102Hz所指示的滚动轴承故障特征频率成分并不突出，难以判断滚动轴承发生了故障损伤；

S2：设置平衡参数α的变化范围为[50,10000]，令α从50开始，以步长s

α

其中，i＝1,2,3,…,N

S3：以α

S4：计算模态分量u

S4-1：针对模态分量u

S4-2：对平方包络谱SES

其中，max[SES

S4-3：利用归一化平方包络谱NS

S4-3-1：对归一化平方包络谱NS

S4-3-2：计算SNS

S4-4：利用归一化平方包络谱NS

S4-4-1：找出归一化平方包络谱NS

S4-4-2：计算位置精度指标PA

PA

其中，f

S4-5：计算模态分量u

其中，p为调整系数，且p＝1，β为平衡系数，β的计算公式为：

S5：从K

S6：根据最终得到的N

S7：根据如图4所示的α

S8：根据如图5所示的根据α

S9：以α

S10：计算最优目标模态分量u

本说明书实施例所述的内容仅仅是对发明构思的实现形式的列举，本发明的保护范围不应当被视为仅限于实施例所陈述的具体形式，本发明的保护范围也及于本领域技术人员根据本发明构思所能够想到的等同技术手段。