陳 鵬, 趙小強(qiáng),2,3

(1. 蘭州理工大學(xué) 電氣工程與信息工程學(xué)院， 蘭州 730050; 2. 甘肅省工業(yè)過程先進(jìn)控制重點(diǎn)實(shí)驗(yàn)室, 蘭州 730050) 3. 蘭州理工大學(xué) 國家級(jí)電氣與控制工程實(shí)驗(yàn)教學(xué)中心, 蘭州 730050)

滾動(dòng)軸承作為旋轉(zhuǎn)類機(jī)械設(shè)備的關(guān)鍵零部件，其運(yùn)行狀態(tài)直接關(guān)系到整個(gè)設(shè)備的安全狀況[1]。在其旋轉(zhuǎn)過程中，由于與其它元件接觸摩擦引起損傷，當(dāng)反復(fù)經(jīng)過損傷點(diǎn)時(shí)，采集的振動(dòng)信號(hào)會(huì)產(chǎn)生周期性脈沖[2]，該脈沖是判斷滾動(dòng)軸承發(fā)生故障的主要依據(jù)。但是軸承的運(yùn)行環(huán)境復(fù)雜，常受到傳輸路徑和強(qiáng)噪聲的干擾，導(dǎo)致早期故障的微弱信號(hào)淹沒在噪聲中而難以提取[3]。因此，如何有效地分析滾動(dòng)軸承振動(dòng)信號(hào)，提取出故障引起的沖擊特征進(jìn)行故障辨識(shí)是整個(gè)故障診斷的關(guān)鍵。

對(duì)于滾動(dòng)軸承的振動(dòng)信號(hào)分析而言，應(yīng)用信號(hào)處理方法對(duì)振動(dòng)信號(hào)進(jìn)行降噪的意義是顯而易見的。近年來，在小波降噪的基礎(chǔ)上，學(xué)者們提出了許多自適應(yīng)信號(hào)處理降噪方法，并應(yīng)用于滾動(dòng)軸承的故障特征提取[4]。如經(jīng)驗(yàn)?zāi)B(tài)分解(empirical mode decomposition，EMD)[5]、集總經(jīng)驗(yàn)?zāi)B(tài)分解(ensemble empirical mode decomposition， EEMD)[6]、互補(bǔ)集合經(jīng)驗(yàn)?zāi)B(tài)分解(complete ensemble empirical mode decomposition,CEEMD)[7]、局部均值分解(local mean decomposition,LMD)[8]、局部特征尺度分解(local characteristic-scale decomposition， LCD)[9]和固有時(shí)間尺度分解(intrinsic time-scale decomposition，ITD)[10]等。這些方法都在不同程度上實(shí)現(xiàn)了滾動(dòng)軸承振動(dòng)信號(hào)的降噪，但是這些方法都是基于遞歸分解原理，沒有從本質(zhì)上解決其模態(tài)混疊和端點(diǎn)效應(yīng)的缺點(diǎn)。而Dragomiretskiy等[11]提出了一種新的信號(hào)自適應(yīng)分解方法——變分模態(tài)分解(variational modal decomposition, VMD)。由于其通過引入變分約束來實(shí)現(xiàn)信號(hào)的分解，克服了上述方法存在的模態(tài)混疊和端點(diǎn)效應(yīng)問題，被廣泛應(yīng)用于旋轉(zhuǎn)機(jī)械的振動(dòng)信號(hào)降噪[12]。然而VMD存在模態(tài)分解個(gè)數(shù)k和二次罰因子α需要進(jìn)行預(yù)設(shè)置且不同參數(shù)值對(duì)于分解效果影響大的問題。對(duì)于此問題，如通過手動(dòng)調(diào)節(jié)較難獲得最優(yōu)參數(shù)。因此，如何自適應(yīng)地實(shí)現(xiàn)VMD最優(yōu)參數(shù)選擇，對(duì)于信號(hào)的分解效果和后續(xù)特征提取至關(guān)重要。對(duì)此，王奉濤等[13]提出了基于能量的模態(tài)個(gè)數(shù)k選擇；Lian等[14]提出了預(yù)設(shè)置k的范圍，在設(shè)置的范圍內(nèi)分解信號(hào)，計(jì)算不同分解模態(tài)個(gè)數(shù)下的能量損失來獲得最佳k。兩位學(xué)者采取預(yù)設(shè)一定的范圍值，對(duì)k進(jìn)行了優(yōu)化，自適應(yīng)性不夠。唐貴基等[15]利用粒子群優(yōu)化算法優(yōu)化組合參數(shù)，建立了基于最小包絡(luò)熵的適應(yīng)度函數(shù)；劉建昌等[16]詳細(xì)分析了組合參數(shù)對(duì)VMD分解效果的影響，并通過遺傳變異粒子群算法對(duì)其進(jìn)行優(yōu)化，適應(yīng)度函數(shù)采用了唐貴基等人提出的最小包絡(luò)熵；Yan等[17]提出了基于最小包絡(luò)熵的CSA優(yōu)化VMD分解方法。學(xué)者們通過智能優(yōu)化算法實(shí)現(xiàn)VMD組合參數(shù)的自適應(yīng)選擇，建立了能夠表征信號(hào)稀疏性的包絡(luò)熵適應(yīng)度函數(shù)，實(shí)現(xiàn)了信號(hào)的分解降噪。可見，通過選擇全局和局部搜索能力較強(qiáng)的優(yōu)化算法，建立合理的適應(yīng)度函數(shù)，可實(shí)現(xiàn)VMD參數(shù)的自適應(yīng)選擇分解降噪而獲得較好的效果。

Mirjalili等[18]在2016年提出了WOA算法，證明具有較強(qiáng)的局部和全局搜索能力。峭度由于可以檢測脈沖的瞬變性而應(yīng)用于選擇模態(tài)分量，但是對(duì)異常噪聲敏感、穩(wěn)定性不高。L-峭度[19]是一種描述信號(hào)分布特性的指標(biāo)，是不同階次統(tǒng)計(jì)的線性組合，相比峭度只是樣本的四階統(tǒng)計(jì)量而言，穩(wěn)定性更高，能夠消除異常值和噪聲信號(hào)的影響，可以更好的檢測周期性脈沖[20]。相關(guān)系數(shù)可以表達(dá)模態(tài)分量和原始信號(hào)的相關(guān)性，但對(duì)噪聲敏感。

基于以上分析，本文構(gòu)建L-峭度與相關(guān)系數(shù)的自適應(yīng)度函數(shù)，利用WOA優(yōu)化VMD參數(shù)實(shí)現(xiàn)振動(dòng)信號(hào)的自適分解降噪，通過構(gòu)建的準(zhǔn)則選取最優(yōu)模態(tài)分量，再通過改進(jìn)閾值降噪技術(shù)對(duì)選取的最優(yōu)模態(tài)分量進(jìn)行進(jìn)一步降噪，實(shí)現(xiàn)振動(dòng)信號(hào)故障特征的有效提取和故障辨識(shí)。

1 VMD算法

VMD的本質(zhì)是通過迭代搜尋變分模型，把信號(hào)x(t)分解為具有一定稀疏性的k個(gè)模態(tài)分量，每個(gè)模態(tài)分量以中心頻率ωk進(jìn)行分解，最終得到約束變分問題為

(1)

式中:{uk}={u1,…,uk}為分解后的k個(gè)模態(tài)分量；{ωk}={ω1,…,ωk}為分解后每個(gè)模態(tài)分量的頻率中心。為了求解式(1)，利用二次罰函數(shù)項(xiàng)和Lagrange乘子將式(1)轉(zhuǎn)化為如下形式的無約束問題

L({uk},{ωk},λ)=

(2)

式中:α為懲罰因子；λ(t)為Lagrange乘子。

(3)

VMD的具體實(shí)現(xiàn)過程如下：

步驟2重復(fù)n←n+1；

(4)

更新ωk

(5)

對(duì)于所有ω，雙上升步長

(6)

直至滿足迭代終止條件

(7)

式中,ε為判別精度，ε>0。

2 優(yōu)化VMD參數(shù)與改進(jìn)閾值降噪

2.1 構(gòu)建L-峭度和相關(guān)系數(shù)準(zhǔn)則

傳統(tǒng)的均值、偏差和峭度等統(tǒng)計(jì)量被用于選擇有效模態(tài)，這些統(tǒng)計(jì)量識(shí)別異常的方法只是剔除標(biāo)準(zhǔn)差和均值之間偏離甚遠(yuǎn)的數(shù)據(jù)，這對(duì)于小樣本數(shù)據(jù)而言是不準(zhǔn)確的，而恰恰旋轉(zhuǎn)機(jī)械設(shè)備的故障數(shù)據(jù)都屬于小樣本數(shù)據(jù)，因此常用的以上統(tǒng)計(jì)量影響樣本數(shù)據(jù)統(tǒng)計(jì)的正確性。而L-峭度是一種描述信號(hào)分布特性的指標(biāo)，是不同階次統(tǒng)計(jì)的線性組合，相比峭度只是樣本的四階統(tǒng)計(jì)量而言，穩(wěn)定性更高，能夠消除異常值和噪聲信號(hào)的影響，可以用來很好地確定最佳模態(tài)分量。L-峭度的表達(dá)式為式(8)

(8)

如果信號(hào)的分布形狀扁平，則λ2的值將趨于很大的值，這意味著分布發(fā)散，其與標(biāo)準(zhǔn)偏差表征一致，而L-峭度將成為一個(gè)很小的值，反之亦然。此外，L-峭度值也由λ4確定，它可以衡量兩個(gè)極值(即EX4:4-EX1:4)和兩個(gè)中心值(即EX3：4-EX2：4)之間的距離。如果分布是平坦的，則信號(hào)分布通常是均勻的，并且λ4接近于零。另一方面，如果分布具有尖峰態(tài)，則兩個(gè)極值之間的距離通常較大，而兩個(gè)中心值之間的距離較小，且λ4較大。因此，L-峭度，也可以代表分布特征。

同時(shí)，相關(guān)系數(shù)的表達(dá)式為(9)，其可以表征分解后模態(tài)分量與原始信號(hào)的相關(guān)程度，但對(duì)噪聲敏感。

(9)

因此，考慮兩者的優(yōu)點(diǎn)構(gòu)建L-峭度和相關(guān)系數(shù)的適應(yīng)度函數(shù)如式(10)，該適應(yīng)度函數(shù)具備篩選出與原始信號(hào)相關(guān)性最大的分量和能夠穩(wěn)定檢測信號(hào)中周期性故障脈沖并有效消除噪聲干擾的能力。

fitness=LKCI=(LKI·CI)max

(10)

2.2 WOA優(yōu)化VMD

設(shè)初始最優(yōu)位置為當(dāng)前獵物的位置，其它獵物在最初以該位置目標(biāo)進(jìn)行移動(dòng)，并不斷利用式(11)和(12)更新位置

D=|C·X*(t)-X(t)|

(11)

X(t+1)=X*(t)-A·D

(12)

A=2a·rand-a

(13)

C=2·rand

(14)

式中：rand表示產(chǎn)生的[0,1]之間隨機(jī)數(shù)，隨著循環(huán)迭代的增加，a作為收縮因子將從2線性減小到0，其表達(dá)式為

a=(2-2t/Tmax)

(15)

式中，Tmax表示設(shè)置的最大運(yùn)行迭代次數(shù)。鯨魚通過式(13)和(14)來實(shí)現(xiàn)獵物搜索過程中包圍區(qū)域的不斷縮小。同時(shí)，每個(gè)鯨魚又通過螺旋式來更新其與目標(biāo)位置的距離，其模擬螺旋式更新的數(shù)學(xué)模型可以表示為如下

X(t+1)=D′·ebl·cos(2πl(wèi))+X*(t)

(16)

式中：D′=|X*(t)-X(t)|表示第i只鯨魚和當(dāng)前最優(yōu)位置之間的距離；b為常量系數(shù)，通常為b=1;l為[1,-1]之間的隨機(jī)數(shù)。

為了收縮包圍和螺旋更新同步進(jìn)行，通過概率實(shí)現(xiàn)對(duì)兩種方式的選擇更新，其數(shù)學(xué)表達(dá)式如下

(17)

當(dāng)|A|≥1時(shí)，將隨機(jī)選擇鯨魚遠(yuǎn)離當(dāng)前最優(yōu)鯨魚來尋找最優(yōu)的鯨魚，可以提高算法的全局搜索能力，該數(shù)學(xué)表達(dá)式如下

D=|C·Xrand-X(t)|

(18)

X(t+1)=Xrand-A·D

(19)

式中：Xrand表示隨機(jī)選取鯨魚的位置向量。

通過WOA對(duì)VMD的參數(shù)k和α進(jìn)行組合優(yōu)化，如圖1所示。其步驟如下

圖1 WOA優(yōu)化VMD流程

(1) WOA初始化參數(shù)k和α的取值范圍為k=[2,15]，α=[200,8 000]以避免不完整的信號(hào)分解。

(2) 對(duì)振動(dòng)信號(hào)進(jìn)行VMD分解。

(3) 計(jì)算每個(gè)組合參數(shù)對(duì)應(yīng)位置的適應(yīng)度值。更新適應(yīng)度值大于當(dāng)前值時(shí)的最佳適應(yīng)值。

(4) 確定是否終止迭代。如果t小于Tmax，讓t=t+ 1，并且鯨魚的位置會(huì)更新。否則，迭代終止，并保存最佳參數(shù)。

2.3 改進(jìn)閾值算法

硬閾值和軟閾值是常用的降噪技術(shù)，在各種應(yīng)用中得到了廣泛應(yīng)用[21]，其數(shù)學(xué)表達(dá)式分別為式(20)和(21)

(20)

(21)

(22)

式中：ki是通過試驗(yàn)設(shè)定;σi是第i分量的噪聲標(biāo)準(zhǔn)偏差；N是信號(hào)的長度。噪聲標(biāo)準(zhǔn)偏差通過式(23)獲得

(23)

硬閾值函數(shù)在閾值處是不連續(xù)的，會(huì)產(chǎn)生嚴(yán)重的波動(dòng)。在軟閾值函數(shù)中，不連續(xù)問題得到了解決，但該函數(shù)在閾值位置是不可微的，這將導(dǎo)致對(duì)于隨機(jī)噪聲信號(hào)不能準(zhǔn)確地估計(jì)。在本文中，引用改進(jìn)的閾值降噪[22]公式為

(24)

其中u=α(|ci(t)|-λi),0<αandμ≤1。

在式(24)中，當(dāng)α接近0時(shí)，閾值式(24)中的函數(shù)具有硬閾值函數(shù)功能，當(dāng)α接近1時(shí)，則該閾值函數(shù)的功能接近于軟閾值功能?？梢?，該函數(shù)具有一定的靈活性，參數(shù)α一般設(shè)置為0.05。

3 設(shè)計(jì)的故障特征提取方法步驟

滾動(dòng)軸承運(yùn)行工況復(fù)雜、環(huán)境噪聲污染大導(dǎo)致早期故障信號(hào)中的故障特征難以提取。本文構(gòu)建基于L-峭度和相關(guān)系數(shù)準(zhǔn)則，利用WOA優(yōu)化VMD參數(shù)實(shí)現(xiàn)振動(dòng)信號(hào)的自適應(yīng)分解降噪，再選取模態(tài)分量中噪聲信息少且與原始信號(hào)最大相關(guān)的信號(hào)進(jìn)行改進(jìn)閾值降噪，突出故障特征頻率。其故障特征提取步驟如圖2所示。

圖2 故障提取步驟

(1) 對(duì)采集得到的振動(dòng)信號(hào)進(jìn)行WOA優(yōu)化VMD自適應(yīng)分解；

(2) 計(jì)算分解之后各分量的LKCI指標(biāo)；

(3) 選擇最大LKCI指標(biāo)分量并對(duì)最大指標(biāo)分量進(jìn)行改進(jìn)閾值降噪；

(4) 對(duì)降噪后的信號(hào)進(jìn)行希爾伯特分析，獲得故障特征頻率。

4 試驗(yàn)驗(yàn)證

4.1 試驗(yàn)1

為了驗(yàn)證所提故障特征提取方法的準(zhǔn)確性和有效性，構(gòu)建一組周期性脈沖序列信號(hào)并添加高斯白噪聲，其仿真信號(hào)模型如下

y(t)=s(t)+n(t)

(25)

式中：s(t)為模擬故障沖擊的周期性沖擊分量;位移常數(shù)A0=0.6;阻尼系數(shù)g=0.04;重復(fù)周期T=0.01 s;固有頻率fn=4 000 Hz;采樣點(diǎn)數(shù)N=4 096;采樣頻率為20 kHZ;n(t)為高斯白噪聲。產(chǎn)生的周期性沖擊信號(hào)如圖3所示，對(duì)該沖擊性信號(hào)進(jìn)行加噪獲得強(qiáng)噪聲干擾的混合仿真信號(hào)如圖4所示。

圖3 周期性沖擊信號(hào)

圖4 仿真信號(hào)波形

對(duì)仿真信號(hào)進(jìn)行直接包絡(luò)如圖5所示，可見在頻率500～1 000 Hz，故障特征頻率出現(xiàn)了模糊，其余特征頻率被噪聲淹沒而沒有得到突出。再對(duì)仿真信號(hào)按照本文建立的適應(yīng)度函數(shù)，通過WOA優(yōu)化VMD進(jìn)行分解，然后進(jìn)行Teager能量算子解調(diào)得到結(jié)果如圖6所示。其中Teager能量算子解調(diào)相比直接對(duì)信號(hào)進(jìn)行包絡(luò)的故障特征更加明顯，但在頻率800～1 000 Hz之間出現(xiàn)了故障特征淹沒。繼續(xù)通過建立最小包絡(luò)熵適應(yīng)度函數(shù)，采用WOA優(yōu)化VMD分解得到包絡(luò)熵最小的最佳模態(tài)分量，然后進(jìn)行改進(jìn)閾值降噪，得到包絡(luò)分析結(jié)果如圖7所示，顯而易見其故障特征頻率和其它噪聲頻率相混合而不易觀察，診斷結(jié)果較差。最后，應(yīng)用本文所提出的方法對(duì)仿真故障振動(dòng)信號(hào)進(jìn)行分析如圖8所示，從中可以明顯的看見所有故障特征頻率都得以突出，取得了較好的效果。由上述分析可見，本文所提出方法的可行性和優(yōu)越性。

圖5 仿真信號(hào)包絡(luò)

圖6 Teager能量算子

圖7 包絡(luò)熵優(yōu)化

圖8 本文方法故障特征提取

4.2 試驗(yàn)2

工程試驗(yàn)數(shù)據(jù)選擇來自美國西儲(chǔ)大學(xué)軸承試驗(yàn)中心在6205-2RSJEMSKF型深溝球滾動(dòng)軸承上采集的振動(dòng)加速度信號(hào)，該軸承上的損傷故障采用電火花加工技術(shù)設(shè)置。本文以故障直徑為0.018 mm、負(fù)載為0、轉(zhuǎn)速為1 730 r/min來模擬滾動(dòng)軸承的早期微弱故障。在采樣頻率為12 kHZ下進(jìn)行軸承的加速度信號(hào)采集。由滾動(dòng)軸承的故障特征計(jì)算理論，可以計(jì)算得到內(nèi)圈故障特征頻率fi=156.1 Hz、外圈故障特征頻率fo=103.4 Hz、轉(zhuǎn)頻為Fr=29.95 Hz。首先，對(duì)內(nèi)圈進(jìn)行故障特征提取，其內(nèi)圈的加速度振動(dòng)信號(hào)時(shí)域圖如圖9所示。

圖9 內(nèi)圈振動(dòng)信號(hào)時(shí)域波形

對(duì)內(nèi)圈振動(dòng)信號(hào)應(yīng)用本文建立的適應(yīng)度函數(shù)通過WOA優(yōu)化VMD進(jìn)行分解，其WOA優(yōu)化VMD的迭代曲線如圖10所示。

圖10 WOA優(yōu)化VMD

根據(jù)最優(yōu)組合參數(shù)進(jìn)行VMD分解，得到分解后的IMF分量如圖11所示。

(a) IMF1

其中兩個(gè)IMF分量的LKCI值如圖12所示。

圖12 IMF分量的LKCI值

選擇最優(yōu)的IMF2分量進(jìn)行Teager能量算子解調(diào)如圖13 所示，可見轉(zhuǎn)頻、故障特征頻率的1倍頻和2倍頻較為突出，3倍頻已較弱。為了進(jìn)一步對(duì)比分析本文方法的有效性，以最小包絡(luò)熵為自適應(yīng)度函數(shù)，建立WOA優(yōu)化VMD選取最佳模態(tài)分量，通過改進(jìn)閾值降噪技術(shù)進(jìn)行降噪得到最后的結(jié)果如圖14所示，從圖14中可以看到只有轉(zhuǎn)頻特征較突出，對(duì)于故障特征頻率和2倍頻實(shí)現(xiàn)了較微弱的識(shí)別，其余倍頻都沒有得到有效的識(shí)別。最后，通過本文提出的方法進(jìn)行分析，其中設(shè)定改進(jìn)閾值降噪中的ki=0.3，得到結(jié)果如圖15所示。由圖15可見，內(nèi)圈故障特征頻率以及2倍、3倍和4倍頻得到了顯著的突出，清楚的判斷出屬于內(nèi)圈故障。綜上所分析，本文提出的方法得到了較優(yōu)的特征提取效果。

圖13 Teager能量算子

圖14 包絡(luò)熵優(yōu)化

圖15 本文方法故障特征提取

對(duì)外圈進(jìn)行故障特征提取分析，外圈故障振動(dòng)信號(hào)的時(shí)域波形如圖16所示。同理，應(yīng)用本文建立的適應(yīng)度函數(shù)，采用WOA優(yōu)化VMD進(jìn)行分解，然后進(jìn)行Teager能量算子解調(diào)如圖17所示，可以觀察到1倍頻到4倍頻得到突出，但在高頻部分的故障特診頻率被淹沒而難以辨識(shí)。其次，建立最小包絡(luò)熵適應(yīng)度函數(shù)，采用WOA優(yōu)化VMD后分解后，選取最佳模態(tài)分量進(jìn)行改進(jìn)閾值降噪后的結(jié)果如圖18所示，可以觀察到較微弱的1倍頻和2倍頻故障特征頻率，特征提取效果較差。

圖16 外圈振動(dòng)信號(hào)

圖17 Teager能量算子

圖18 包絡(luò)熵優(yōu)化

最后，通過本文所提的方法進(jìn)行分析，其結(jié)果如圖19所示，故障特征頻率在低頻和高頻部分都明顯突出，提取效果更佳。

圖19 本文方法故障特征提取

5 結(jié)論

針對(duì)復(fù)雜工況和強(qiáng)背景噪聲干擾下，滾動(dòng)軸承早期故障信號(hào)微弱而故障特征頻率難以提取的問題，提出了優(yōu)化VMD與改進(jìn)閾值降噪的滾動(dòng)軸承故障特征提取方法，主要結(jié)論如下：

(1) 該方法對(duì)VMD組合參數(shù)進(jìn)行了優(yōu)化，建立了L-峭度和相關(guān)系數(shù)的最優(yōu)模態(tài)分量選取準(zhǔn)則，結(jié)合改進(jìn)的閾值降噪對(duì)選取的最佳分量進(jìn)行進(jìn)一步的降噪，實(shí)現(xiàn)了故障特征頻率的有效提取，通過仿真和工程數(shù)據(jù)集驗(yàn)證了有效性。

(2) 通過WOA優(yōu)化VMD對(duì)故障信號(hào)進(jìn)行降噪，其建立的L-峭度和相關(guān)系數(shù)準(zhǔn)則能夠最大限度保留故障信息和去除噪聲干擾。并通過與最小包絡(luò)熵的準(zhǔn)則優(yōu)化結(jié)合改進(jìn)閾值降噪方法進(jìn)行了對(duì)比分析，驗(yàn)證了提出準(zhǔn)則的優(yōu)越性。

(3) 通過希爾伯特分析表明，相比于WOA優(yōu)化VMD降噪方法結(jié)合Teager能量算子的故障特征提取方法，本文方法降噪后的故障特征更加突出，具有明顯優(yōu)勢。