張毅,郝高巖,劉璇
(1.北京市軌道交通建設管理有限公司,北京 100068;2.城市軌道交通全自動運行系統與安全監控北京市重點實驗室, 北京 100068;3.北京博華信智科技股份有限公司,北京 100029;4.洛陽軸承研究所有限公司,河南 洛陽 471039)
自動扶梯電動機、齒輪箱、驅動軸等關鍵部件的軸承故障時有發生,屬于典型的旋轉機械類故障,通常在關鍵部位安裝振動傳感器獲取振動數據并通過帶通濾波和包絡解調等方法識別故障[1-2]。對于軸承等關鍵部件的早期故障,故障沖擊特征易被淹沒在其他高幅值的設備運轉噪聲中,且由于信號的隨機性和非平穩性,常規的帶通濾波和包絡解調方法僅為近似平均方法,無法真實刻畫實際信號特征。
目前,信號的自適應分解技術由于擺脫了傅里葉分析、小波分析[3]等固定基函數的缺陷,可針對原始非平穩信號自適應選擇分解模式,經驗模態分解(EMD)[4-5]方法的分解效果高度依賴極值點的查找、極值點的擬合方法以及迭代的結束條件,且缺乏一定的數學理論基礎,這些條件決定了EMD方法的魯棒性不佳。
變分模態分解(Variational Mode Decomposition,VMD)方法克服了EMD方法的不足[7],在軸承故障診斷中應用廣泛:文獻[8-12]利用包絡熵對變分模態分解參數進行優化,但包絡熵計算復雜,且沒有考慮信號相關性;文獻[13]對分解層數和懲罰因子的影響關系進行分析,利用分量信號的能量之和占原信號能量的比值判斷分解效果,從而找出變分模態分解的最佳分解層數,但其無法準確確認懲罰因子;文獻[14]僅利用峭度指標作為尋優適應度函數,導致出現了過分解問題。綜上分析,本文提出新的適應度函數尋優算法對變分模態分解的分解層數和懲罰因子進行尋優,并通過自動扶梯電機軸承的真實故障數據進行驗證分析。
變分模態分解可理解為將原信號分解為具有一定中心頻率和帶寬的子信號,其主要步驟為[9]:
1)對每一個模態分量uk,求得實信號的解析信號。
2)將uk與具有中心頻率ωk的復指數信號相乘,把各模態分量的頻率搬移到基帶。
3)通過調制信號的高斯平滑度方法估計各模態分量的帶寬,即梯度的L2范數。
以上過程轉換為約束變分問題,即
(1)
式中:uk為分解得到的第k個模態分量;ωk為第k個模態分量的中心頻率;δ(t)為脈沖函數;f為原信號。

(2)
(3)
(4)
式中:α為帶寬懲罰因子;τ為雙上升步長;λ為拉格朗日乘子;ε為誤差因子,f(ω)為信號頻域。
重復(2)—(4)式進行更新,滿足迭代停止條件時結束循環[7],得到K個IMF分量。
蜻蜓優化算法(Dragonfly Algorithm,DA)是一種新型智能優化算法[15],該算法模擬蜻蜓的生活習性并將其歸納為5類行為方式:分離、對齊、內聚、尋找食物和遠離天敵。
根據蜻蜓行為方式進行蜻蜓位置的更新迭代[16],計算公式為
Xt+1=Xt+ΔXt+1,
(5)
ΔXt+1=(sSi+aAi+cCi+fFi+eEi)+wΔXt,
式中:ΔXt+1為位置更新步長;t為當前迭代次數;s為分離度權重;Si為第i個個體的分離度;a為對齊度權重;Ai為第i個個體的對齊度;c為內聚度權重;Ci為第i個個體的內聚度;f為食物因子;Fi為食物位置對第i個個體的吸引力;e為天敵因子;Ei為天敵位置對第i個個體的排斥力;w為慣性權重。
引入Lévy飛行隨機游走方法進行蜻蜓位置更新,即
Xt+1=Xt+Lévy(d)ΔXt,
(6)

(7)
Γ(x)=(x-1)!,
式中:d為維度;r1,r2為[0,1]范圍內的隨機數;β為常數,此處取值為1.5。
變分模態分解的實際使用中,需要人為確定分解層數K和懲罰因子α,因此,本文結合蜻蜓優化方法構造新的適應度函數,提出自適應的變分模態分解方法,主要步驟為:
1)對蜻蜓優化算法參數進行初始化,初始化群體個數為20,迭代次數為10;蜻蜓位置初始化,初始化K∈[2,10],α∈[500,10 000],兩者在該范圍內隨機取值;初始化位置變化步長。
2)計算各權重值,根據初始值隨機初始化初始參數s,f,w,a,c,e。
3)采用相關系數與峭度的乘積構造適應度函數,即
Ti=max(ciki);s.t.ki>3.2,
(8)
式中:ci為第i個分量的相關系數;ki為第i個分量的峭度,為避免相關性很大而信噪比較小的情況,僅考慮峭度大于3.2的分量。
4)依據(6)式更新蜻蜓位置,迭代蜻蜓位置及步長,滿足最大迭代次數后停止程序。
5)采用最優的K和α進行變分模態分解,對分解后適應度最大的分量進行希爾伯特包絡分析并計算包絡波形的頻譜,觀察包絡譜中是否存在明顯的軸承故障特征頻率。
采用凱斯西儲大學軸承故障數據進行算法驗證,使用1 750 r/min轉速工況下驅動端6205-2RS深溝球軸承的數據,其外圈故障特征頻率為104.7 Hz,內圈故障特征頻率為157.5 Hz。為驗證算法在實際工況下的有效性,分別在外圈、內圈原始信號中加入-15,-10 dB的白噪聲信號。
軸承外圈故障振動信號的時域波形及其包絡譜如圖1所示,故障特征頻率被淹沒在噪聲中,包絡譜中也無法有效提取軸承故障特征頻率。

圖1 軸承外圈故障振動信號的時域波形及其包絡譜Fig.1 Time domain waveform and envelope spectrum of bearing outer ring fault vibration signal
對軸承外圈故障振動信號進行變分模態分解,采用新構造的適應度函數進行蜻蜓尋優,確定分解參數K=6,α=5 200,分解所得IMF分量如圖2所示,各分量的相關系數、峭度及其乘積見表1。分析可知IMF4分量的適應度值最大,且其時域波形中有一些沖擊成分,因此選擇IMF4分量進行包絡解調,從圖3所示包絡譜中可明顯觀察到軸承外圈故障特征頻率及其2倍頻。

圖2 軸承外圈故障振動信號變分模態分解所得分量(K=6,α=5 200)Fig.2 VMD components of bearing outer ring fault vibration signal(K=6,α=5 200)

表1 軸承外圈故障振動信號各IMF分量的相關系數及峭度(K=6,α=5 200)Tab.1 Correlation coefficient and kurtosis of IMF of bearing outer ring fault vibration signal(K=6,α=5 200)

圖3 軸承外圈故障振動信號IMF4分量的包絡譜(K=6,α=5 200)Fig.3 IMF4 envelope spectrum of bearing outer ring fault vibration signal(K=6,α=5 200)
軸承內圈故障振動信號的時域波形及其包絡譜如圖4所示,同樣無法發現軸承內圈故障特征頻率。采用新構造的適應度函數進行蜻蜓尋優可得K=4,α=8 820,變分模態分解所得IMF分量如圖5所示,各分量的相關系數、峭度及其乘積見表2。分析可知IMF3分量的適應度值最大,其包絡譜如圖6所示,可明顯觀察到軸承內圈故障特征頻率及其2倍頻(315.7 Hz)、3倍頻(473.1 Hz)。

圖4 軸承內圈故障振動信號的時域波形及其包絡譜Fig.4 Time domain waveform and envelope spectrum of bearing inner ring fault vibration signal

圖5 軸承內圈故障振動信號變分模態分解所得分量(K=4,α=8 820)Fig.5 VMD components of bearing inner ring fault vibration signal(K=4,α=8 820)

表2 軸承內圈故障振動信號各IMF分量的相關系數及峭度(K=4,α=8 820)Tab.2 Correlation coefficient and kurtosis of IMF of bearing inner ring fault vibration signal(K=4,α=8 820)

圖6 軸承內圈故障振動信號IMF3分量的包絡譜Fig.6 IMF3 envelope spectrum of bearing inner ring fault vibration signal
為驗證本文所提適應度函數的優越性,僅采用峭度作為適應度進行尋優,得到外圈故障(K=9,α=4 880)、內圈故障(K=9,α=5 500)的最佳參數并選擇峭度最大的分量進行包絡譜分析,結果如圖7所示。由圖7可知:可以觀察到軸承外圈故障特征頻率104.7 Hz但幅值相對較低,而其2倍頻則被噪聲淹沒;軸承內圈故障特征頻率157.5Hz相對明顯,但同樣觀察不到其倍頻特征;本文所提適應度函數兼顧了軸承故障沖擊信號特征及與原信號的相關性,能更好的反映軸承的故障特征。

圖7 峭度作為適應度函數時的故障特征提取結果Fig.7 Fault feature extraction results of kurtosis as fitness function
對某自動扶梯電機軸承故障進行診斷分析,自動扶梯類型為單驅動,運行速度0.51 m/s,電動機實際轉速778 r/min,轉頻12.97 Hz。軸承型號為6310,在運行轉速下軸承外圈故障特征頻率為39.516 Hz,內圈故障特征頻率為64.217 Hz。
采集電機軸承機殼處的振動加速度,采樣頻率12 800 Hz,采樣點數16 384,所采集振動信號的時域波形及其頻譜、包絡譜如圖8所示,從圖中無法有效判別軸承狀態,故障特征頻率淹沒在噪聲和變頻器干擾(78 Hz為間隔的高頻成分)中。

圖8 自動扶梯電機軸承振動信號的時域波形及其頻譜、包絡譜Fig.8 Time domain waveform,spectrum and envelope spectrum of vibration signal of escalator motor bearing
采用改進的變分模態分解方法,以新構造的適應度函數進行蜻蜓尋優獲得最優的分解參數為K=4,α=5 800。采用最優參數進行變分模態分解,所得分量如圖9所示,各分量的相關系數、峭度及其乘積見表3。分析可知:IMF2分量的適應度值最大且其時域波形中具有明顯的周期性軸承故障沖擊特征。

圖9 自動扶梯電機軸承振動信號變分模態分解所得分量Fig.9 VMD components of vibration signal of escalator motor bearing

表3 自動扶梯電機軸承各IMF分量的相關系數及峭度Tab.3 Correlation coefficient and kurtosis of IMF of vibration signal of escalator motor bearing
對IMF2分量進行希爾伯特包絡解調分析,結果如圖10所示,從包絡譜中可明顯觀察到40.63,80.47,121.1,160.9,201.6 Hz頻率處存在明顯譜峰,與自動扶梯電機軸承外圈故障特征頻率及其倍頻相近,因此判定自動扶梯電機軸承存在外圈故障,與拆機檢查結果一致,說明本文方法有效可行。

圖10 自動扶梯電機軸承IMF2分量的包絡譜Fig.10 IMF2 envelope spectrum of vibration signal of escalator motor bearing
結合蜻蜓優化算法,以相關系數與峭度的乘積作為適應度函數,對變分模態分解的分解層數和懲罰因子進行尋優,新構造的適應度函數兼顧了峭度指標的沖擊敏感性和IMF分量與原信號的相關性,避免了僅通過峭度為適應度函數易出現的過分解問題。
試驗臺數據以及自動扶梯電機軸承的分析結果表明,改進變分模態分解與包絡譜分析相結合的方法能夠在低信噪比情況下充分提取軸承故障特征信息,在軸承故障診斷中具有較好的應用前景。