李翠省， 廖英英， 劉永強

(1.石家莊鐵道大學 交通運輸學院，石家莊 050043; 2.石家莊鐵道大學 省部共建交通工程結構力學行為與系統(tǒng)安全國家重點實驗室，石家莊 050043)

輪對軸承是列車走行部中關鍵的旋轉部件之一，起著支撐、傳動、運動轉換等重要作用，其運行狀態(tài)直接關系到列車的運行安全[1-2]。因此，對輪對軸承進行故障診斷研究具有重要的意義。

輪對軸承長期處于高速重載的運行環(huán)境中，承受著多種交變載荷的作用，同時還受到踏面磨耗、軌道不平順等外部復雜激勵的影響。所以軸承元件的工作表面極易出現(xiàn)壓痕、點蝕、剝落等局部損傷[3]。這些故障沖擊響應被強烈的背景噪聲所淹沒，通常無法直接識別。因此，在復雜干擾下準確分離與識別出故障特征，是輪對軸承故障診斷的重點和難點。

針對這一難點，學者們提出了多種分析方法，比如共振解調法，短時傅里葉變換，小波變換等，其中重要一類就是信號分解方法[4]。經(jīng)驗模態(tài)分解(empirical mode decomposition，EMD)[5]是一種自適應信號處理方法，可將非平穩(wěn)信號分解成不同頻段的模態(tài)分量，在信號處理中應用較廣，但EMD存在端點效應、模態(tài)混疊等問題[6]。集成經(jīng)驗模式分解(ensemble empirical mode decomposition, EEMD)是EMD的一種改進算法，能有效抑制模態(tài)混疊現(xiàn)象[7]。局部均值分解(local mean decomposition，LMD)[8]改善了EMD中過包絡或欠包絡的問題，但其本質與EMD、EEMD一樣都屬于遞歸式模態(tài)分解，終究避免不了端點效應和模態(tài)混疊現(xiàn)象。Dragomiretskiy等[9]根據(jù)約束性變分問題，提出了一種非遞歸式模態(tài)分解方法，即變分模態(tài)分解(variational modal decomposition，VMD)，可有效減少模態(tài)混疊現(xiàn)象的發(fā)生。該算法可以將信號在頻域內靈活分割，能有效挖掘信號中被噪聲淹沒的故障特征信息。一經(jīng)提出，受到國內外眾多學者的關注[10-11]。

學者們研究發(fā)現(xiàn)，VMD對信號的處理效果受其自身參數(shù)設置的影響。文獻[12]通過中心頻率觀察法確定分解模態(tài)數(shù)，實現(xiàn)了滾動軸承故障特征提取；文獻[13]利用奇異值最佳有效秩階次來確定分量個數(shù)，用于風電齒輪箱故障分析。上述文獻在應用VMD時，僅考慮了模態(tài)數(shù)K的影響，而忽略了二次懲罰因子α對模態(tài)帶寬的影響[14]。目前流行方式之一是通過引入智能優(yōu)化算法，來自適應地確定[K,α]參數(shù)組合。如文獻[15]采用混合灰狼進化算法，以局部最小包絡熵為目標函數(shù)優(yōu)化VMD參數(shù)；文獻[16]提出了基于相關峭度和蝗蟲優(yōu)化算法的參數(shù)自適應VMD方法。不難發(fā)現(xiàn)，基于智能優(yōu)化算法的參數(shù)自適應VMD方法的關鍵在于確定合適的目標函數(shù)，引導全局尋優(yōu)。倘若目標函數(shù)確定不恰當，將導致尋優(yōu)失敗。如上文提到的包絡熵指標，其對異常尖脈沖比較敏感，當信號中存在該項干擾時，優(yōu)化后的VMD將不能準確提取出軸承的故障特征信息。為提升目標函數(shù)的抗噪性，學者們經(jīng)過不斷探索，提出了一些先進指標。文獻[17]提出了加權峭度指標(weighted kurtosis index，KCI)，文獻[18]提出了全體峭度指標(ensemble kurtosis index，EKI)，文獻[19]提出了融合沖擊指標(syncretic impact index，SII)，均取得了不錯的效果。然而，高速列車輪對軸承工作環(huán)境復雜，振動信號中時常伴有沖擊性噪聲和循環(huán)平穩(wěn)性噪聲。筆者研究發(fā)現(xiàn)，在這些強干擾存在的環(huán)境中上述三種先進指標依然會受到影響。

通過上述分析可知，參數(shù)自適應VMD方法的核心在于確定合適的目標函數(shù)。而在強噪聲的環(huán)境下，現(xiàn)有的指標又很難避開這些干擾的影響。為此，本文提出了一種基于EEMD預處理的參數(shù)自適應VMD方法。首先將原始信號進行EEMD分解，選取包絡峭度值大于原始信號包絡峭度值的分量進行重構，生成新的信號。其次再以局部最大包絡譜峭度作為目標函數(shù)，利用粒子群優(yōu)化算法對VMD進行參數(shù)優(yōu)化。最后使用優(yōu)化后的VMD對新信號進行分解，選取包絡譜峭度值最大的分量進行包絡解調分析，可以準確地提取到輪對軸承的故障特征信息。

1 基于EEMD預處理的參數(shù)自適應VMD方法

1.1 集成經(jīng)驗模態(tài)分解(EEMD)

EEMD是一種基于信號本身時間尺度特征進行自適應分解的方法，其本質是在原信號中加入高斯白噪聲，利用白噪聲頻率均勻分布的統(tǒng)計特性，經(jīng)多次EMD分解，達到抑制模態(tài)混疊的目的。

EEMD具體分解過程如下[20]

(1) 在信號x(t)中加入隨機高斯白噪聲

xm(t)=x(t)+knm(t)

(1)

式中,k為所加噪聲幅值系數(shù)。

(2) 利用EMD方法對信號xm(t)進行分解；

(3) 每次加入不同的高斯白噪聲序列，重復步驟(1)、(2)，直至達到終止條件；

(4) 將分解得到的各個IMF分量的均值記為最終分解結果

(2)

式中：N為加噪聲的次數(shù);ci,m表示第m次EMD分解得到的第i個IMF分量。

1.2 VMD算法

VMD是一種完全非遞歸的信號分解模式，其實質是構造一個受約束的變分模型，通過迭代搜尋最優(yōu)解。VMD可將一個非平穩(wěn)的輸入信號y分解成K個具有特定稀疏性的模態(tài)分量yk，并能確定出每個分量的中心頻率ωk和帶寬。

算法中通過計算調制信號梯度的L2范數(shù)來估計各分量的帶寬，其對應的約束變分模型表達式如下

(3)

式中：{yk}表示分量信號;{ωk}表示相應的中心頻率的集合;K為分解層數(shù)。

引入二次懲罰因子和Lagrange乘子將式(3)轉化為無約束優(yōu)化模型，表達式如下

(4)

式中：α為二次懲罰因子;λ為Lagrange乘子。

對式(4)采用交替乘子算法求其最優(yōu)解，結果如下

(5)

(6)

(3) 根據(jù)式(7)更新

(7)

(4) 根據(jù)式(8)判斷是否收斂，若滿足收斂條件分解過程結束；否則迭代次數(shù)n=n+1，返回步驟(2)繼續(xù)分解。

(8)

1.3 改進粒子群優(yōu)化算法

粒子群算法擁有進化計算和群體智能的特點，具有良好的全局優(yōu)化能力。在一個D維空間中，由m個粒子組成的種群Z={Z1，Z2，…，Zm}，每個粒子所處的位置Zi={zi1,zi2，…，zin}都表示一個解，且會根據(jù)目標函數(shù)搜索新解，每個粒子的速度為Vi={vi1,vi2,…,vin}。每次迭代中，粒子跟蹤個體最好解Pid和種群最優(yōu)解Pgd，根據(jù)式(9)來更新自身速度和位置

(9)

式中：ω為慣性權重;d=1,2,…,D;i=1,2,…,m;k為當前跌迭代次數(shù);c1和c2為學習因子;η為介于[0,1]間的隨機數(shù)。

在標準粒子群算法中，隨著粒子速度的不斷更新會對ω進行動態(tài)調整，用于權衡全局和局部搜索能力，線性遞減權值策略如下

(10)

式中：Tmax為最大迭代次數(shù)；t表示當前迭代次數(shù)；ωmax為最大慣性權重；ωmin為最小慣性權重。

粒子群的性能在很大程度上依賴于參數(shù)ω、c1和c2的設置，文獻[21]在上述標準粒子群的基礎上，提出了動態(tài)學習因子，成功用于VMD參數(shù)優(yōu)化，公式如下

(11)

式中,R1、R2、R3和R4為初始系數(shù)值。

文獻[22]對慣性權重進行了深入研究，發(fā)現(xiàn)凹函數(shù)遞減策略在不影響收斂精度的前提下，能較大幅度地提升粒子群的收斂速度，整體效果優(yōu)于線性遞減策略。因此，為進一步提升粒子群的優(yōu)化性能，本文將凹函數(shù)遞減權值代替式(10)，構成新的改進粒子群優(yōu)化算法，具體公式如下

(12)

根據(jù)文獻[23],設置改進PSO的具體參數(shù)，如表1所示。

表1 粒子群參數(shù)

1.4 故障指標性能分析

滾動軸承的故障信號在時域具有沖擊性，在頻域具有循環(huán)平穩(wěn)性，在進行故障特征提取時需同時考慮這兩大特征。顧曉輝等指出振動信號包絡的稀疏性可表征軸承故障的沖擊性，信號包絡譜的稀疏性可表征軸承故障的循環(huán)平穩(wěn)性，對于稀疏性的度量，峭度為最常用指標之一。

信號包絡的峭度如式(13)所示

(13)

式中：Ex為信號經(jīng)Hilbert解調后得到的包絡信號;μE為Ex的均值;σE為Ex的標準差。

信號包絡譜的峭度如式(14)所示

(14)

式中：ESx=DFT[Ex];μES為ESx的均值;σES為ESx的標準差。

在對軸承故障進行評價時，單一指標易受沖擊性噪聲或循環(huán)平穩(wěn)性噪聲的影響。為此，學者們經(jīng)過不懈努力，提出了一些新的復合型指標。Zhang等提出了加權峭度指標IKC，公式如下

IKC=KI·|C|

(15)

式中：KI代表信號的峭度;|C|代表相關系數(shù)。

Miao等提出了全體峭度指標IEK，公式如下

IEK=KES·K

(16)

式中：K代表信號的峭度;KES代表信號的包絡譜峭度。

He等提出了融合沖擊指標ISI，公式如下

ISI=KEPS·|CC|

(17)

式中：KEPS代表信號的包絡功率譜峭度；|CC|代表相關系數(shù)。

為了研究這些指標的性能，用第3章的仿真信號進行分析。仿真信號是由故障沖擊、高幅值隨機脈沖、周期性干擾和高斯白噪聲組成。用上述五種指標對這四種成分進行敏感性分析，如圖1所示。我們不難發(fā)現(xiàn)，KE和IKC對高幅值隨機脈沖敏感，KES、IEK和ISI對周期性干擾敏感。在強干擾環(huán)境下，現(xiàn)有的單一或復合型指標都很難統(tǒng)一并平衡軸承故障響應的沖擊性和循環(huán)平穩(wěn)性特征。

圖1 五種指標對不同成分的敏感性Fig.1 Sensitivity of five indexes to different components

2 方法步驟

參數(shù)自適應VMD方法有很強的故障特征提取能力，但診斷效果嚴重依賴于優(yōu)化算法中目標函數(shù)的設置。在沖擊性噪聲和循環(huán)平穩(wěn)性噪聲存在的環(huán)境下，現(xiàn)有的故障評價指標又很難避開這些干擾的影響。針對該問題，本文同時考慮到EEMD具有優(yōu)秀的降噪能力和自適應性[24]，為此提出了基于EEMD預處理的參數(shù)自適應VMD方法。首先將原信號進行EEMD分解，以包絡峭度為評價指標，進行信號重構，可排除循環(huán)平穩(wěn)性噪聲的影響，但新信號中仍會存在沖擊性噪聲；其次再以包絡譜峭度為目標函數(shù)，優(yōu)化VMD參數(shù)，可消除沖擊性噪聲的干擾，最終達到精準提取故障特征信息的目的。詳細步驟如下：

步驟1將采集到的振動信號進行EEMD分解，計算原信號以及各分量的包絡峭度值，選取峭度值大于原始信號峭度值的分量進行重構，生成新的振動信號。

步驟2給定VMD中二次懲罰因子α的尋優(yōu)范圍[100,5 000]，模態(tài)分解個數(shù)K的尋優(yōu)范圍[2,10]。

步驟3根據(jù)表1設置PSO中的參數(shù)，以局部最大包絡譜峭度為目標函數(shù)，搜尋VMD的最佳參數(shù)組合。

步驟4利用優(yōu)化后的VMD對新信號進行分解，生成K個分量，選取包絡譜峭度值最大的分量，進行包絡解調分析，與理論故障特征頻率對比，判斷軸承故障。

具體診斷流程如圖2所示。

圖2 輪對軸承故障診斷流程圖Fig.2 The flow chart of wheelset bearing fault diagnosis

3 仿真信號分析

為了驗證本文方法的正確性，采用文獻[25]中介紹的模擬信號模型，來模擬復雜工況下真實的振動信號。信號中包含四部分：軸承故障引起的重復性瞬態(tài)，電磁干擾引起的高幅值隨機脈沖，軸旋轉和齒輪嚙合產生的周期性干擾，以及高斯白噪聲。公式如下

(18)

式中：fr代表轉頻;Z代表齒輪的齒數(shù);Td為故障周期;τi為滾動體的隨機滾動誤差;S(t)為脈沖響應函數(shù)。S(t)的表達式如下

(19)

在式(18)和(19)中用到的具體參數(shù)如表2所示，n(t)為均值為0，方差為0.25的高斯白噪聲。

表2 仿真信號中的參數(shù)

采樣頻率Fs為10 000 Hz,采樣時間為1 s。仿真信號各組成成分如圖3所示，圖3(a)為軸承故障脈沖信號，圖3(b)為高幅值隨機脈沖信號，圖3(c)為周期性干擾，圖3(d)為高斯白噪聲。

(a) 故障沖擊

(b) 隨機沖擊

(c) 周期干擾

(d) 高斯白噪聲圖3 仿真信號的四種成分Fig.3 The four components of the simulation signal

仿真信號的時域波形、頻譜以及包絡譜如圖4所示，在時域波形圖中可以看出故障沖擊信號已經(jīng)被強噪聲完全掩蓋，其中高幅值隨機脈沖比較明顯；頻譜中轉頻以及齒輪嚙合頻率譜線比較突出；在包絡譜中難以辨別出故障特征頻率，表明傳統(tǒng)的譜分析方法失效。

(a) 時域波形

(b) 頻譜

(c) 包絡譜圖4 仿真信號Fig.4 The simulation signal

應用本文所提方法進行分析。首先，對原始信號進行EEMD分解，生成12個IMF分量和1個余量r。計算原始信號以及各分量的包絡峭度值，如圖5所示，選取IMF1和IMF6進行重構，生成新的信號，其時域波形和包絡譜如圖6所示?？梢园l(fā)現(xiàn)，循環(huán)平穩(wěn)性噪聲被有效地剔除了，沖擊性噪聲依然存在，包絡譜中可以找到故障特征頻率的基頻。

圖5 各分量以及原始信號的包絡峭度值Fig.5 Envelope kurtosis values of each component and original signal

(a) 時域波形

(b) 包絡譜圖6 新信號的時域波形以及包絡譜Fig.6 The time domain waveform and envelope spectrum of the new signal

其次，以局部最大包絡譜峭度為目標函數(shù)，利用PSO對VMD進行參數(shù)優(yōu)化，用優(yōu)化后的VMD對新信號分解，選取包絡譜峭度值最大的分量進行包絡解調分析，結果如圖7所示。在時域波形上可以看到清晰的故障沖擊，高幅值隨機脈沖被有效剔除；在包絡譜上可以清晰地發(fā)現(xiàn)故障特征頻率的基頻及其倍頻,譜線都很突出。結果表明，本文方法可以在復雜干擾下，準確提取出軸承的故障特征信息。

(a) 時域波形

(b) 包絡譜圖7 本文所提方法的最終結果Fig.7 The final result of the method presented in this paper

為說明本文方法的優(yōu)越性，直接對原始信號進行參數(shù)自適應VMD方法分析。以KE為評價指標，最優(yōu)分量如圖8所示，只能發(fā)現(xiàn)高幅值沖擊，找不到任何故障信息；以KES為評價指標，最優(yōu)分量如圖9所示，搜尋到的是齒輪嚙合信息。用目前提出的先進復合指標進行尋優(yōu)，分別以EKI、KCI和SII為評價指標，最終分析結果如圖10、圖11和圖12所示，由于受沖擊性噪聲或循環(huán)平穩(wěn)性噪聲的影響，這三種指標引導的參數(shù)自適應VMD方法依然不能找到軸承的故障信息。

(a) 時域波形

(b) 包絡譜圖8 以KE為目標函數(shù)的VMD分析結果Fig.8 VMD analysis results with KE as the objective function

(a) 時域波形

(b) 包絡譜圖9 以KES為目標函數(shù)的VMD分析結果Fig.9 VMD analysis results with KES as the objective function

(a) 時域波形

(b) 包絡譜圖10 以EKI為目標函數(shù)的VMD分析結果Fig.10 VMD analysis results with EKI as the objective function

(a) 時域波形

(b) 包絡譜圖11 以KCI為目標函數(shù)的VMD分析結果Fig.11 VMD analysis results with KCI as the objective function

(a) 時域波形

(b) 包絡譜圖12 以SII為目標函數(shù)的VMD分析結果Fig.12 VMD analysis results with SII as the objective function

4 試驗信號分析

為進一步驗證本文所提方法的有效性，將其應用于高速列車輪對軸承的故障診斷分析。鐵路軸承綜合試驗臺以及測試軸承如圖13所示。測試軸承主要結構參數(shù)如表3所示，軸承的外圈上存在長5 mm,寬1 mm,深0.7 mm的損傷。與它的幾何形狀相比，缺陷尺寸相對較小。

(a) 鐵路軸承綜合試驗臺

(b) 鐵路軸承綜合試驗臺結構圖

(c) 測試軸承圖13 鐵路軸承綜合試驗臺及故障軸承Fig.13 Railway bearing comprehensive test-bed and fault bearing

表3 測試軸承主要參數(shù)

試驗時，為模擬輪對軸承真實的運行工況，施加平均徑向力為80 kN，最大軸向力為±40 kN，激勵頻率為0.2～20 Hz的動態(tài)荷載。

采樣頻率為51.2 kHz，轉頻為28.5 Hz。根據(jù)軸承尺寸和試驗工況，計算得到外圈故障頻率fBPFO=208.1 Hz。應用加速度傳感器采集60 s的振動數(shù)據(jù)，截取其中0.5 s的數(shù)據(jù)進行分析，其時域波形和包絡譜如圖14所示。由于受到背景噪聲的干擾，從時域信號中無法找到較明顯的故障沖擊，包絡譜中亦看不到故障特征頻率譜線，表明傳統(tǒng)的譜分析方法失效。

應用本文所提方法對試驗數(shù)據(jù)進行分析。首先，對原始信號進行EEMD分解，生成13個IMF分量和1個余量r。計算原始信號以及各分量的包絡峭度值，如圖15所示，選取IMF1、IMF2、IMF4和IMF5進行重構，生成新的信號，其時域波形和包絡譜如圖16所示。從時域圖上可以看出噪聲含量有所降低，但包絡譜上依然找不到故障特征頻率的譜線。

(a) 時域波形

(b) 包絡譜圖14 試驗信號Fig.14 Test signal

圖15 各分量以及原始信號的包絡峭度值Fig.15 Envelope kurtosis values of each component and original signal

(a) 時域波形

(b) 包絡譜圖16 新信號的時域波形以及包絡譜Fig.16 The time domain waveform and envelope spectrum of the new signal

其次，以局部最大包絡譜峭度為目標函數(shù)，利用PSO對VMD進行參數(shù)優(yōu)化，用優(yōu)化后的VMD對新信號分解，選取包絡譜峭度值最大的分量進行包絡解調分析，結果如圖17所示。在時域波形上可以看到清晰的故障沖擊；在包絡譜上可以清晰地發(fā)現(xiàn)故障特征頻率的基頻及其二倍頻,譜線都很突出。結果表明，本文方法可以在復雜干擾下，準確提取出輪對軸承的故障特征信息。

(a) 時域波形

(b) 包絡譜圖17 本文所提方法的最終結果Fig.17 The final result of the method presented in this paper

為說明本文方法的優(yōu)越性，直接對原始信號進行參數(shù)自適應VMD方法分析。以KE為評價指標，最優(yōu)分量如圖18所示，時域波形上有高幅值沖擊，包絡譜上找不到任何外圈故障信息；以KES為評價指標，最優(yōu)分量如圖19所示，包絡譜上亦發(fā)現(xiàn)不了外圈故障信息。用目前提出的先進復合指標進行尋優(yōu)，分別以EKI、KCI和SII為評價指標，最終分析結果如圖20、圖21和圖22所示，由于受到強噪聲的干擾，這三種指標引導的VMD依然不能找到外圈的故障信息。

(a) 時域波形

(b) 包絡譜圖18 以KE為目標函數(shù)的VMD分析結果Fig.18 VMD analysis results with KE as the objective function

(a) 時域波形

(b) 包絡譜圖19 以KES為目標函數(shù)的VMD分析結果Fig.19 VMD analysis results with KES as the objective function

(a) 時域波形

(b) 包絡譜圖20 以EKI為目標函數(shù)的VMD分析結果Fig.20 VMD analysis results with EKI as the objective function

(a) 時域波形

(b) 包絡譜圖21 以KCI為目標函數(shù)的VMD分析結果Fig.21 VMD analysis results with KCI as the objective function

(a) 時域波形

(b) 包絡譜圖22 以SII為目標函數(shù)的VMD分析結果Fig.22 VMD analysis results with SII as the objective function

5 結 論

(1) 本文在文獻[21]使用的粒子群基礎上，將凹函數(shù)遞減權值代替線性遞減權值，進一步提升了PSO的全局搜索能力。

(2) VMD能夠將信號在頻域內進行精細分析，但效果受到K和α的影響。利用智能優(yōu)化算法可自適應確定這兩個關鍵參數(shù)，但優(yōu)化效果嚴重依賴于目標函數(shù)。在強噪聲干擾下，現(xiàn)有作為目標函數(shù)的故障評價指標又易受沖擊性噪聲或循環(huán)平穩(wěn)性噪聲的影響。

(3) 針對傳統(tǒng)參數(shù)自適應VMD方法易受噪聲影響的問題，提出了基于EEMD預處理的參數(shù)自適應VMD方法。首先信號經(jīng)EEMD分解，以包絡峭度為指標，可以消除循環(huán)平穩(wěn)性噪聲的影響。其次利用VMD優(yōu)秀的故障特征提取能力，以包絡譜峭度為指標可以避開沖擊性噪聲的影響。通過仿真和試驗分析，驗證了本文所提方法的有效性和魯棒性。結果表明，該種方法在強干擾環(huán)境下依然能夠實現(xiàn)輪對軸承故障的精準診斷，效果優(yōu)于傳統(tǒng)參數(shù)自適應VMD方法。