基于VMD的Volterra模型奇异值熵的转子故障诊断方法*

杨恭勇 丁潇男 王珺琦 魏迎东 周小龙

(①东北电力大学工程训练教学中心，吉林 吉林 132012；②东北电力大学机械工程学院，吉林 吉林132012；③江苏川玛工业科技有限公司，江苏 昆山 215300；④北华大学机械工程学院，吉林 吉林 132021)

旋转机械是航空航天、铁路交通等众多行业的关键性设备，转子更是其核心部件之一[1]。受工作环境复杂性的影响，转子是旋转机械设备中的易损零件。据统计，由于转子故障导致的旋转机械故障占比50%以上[2]。当转子出现故障时，其振动信号表现出非平稳特性，传统时频分析方法无法实现故障特征的精准提取。因此，获取可有效表征转子状态的敏感故障特征已成为该领域研究的热点与难点。

当前，针对非平稳性转子故障信号特征提取的研究，主要以小波变换和经验模态分解[3](empirical mode decomposition，EMD)方法为主。孙嘉兵等[4]将小波分解方法同灰色相似关联度相结合并成功应用于转子系统的故障诊断。对于转子碰摩故障特征微弱难以识别的问题，王雷飞等[5]通过频谱校正和复合小波包变换相结合的方法，实现了碰摩故障特征的有效提取。金志浩等[6]以转子声振信号为研究对象，通过小波分解获取信号不同尺度能量百分比，并通过最小二乘支持向量机实现了不同材料转子碰摩故障的准确诊断。周玉平等[7]将ANSYS技术同EMD方法相结合，提出了一种可有效诊断滑动轴承转子系统裂纹故障的新方法。童靳于等[8]在EMD基础上，提出一种极点加权模态分解方法，并通过此方法实现了转子碰摩故障的有效诊断。

相较小波变换等线性处理方法，EMD具有自适应性特点，可根据信号特性将其分解成多个包含单一频率成分的固有模态函数(intrinsic mode function，IMF)分量，适用性更好，已成为该领域的主要研究方法之一。但由于EMD自身算法致使分解具有不稳定性，存在模态混叠问题[9]，这将导致分解结果中出现虚假分量，各IMF分量无法有效表征信号特征，从而影响信号特征提取的准确性。变分模态分解[10](variational mode decomposition，VMD)拥有坚实的理论基础，是一种非递归式自适应信号处理方法，在信号分解过程中可有效避免EMD分解时产生的模态混叠问题，保证信号特征提取的可靠性。

已有研究表明，自回归(auto regressive，AR)模型[11]的参数可有效反映系统状态的变化规律，但AR模型对非平稳信号的分析效果并不理想，模型受信号采样频率的影响显著[12]。在非线性系统建模研究中，Volterra模型的应用最为广泛，其计算效率高，模型参数也承载了系统状态的信息[13]。奇异值熵在信号信息量评估方面有显著优势且不受信号采样频率的影响[14]，若将其与Volterra模型相结合，在充分利用奇异值熵的信号信息评估优势的同时，也可有效避免采用时间对模型预测参数准确性的影响，从而强化信号故障特征的提取。

对于转子的故障诊断，常以其振动信号的频谱或包络谱为分析对象，但转子的多种故障特征频率都与其转频有关[15]，导致故障特征间存在较强的相似性，难以精确识别。而检测故障特征的信号中往往包含各种复杂的模糊联系[16]。因此，可采用模糊聚类方法对转子信号的故障类型进行识别。目前基于目标函数的模糊聚类方法最为常用，其中模糊C均值(fuzzy c-means，FCM)聚类算法理论具有最好的完备性。

鉴于上述分析，本文提出一种基于VMD的Volterra模型奇异值熵和FCM相结合的转子故障特征提取方法。通过实测信号的分析，证明了该方法的可行性与有效性。

1 VMD原理

通过预设尺度参数K的设置，VMD可将信号分解为K个中心频率是ωk的IMF分量。则可得到变分约束问题：

(1)

式中：∂t为对函数求时间t的偏导数；δ(t)为单位脉冲函数，uk(t)为第k个IMF分量。

求解上述问题，由此引入增广拉格朗日函数ζ，将约束问题转化为非约束问题：

(2)

式中：α为惩罚参数，以保证信号的重构精度；< >表示向量内积。

则IMF分量uk及其中心频率ωk可表示为：

(3)

(4)

VMD具体实现过程如下：

(2)执行循环n=n+ 1。

(3)根据式(3)、式(4)更新uk和ωk。

1.2 基于VMD的Volterra奇异值熵构建特征向量

设X(n) = [x(1),x(2),…,x(n)]为采集到的转子振动信号，U(n) = [u(1),u(2),…,u(n)]是敏感IMF，对其相空间进行重构，其中重构方法为延迟坐标法[13]，则：

U′(n)=[u(n),u(n-),···,u(n-(m-1))]

(5)

式中：m和分别为嵌入维数和时间延迟。

以U′(n)为输入，输出为y(n) =u(n+1)，则其Volterra级数展开式为：

(6)

式中：

(7)

其中：hk(i1,…,ik)是k阶Volterra核；q是Volterra展开级数；b是记忆长度。由于Volterra级数是无穷级数，故实际应用中以二阶Volterra级数为主，因此，本文选择二阶Volterra级数对敏感IMF进行预测，即：

(8)

令

W(n)=[h0,h1(0),h1(1),···,h2(0,0),h2(0,1),···,h1(b-1,b-1)]T

(9)

Z(n)=[1,u(n),u(n-),···,u(n-(b-1)),u2(n),u(n),u(n-),···,u2(n-(b-1))]T

(10)

则式(10)可表示为：

u(n+1)=ZT(n)W(n)

(11)

采用归一化最小均方自适应算法对上式进行求解，获取表征信号特性的模型参数。即由W(n)组成状态特征向量用以表征IMF分量u(t)的特征。

按上述方法，求解每个敏感IMF分量(设有k个)的状态特征向量并组成初始特征矩阵A：

A=[W1W2，…，Wk]T

(12)

经上述分析，转子振动信号X(n)的特征可由初始特征矩阵A所描述。

对A进行奇异值分解，获得奇异值。当转子出现不同故障时，敏感IMF分量的奇异值会产生相应改变，在此构造奇异值熵以定量描述这种变化[2]。

设矩阵A经奇异值分解得到的奇异值为p= {p1,p2,…,pk}，并对每个分量归一化，可得:

(13)

由信息熵定义可得Volterra模型奇异值熵计算公式为:

(14)

2 基于VMD的Volterra模型奇异值熵的转子故障特征提取

2.1 信号采集

为验证所提方法的有效性，在ZT-3转子实验台上进行故障模拟，提取故障数据并进行故障诊断分析。图1为转子振动信号采集装置。通过调速器调节实验台转速，实验台采用直流并励电动机驱动，电机额定电流为2.5 A，输出功率250 W；由输出端安装的光电传感器测得转速；转子加速度信号由AI005型加速度传感器获取，加速度信号通过MJ5936型动态信号测试器进行处理；并通过计算机获取实时测得的转子加速度信号。

试验时模拟正常、不对中、不平衡、动静碰摩和轴承座松动等5种状态。信号采集过程中，电动机转速2 700 r/min，采样频率2 000 Hz。不同状态下采集到转子振动信号的原始时域波形，如图2所示。由图2可知，各振动信号的时域波形虽有一定差异，但以实现转子工作状态和故障类型的准确诊断。

2.2 VMD关键参数的确定

当采用VMD方法对信号分解时，预设尺度数K和惩罚参数α是影响分解精度的重要参数[17]。

由VMD算法可知，经VMD所得各IMF分量的中心频率数值由低至高，分布合理。当K取得最优值后时，第K个IMF分量的中心频率取值最大，随着K的增加，其数值也不会明显增大。因此，本文以中心频率最大值法确定K的最优值。

采用VMD对图2b中转子不对中故障信号进分解，不同K值下各IMF分量的中心频率如表1所示。

表1 不同K值对应的各IMF分量中心频率

从表1中可以看出，在预设尺度数K= 5时，IMF分量中心频率取得最大值，并随着K值的增大，中心频率的最大未出现较大波动，表示此时VMD的分解效果最佳。因此，预设尺度数K取5。

惩罚参数α主要用于控制IMF分量的带宽。由于VMD算法具有较好的噪声鲁棒性，当信号经VMD分解后，信号内的干扰成分应得到一定滤除，使重构信号内表征信号特征的冲击成分增多，为准确刻画信号的复杂程度，在此选用多尺度模糊熵值作为惩罚参数α的选取评价参数。

计算10组转子不对中故障信号在预设尺度数K= 5的条件下，惩罚参数α在不同取值范围下多尺度模糊熵值的均值，结果示于图3。计算过程中，多尺度模糊熵值为嵌入维数g=2、相似容限r= 0.15·SD(待分解信号的标准差)、尺度因子h= 1, 2，…，10时的模糊熵均值。

由图3可知，当惩罚参数α=5 400，经VMD分解后重构信号的多尺度模糊熵值最小，由此说明重构信号内同故障特征相关的冲击成分所含最多，呈现较强的规则性和自相似性，故信号分解过程中，取惩罚参数α=5 400。

2.3 分解结果分析

按上节参数选择方法对转子不对中故障信号进行分解，结果如图4所示。由图4可知，分解结果较为合理，各IMF分量主要集中在其中心频率附近，说明该方法有效抑制了模态混叠问题。

为选择对故障特征敏感的IMF分量，参照文献[18]计算图4a中各IMF分量的能量熵增量Δqi，结果如表2所示。

表2 不对中故障信号各IMF分量的能量熵增量

由表2可知，VMD分解得到的IMF1～IMF3分量对于不对中故障敏感，选取它们作为敏感IMF分量。按上述方法对图2中转子不同状态下振动信号进行分析，为保证奇异值特征向量的一致性，根据能量熵增量数值选取3个对于转子状态最为敏感的IMF分量建立Volterra自适应预测模型并计算其奇异值熵，结果如表3所示。

表3 不同状态下转子振动信号的Volterra模型奇异值熵

通过对表3的分析可知，根据敏感IMF分量所求得的Volterra模型奇异值熵对转子故障非常敏感，故障类型不同，其Volterra模型奇异值熵的数值间差距较大。表4为不同采样频率情况下敏感IMF分量的Volterra模型奇异值熵。

由表4可知，同一故障类型在不同采样频率下其Volterra模型奇异值熵的数值差别较小，说明Volterra模型奇异值熵对于转子故障十分敏感，从而保证了故障特征提取的可靠性。

表4 不同采样频率下基于VMD的Volterra模型奇异值熵

2.4 基于FCM的转子故障诊断

选择电动机转速为2 700 r/min，采样频率为2 000 Hz条件下的转子正常状态及不对中、不平衡、动静碰摩、轴承座松动的故障数据进行分析。每种状态下采集3组数据，每组数据的采样时间为15 s。每个样本的截选时间长度为1 s，每种状态下从其中所采集到2组数据中分别截选10段，5种状态共计截取100段数据作为标准样本的原始数据数据，此外，从每种状态下所采集到的另1组数据中分别截取10段数据，5种状态共计截取50段数据作为检测样本的原始数据。

采用VMD方法对各标准样本数据进行分解，并根据能量熵增量准数值取出对于转子状态特征最为敏感的3个IMF分量，对其建立Volterra模型，以Volterra模型参数作为初始特征向量矩阵，求解其奇异值熵并将奇异值归一化处理，形成奇异值特征向量矩阵。

在特征向量计算过程中，先求解出转子5种状态下各20个标准样本的奇异值特征向量，以样本均值作为FCM算法的初始聚类中心；再对每种状态下各10个检测样本的奇异值特征向量进行求解，共获取50个检测样本的特征向量。表5为VMD分解后各敏感IMF分量所求得的样本初始聚类中心和部分检测样本的特征向量。

由表5可知，经VMD建立的Volterra模型奇异值特征向量得到的同类样本间波动性较小。由此表明，以本文所提方法求解出的上述参数作为特征向量对转子工作状态和故障类型诊断具有较好的可分性和诊断可靠性。

表5 VMD处理后得到的初始聚类中心及部分检测样本

基于VMD的Volterra模型奇异值熵得到的50个检测样本的FCM分类识别结果如图5所示，在FCM算法中加权指数m= 2，迭代停止阈值为10-6。图中类别1、2、3、4、5分别表示转子正常状态、不对中、不平衡、轴承座松动和动静碰摩故障。

由图5可知，本文所提方法对50个待检测样本均作出了正确的诊断，由此可以验证，本文所提方法的有效性。

为比较分析，采用集合经验模态分解(ensemble empirical mode decomposition，EEMD)方法对上述信号进行分解，根据敏感IMF选择算法，选取3个对于转子状态最为敏感的IMF分量建立Volterra预测模型，计算其奇异值特征向量并得到初始聚类中心，采用FCM算法对检测样本进行分类识别，结果如图6所示。

由图6可知，采用EEMD方法对转子信号进行处理后，不平衡和动静碰摩故障全部识别正确，但正常状态、不对中和轴承座松动故障的识别结果并不理想，共有5个检测样本出现了错误，平均识别率为90%，低于采用VMD方法的识别结果。究其原因，由于EEMD算法对模态混叠虽有一定抑制作用但仍无法避免，影响特征向量构建的准确性。由此表明，相较于EEMD方法，VMD方法可更为有效地提取出信号各频带的信息，保证后续方法可更加准确地提取出转子故障特征，实现转子工作状态和故障类型的有效识别。

3 结语

由于转子故障信号的非线性特征及其故障特征信息无法有效提取的问题，将VMD和Volterra模型结合，以奇异值熵构建奇异值特征向量的特征提取方法。该方法采用VMD对转子振动信号进行分解，以各IMF的能量熵增量数值选取对故障特征敏感的IMF分量对其建立Volterra模型，以获取模型参数向量组成初始特征向量，对其进行奇异值分解并获得奇异值特征向量矩阵，并采用FCM算法对转子故障类型进行分类识别。通过对实测信号的分析，表明本文所提方法可准确反映转子的故障特征，同EEMD方法相比，基于VMD的信号分解方法具有更有效的信号特征提取能力，诊断效果更好。