基于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的信號分解方法具有更有效的信號特征提取能力，診斷效果更好。