基于匹配穩態隨機共振的軸承故障診斷方法

池 闊，康建設，張星輝，楊志遠，趙 斐，2

(1.陸軍工程大學石家莊校區，石家莊 050000；2.東北大學工商管理學院，沈陽 110819)

0 引言

軸承廣泛用于旋轉機械，用于支撐轉軸等旋轉部件。其故障常導致設備停機，造成巨大的經濟損失，甚至導致人員傷亡。越早發現軸承故障，就能越早維修設備，避免不必要的損失。軸承故障誘導產生的沖擊信號相對微弱、噪聲背景較強、振動傳播路徑較遠等原因導致很難及時準確地診斷軸承故障狀態。為準確診斷軸承故障，相關研究逐步展開。Laha[1]將圖像處理中常用的局部均值降噪方法拓展到1維信號處理，并通過最大化軸承振動信號的峭度值，確定改進的局部均值降噪方法的各參數。張星輝等[2]采用窄帶干擾消除方法，提取軸承振動信號中的周期沖擊成分，以判斷軸承健康狀態。Bessous等[3]采用離散小波變換分析電機電流變化特征，以診斷電機中的軸承健康狀態。考慮到軸承故障誘導所產生的沖擊十分微弱，需要進一步改進現有研究方法，提高弱故障條件下的軸承故障診斷的準確性。

隨機共振(Stochastic Resonance，SR)是一種利用噪聲增強微弱信號的特殊物理現象。隨機共振認為在特殊條件下噪聲信號不僅無害反而有利，適用于強噪聲背景下的弱信號探測。隨機共振提出后，就開始被應用于各個領域，如信號傳輸[4]、能量收集[5]、圖像處理[6]等。根據絕熱近似理論，傳統的隨機共振只能處理小參數信號(噪聲強度D<<1、驅動信號幅值A<<1、驅動信號角頻率ω0<<克萊莫斯躍遷率)。然而，工程實際中所采集的信號通常為大參數信號，不能滿足小參數限制條件。為解決該問題，大參數隨機共振方法應運而生，如變步長隨機共振[7]、歸一化隨機共振[8]、基于頻域信息交換和變尺度的隨機共振[9]等。這些方法通常將大參數信號轉化為小參數信號或改變非線性系統參數，以滿足小參數限制條件。

自從大參數隨機共振提出后，隨機共振逐步應用于機械故障診斷。陸思良等[10]從經典隨機共振、改進隨機共振和隨機共振在旋轉機械的應用等3個方面，對隨機共振在機械故障診斷領域中的應用研究展開綜述。李繼猛等[11]將時延反饋單穩態隨機共振和自適應最小熵反卷積相結合，用于軸承的故障診斷。池闊等[12]提出基于布谷鳥搜索算法的自適應雙穩態隨機共振方法，用于軸承的故障診斷。時培明等[13]提出基于時延反饋三穩態隨機共振方法，用于旋轉機械的故障診斷。然而，目前基于隨機共振的故障診斷方法主要集中在單穩態、雙穩態或三穩態隨機共振。這些隨機共振方法的勢函數結構單一且所含勢阱數量固定不變，既不適宜從復雜多樣的機械振動信號中提取微弱的故障特征信號，也不利于實現振動信號和勢函數的最佳協同效果。為了豐富勢函數的結構，構造了一種新勢函數，即匹配穩態勢函數。匹配穩態勢函數的結構多樣，所含的勢阱數量可變且由勢函數參數C決定。為使勢函數與機械振動信號達到最佳協同效果，以最大化匹配穩態隨機共振輸出信噪比為目標構造目標函數，并采用布谷鳥搜索算法搜索最佳隨機共振參數和輸出信噪比。相比單穩態、雙穩態和三穩態隨機共振，所提匹配穩態隨機共振具有如下優勢：1)匹配穩態隨機共振勢函數的勢阱數量能夠根據所輸入的機械振動信號進行最優匹配，有利于實現從復雜多樣的機械振動信號中增強微弱的軸承故障特征頻率信號；2)采用布谷鳥搜索算法自適應地改變匹配穩態隨機共振參數，有利于使機械振動信號和勢函數協同效果達到最優。

本文在分析匹配穩態勢函數的基礎上，提出基于匹配穩態隨機共振的軸承故障診斷方法，建立相應的故障診斷框架。通過數值仿真，分析匹配穩態隨機共振各參數對輸出信噪比的影響，對比穩態匹配隨機共振和雙穩態隨機共振的抗噪魯棒性。通過軸承內圈故障案例和滾動體故障案例，對比雙穩態隨機共振和匹配穩態隨機共振的應用效果，驗證所提方法可行性和有效性。

1 匹配穩態隨機共振

隨機共振(Stochastic Resonance，SR)可表述為：在非線性系統中，微弱信號和噪聲信號協同作用下，使得微弱信號顯著增強的現象，可用如下朗之萬方程(Langevin Equation，LE)進行描述：

(1)

其中，x為粒子運動軌跡；S(t)=A0sin(2πfdt)為微弱驅動信號(驅動頻率為fd，幅值為A0)；Γ(t)=(2D)1/2ε(t)為高斯白噪聲(Gaussian White Noise，GWN)，D為噪聲強度，ε(t)為標準GWN(均值為0，方差為1)；U(x)為非線性勢函數。

圖1 雙穩態勢函數形狀示意圖Fig.1 Shape of the bi-stable potential

1.1 雙穩態隨機共振

雙穩態隨機共振(Bi-Stable Stochastic Resonance，BSR)是最經典的隨機共振，其勢函數為雙穩態勢函數UB(x)，可表示為

(2)

其中，a和b為雙穩態勢函數的參數，且a∈R+、b∈R+。UB(x)形狀如圖1所示。顯然，UB(x)具有兩個勢阱(極小值)和一個勢壘(極大值)，且xm=(a/b)1/2、ΔU=a2/(4b)、xb=0。

那么，將式(2)代入式(1)，得雙穩態隨機共振的LE為

(3)

當驅動信號幅值A0較小時，BSR的輸出為[14]

(4)

BSR輸出幅值越大，則隨機共振對驅動信號S(t)的提升效果越好。工程實際中所采集的信號為離散信號。對于離散信號，式(3)可采用Runge-Kutta法求解[12]。令(a,b)=1、A0=0.3，采用5階Runge-Kutta法求解式(3)，得到不同驅動頻率fd下的BSR輸出幅值隨噪聲強度D變化關系，如圖2所示。隨噪聲強度增加，BSR輸出幅值先增加再達到極大值后下降，且極大值遠遠大于驅動信號幅值A0=0.3。因此，當非線性系統固定時，合適的噪聲強度能夠極大地提升BSR輸出幅值。這意味著，在一些特殊條件下，噪聲不僅無害反而有益。同時，隨驅動頻率fd降低，BSR輸出的最大幅值明顯增加，故BSR更適合提升低頻信號。

1.2 匹配穩態隨機共振

1.2.1 匹配穩態隨機共振

匹配穩態隨機共振(Matched-stable Stochastic Resonance，MSR)的勢函數UM(x)為

(5)

其中，c∈(-∞,xm]為穩態決定參數；a和b含義與式(2)相同；O為無窮小量。顯然，匹配穩態勢函數UM(x)由UM1(x)和UM2(x)共同組成。UM2(x)為圓心為(0,B)、半徑為R的過渡小圓弧。圓弧UM2(x)與UM1(x)在|x|=O處相交，且交點處切線相同，即：

(6)

匹配勢函數UM(x)的勢阱數量由參數c決定。若c∈(-∞, -xm]，則UM(x)僅有1個勢阱；若c∈(-xm, 0]，則UM(x)有2個勢阱；若c∈(0,xm]，則UM(x)有3個勢阱。當(a,b)=(1, 1)、c=[-1, 0, 1]時，UM(x)的形狀如圖3所示。

圖2 不同驅動頻率fd下的BSR輸出幅值隨噪聲強度D變化關系Fig.2 Amplitude of BSR output vs. noise intensity D under different driving frequencies fd

圖3 匹配穩態勢函數的形狀Fig.3 Shapes of the matched-stable potentials

將式(5)代入式(1)，得匹配穩態隨機共振的LE：

(7)

其中，

(8)

令z=(b/a)1/2x、τ=at，則式(7)可轉換為

(9)

其中，C=c(b/a)1/2∈(-∞, 1]為穩態決定參數；K=(b/a)1/2為幅值增益。經變換后，參數a和b均轉換為常數1，參數C取值上限變為1，故稱該變換為類歸一化變換。經類歸一化變換后，驅動信號頻率縮小了a倍，幅值擴大了K倍。因此，類歸一化變換能夠將驅動信號轉化為小參數信號，以滿足隨機共振的小參數限制。

1.2.2 數值實施方法和隨機共振輸出評價指標

1) 數值實施方法

工程采集的信號多為離散信號s(n)，且為驅動信號S(n)和噪聲信號Γ(n)的疊加。令s=S+Γ，采用5階Runge-Kutta法求解式(9):

(10)

其中，H=a/fs表示積分步長，fs為采樣頻率。顯然，當輸入信號s固定時，參數(C,H,K)的取值嚴重影響隨機共振的輸出z。

2) 隨機共振輸出評價指標

目前，用隨機共振輸出評價指標有很多。這些評價指標可分為兩類：第一類為需預知驅動頻率的評價指標，如信噪比[13]、局部信噪比[12]、加權信噪比[15]等；第二類為不需預知驅動頻率的評價指標，如近似熵[7]、加權功率譜峭度[16]、MPSK指標[11]、綜合定量指標[17]等。本文采用最常見的信噪比(Signal-to-Noise Ratio，SNR)作為隨機共振輸出評價指標，其定義為

(11)

其中，Ad和An分別為驅動信號和噪聲信號的功率。SNR越大，隨機共振輸出z越好。因此，可建立搜索最優參數(C,H,K)的目標函數

(C,H,K)=arg max[SNR(z)]

(12)

1.3 布谷鳥搜索算法

布谷鳥搜索算法(Cuckoo Search，CS)是Yang[18]提出的一種啟發式群智能優化算法，具有良好的全局尋優能力，已應用于參數估計[19]、圖像降噪[20]、信號處理[12]等研究。本文將采用CS算法搜索式(12)的最優參數和最優值。下面簡要介紹CS。

CS模擬了布谷鳥寄生繁衍策略。它將參數取值范圍視為整個搜索空間，將參數的取值(即解決方案)視為布谷鳥蛋(或鳥巢)，將目標函數值視為布谷鳥蛋的適應度。對于最大化問題，目標函數值(適應度)越大，解決方案(布谷鳥蛋)就越好。為簡單有效地模擬布谷鳥寄生繁衍策略，CS假定了3種理想條件：

1)布谷鳥每次僅產下一個蛋，隨機置于一個鳥巢中；

2)最優質的布谷鳥蛋將保留至下一代；

3)鳥巢(布谷鳥蛋)的數量是固定的，布谷鳥蛋以概率Pa∈(0,1)被寄主鳥發現。一旦布谷鳥蛋被發現，該布谷鳥蛋將被寄主鳥丟棄，新的布谷鳥蛋將隨機生成。

CS包含了3項重要操作，即候選種群產生、擇優選擇和隨機遷移。下面分別介紹這3項重要操作。

1)候選種群產生

(13)

2)擇優選擇

擇優選擇通過比較父體與其候選子代的適應度，保留適應度到下一代。對于最大化問題，擇優選擇可表示為

(14)

其中，f(·)為適應度函數。擇優選擇采用了貪婪策略，是第2條理想假設的實現。該操作不僅有效防止迭代過程中反生退化現象，而且加快適應度收斂速度。

3)隨機遷移

解決方案以概率Pa∈(0,1)被隨機遷移。隨機遷移操作可表示為

(15)

其中，xi,j和yi,j分別為遷移前和遷移后第i個解決方案的第j個參數；xp,j和xq,j分別為為遷移前第p個和第q個解決方案的第j個參數，p和q均為隨機正整數；w和o均為服從(0, 1)均勻分布的隨機數。隨機遷移體現了第3條理想假設，增加了種群多樣性，降低了種群早熟概率。

圖4 基于匹配穩態隨機共振的軸承故障診斷框架Fig.4 Bearing fault diagnosis framework based on MSR

1.4 基于匹配穩態隨機共振的軸承故障診斷框架

建立基于穩態匹配隨機共振的軸承故障診斷框架如圖4所示，具體步驟如下：

1)信號的采集和預處理。采集軸承勻速運轉條件下的振動加速度信號s(n)，用高通濾波器濾除所采信號s(n)中的低頻成分，用希爾伯特變換求取包絡信號se(n)，計算軸承故障特征頻率，采用譜編輯技術消除包絡信號se(n)中低于故障特征頻率且功率譜幅值較高的頻率成分。

2)參數初始化。設置參數(C,H,K)的取值范圍、種群數量n、最大迭代次數Gmax、n個解決方案初始值(記為A)，步長尺度R、常數β、隨機遷移概率Pa等參數，令迭代計數G=1。

3)匹配穩態隨機共振和信噪比計算。根據式(10)求解各初始解決方案A的匹配穩態隨機共振輸出z，然后根據式(11)和軸承故障特征頻率計算各輸出的信噪比。

4)候選種群產生、隨機遷移和邊界約束。根據式(13)和解決方案A產生候選種群，根據式(15)對候選種群進行隨機遷移，得到候選解決方案B。若解決方案B中的某參數超出所設置參數取值范圍，則隨機生成取值范圍內參數值，替代B中該參數。

5)匹配穩態隨機共振和信噪比計算。根據式(10)求解各初始解決方案B的匹配穩態隨機共振輸出z，然后根據式(11)和軸承故障特征頻率計算各輸出的信噪比。

6)擇優選擇。根據式(14)逐一對比解決方案A和B的信噪比，保留信噪比較高的解決方案(記為A)，記錄當前最佳解決方案、最佳隨機共振輸出和最佳信噪比，令G=G+1。

7)判斷是否結束。若G≤Gmax，則跳至步驟4)；否則，輸出最佳解決方案、最佳隨機共振輸出和最佳信噪比。

8)信號分析和軸承故障診斷。分析最優隨機共振輸出的功率譜，判斷軸承健康狀態。

2 參數分析和抗噪魯棒性分析

本節通過仿真試驗，分析MSR參數對輸出的影響和MSR的抗噪魯棒性。在所有仿真試驗中，驅動信號S(t)為正弦信號，噪聲信號Γ(t)為高斯白噪聲，采樣頻率fs為104Hz，信號長度N為3 000。

2.1 參數分析

參數分析對確定MSR各參數對MSR輸出影響具有重要意義。下面通過數值仿真分析各參數對MSR輸出SNR的影響。令正弦信號幅值A0為1，頻率fd為50Hz；高斯白噪聲強度D為1；輸入信號的理論SNR為-6.021dB。為降低輸出SNR隨機性的影響，各仿真重復50次，計算每次MSR輸出SNR，取均值作為最終結果。繪制MSR輸出SNR隨各參數的變化曲線，如圖5所示。

圖5 MSR各參數對輸出SNR的影響分析Fig.5 Influence analysis of every MSR parameter on output SNR

由圖5a可知，(1)當參數H=0.002時，隨著參數C增加，MSR輸出SNR先增加、后降低、再增加、最后降低。此時，在單勢阱和三勢阱狀態下，MSR輸出SNR存在極大值；在雙勢阱狀態下，存在極小值。(2)當參數H=0.06時，隨著參數C增加，MSR輸出SNR先增加、后降低。在雙勢阱狀態下，存在極大值。因此，參數C嚴重影響MSR輸出SNR，但同時受到其他參數的制約。由圖5b(或5c)可知，隨著參數H(或K)的增加，MSR輸出SNR先增加至極大值、后降低。因此，參數H和K也嚴重影響MSR輸出SNR。綜上，想要達到最優MSR輸出SNR，應當同時調節參數(C,H,K)。

2.2 抗噪魯棒性

抗噪魯棒性反映了MSR抵抗不同強度噪聲的能力。令正弦信號幅值A0為1，頻率fd為50Hz；高斯白噪聲強度D從0.4遞增至10(步長為0.4)。采用布谷鳥搜索算法搜索滿足式(12)的最優MSR參數，記錄MSR輸出的最大SNR。每項試驗重復40次，取輸出最大SNR的均值，作為最終結果。同時分析BSR的抗噪魯棒性，作為MSR的對比。繪制抗噪魯棒性分析結果，如圖6所示。

圖6 MSR和BSR的抗噪魯棒性對比Fig.6 Comparison of the anti-noise robustness between MSR and BSR

由圖6a可知，隨噪聲強度D增加，BSR和MSR的輸出SNR逐漸降低，即抗噪能力逐漸減弱。BSR和MSR的輸出SNR均遠遠大于輸入SNR，故SR能夠顯著提升微弱信號。由圖6b可知，MSR輸出SNR總是大于BSR輸出SNR，故MSR抗噪魯棒性優于BSR。值得注意的是，當噪聲強度D大于2后，隨噪聲強度增加，ΔSNR呈現遞增趨勢。也就是說，相比于BSR，MSR更適合探測強噪聲背景下的微弱信號。

圖7 軸承故障模擬試驗臺Fig.7 Bearing fault test rig

3 軸承故障案例驗證

3.1 軸承預植故障試驗

軸承預植故障試驗臺如圖7所示。該試驗臺由三部分組成，即動力與控制部分、軸承故障模擬部分和振動數據采集部分(未在圖7中顯示)。動力與控制部分由電機(提供動力)、電機調頻器(控制電機輸出轉速)和轉速顯示器(顯示電機實時轉速)組成。軸承故障模擬部分由兩深溝球軸承套件(型號ER-12K，主要尺寸如表1所示)、飛輪(提供徑向負載)和光軸組成。振動數據采集部分由振動加速度傳感器、數據采集卡和數據采集軟件組成。

分別進行內圈預植故障試驗(內圈0.5mm寬深溝槽)和滾動體預植故障試驗(單個滾動體0.5mm深溝槽)。預植故障軸承安裝于軸承1位置。傳感器安裝在軸承1附近，采集豎直方向振動加速度數據。每次試驗設置參數如下：電機轉速fr為20 r/s，采樣頻率fs為12.8 kHz，總采樣時間t為1 s。根據電機轉速和軸承主要齒輪，分別計算軸承內圈故障和滾動體故障的特征頻率，如表1所示。

表1 ER-12K軸承主要尺寸和故障特征頻率Tab.1 Main dimensions and fault characteristic frequencies of Bearing ER-12K

3.2 內圈故障案例

首先分析軸承內圈故障信號。采用20階Butterworth高通濾波器(截止頻率為1kHz)濾除低頻信號，采用希爾伯特變換求取包絡信號(消除直流分量)。其波形和功率譜如圖8a和b所示。經過共振解調后，內圈故障特征頻率fic得到增強。由圖8b可知，包絡信號中包含了較強的旋轉頻率諧波2fr和3fr。為降低諧波影響，采用譜編輯技術消除諧波2fr和3fr，即將信號中諧波2fr和3fr對應的傅里葉變換幅值置零再進行傅里葉逆變換，所得信號及其功率譜如圖8c和d所示。將消除諧波后的信號作為驅動信號和噪聲信號的疊加分別輸入到BSR系統和MSR系統中，以最大化SR輸出SNR為目標，用CS搜索最優參數，分別得到最優BSR輸出和最優MSR輸出，如圖8e～h所示。CS參數設置如表2所示。

表2 BSR和MSR系統下的CS參數設置Tab.2 CS parameter set of BSR and MSR

對比圖8d、f和h可知，BSR和MSR的輸出SNR均遠遠大于輸入SNR(-9.237 7 dB)，故BSR和MSR都能夠提升微弱信號。分析圖8f可知，在BSR輸出中，高頻成分得到顯著壓縮，而低頻成分(低于fic)得到了增強，這驗證了BSR的低通濾波的效果。對比圖8f和h可知，1)MSR輸出SNR大于BSR輸出SNR，故MSR性能優于BSR；2)因C=-0.994≈-1，故匹配穩態勢函數可近似認為僅包含一個勢阱，MSR轉換為單穩態隨機共振；3)當MSR轉化為單穩態隨機共振時，MSR輸出的高頻成分得到壓縮，低頻成分略微提高，但小于BSR的低頻提高程度。

圖8 軸承內圈故障信號分析結果Fig.8 Analyzed results of the bearing inner race fault signal

3.3 滾動體故障案例

相比于軸承內圈故障信號，滾動體故障信號就比較微弱了。下面對滾動體故障信號進行分析。采用20階Butterworth高通濾波器(截止頻率為1kHz)濾除低頻信號，采用希爾伯特變換求取包絡信號(消除直流分量)，其波形和功率譜如圖9a和b所示。顯然，滾動體故障頻率fbc十分微弱，幾乎淹沒在強大的噪聲背景中。將包絡信號作為驅動信號和噪聲信號的疊加，分別輸入到BSR系統和MSR系統中，以最大化SR輸出SNR為目標，用CS搜索最優參數，分別得到最優BSR輸出和最優MSR輸出，如圖9c～f所示。CS的參數設置與內圈故障案例相同，如表2所示。

對比圖9b、d和f可知，BSR和MSR的輸出中的頻率fbc成分顯著提高，BSR和MSR的輸出SNR均遠遠大于輸入SNR(-26.487 1dB)，故BSR和MSR都能夠顯著提升微弱信號。對比圖9d和f可知，1)MSR輸出SNR大于BSR輸出SNR，故MSR性能優于BSR；2)因C=0.950∈(0, 1]，故匹配穩態勢函數包含3個勢阱，MSR轉換為三穩態隨機共振；3)當MSR轉化為三穩態隨機共振時，MSR輸出的高頻成分得到壓縮，低頻成分略微提高。

綜合上述兩案例，可得出下列結論：1)MSR能夠在壓縮高頻成分的同時，顯著提升故障特征頻率成分，且效果優于BSR；2)基于MSR的軸承故障診斷方法能夠有效提升軸承故障特征頻率，便于判斷軸承健康狀態。

圖9 軸承滾動體故障信號分析結果Fig.9 Analyzed results of the bearing ball fault signal

4 結論

本文提出了基于匹配穩態隨機共振的軸承故障診斷方法。通過數值仿真分析了匹配穩態隨機共振的參數對輸出的影響及其抗噪魯棒性，得出以下結論：

1)匹配穩態隨機共振的參數C、H和K都能夠顯著影響匹配穩態隨機共振的輸出信噪比，需同時調節3個參數，才能獲得最優匹配穩態隨機共振輸出和最大輸出信噪比；

2)匹配穩態隨機共振的抗噪魯棒性優于雙穩態隨機共振；

3)相比于雙穩態隨機共振，匹配穩態隨機共振更適合探測強噪聲背景下的微弱信號。

采用軸承內圈故障案例和滾動體故障案例對所提故障診斷方法進行驗證，得到以下結論：

1)所提基于匹配穩態隨機共振的軸承故障診斷方法能夠有效提升故障信號中的軸承故障特征頻率，實現軸承的故障診斷，且效果優于基于雙穩態隨機共振的軸承故障診斷方法；

2)與雙穩態隨機共振類似，匹配穩態隨機共振也能夠在壓縮高頻成分的同時，顯著提高故障特征頻率成分。

但本文所提方法僅提升了軸承故障特征頻率而忽略了其諧波。考慮到軸承故障特征頻率諧波能夠在一定程度上表征軸承的健康狀態，下一步將研究基于隨機共振的多頻微弱信號提升方法，以同時提高軸承故障特征頻率及其諧波。