俞森， 馬潔

(北京信息科技大學(xué)機(jī)電工程學(xué)院， 北京 100192)

0 前言

軸承被譽(yù)為“工業(yè)的關(guān)節(jié)”，在汽車、高鐵、風(fēng)力發(fā)電、冶金、化工、船舶以及大型制造業(yè)中廣泛使用。但由于其在高壓、碰撞、高溫等復(fù)雜的工作環(huán)境極易發(fā)生損傷，在生產(chǎn)生活中，約1/3的機(jī)械設(shè)備損壞由軸承故障引起[1]。因此，如何更早地發(fā)現(xiàn)軸承故障并進(jìn)行維修保養(yǎng)，對減少事故發(fā)生和降低經(jīng)濟(jì)損失非常重要。軸承故障診斷主要涉及2個任務(wù)：特征提取和故障分類[2-3]。

基于振動信號分析的故障診斷方法能有效診斷機(jī)械設(shè)備缺陷和破損[4]。受工作環(huán)境影響，振動信號特征往往被噪聲覆蓋[5]，因此，噪聲處理是滾動軸承故障診斷中的關(guān)鍵問題之一。1980年，MARR 和HILDRETH[6]首次提出拉普拉斯高斯濾波器，被廣泛應(yīng)用于圖像邊緣檢測問題中，但在振動信號中并沒有得到好的應(yīng)用；SAAD等[7]將拉普拉斯濾波器用于地震波去噪，計算地震波到達(dá)時間；LIU等[8]將teager能量算子和拉普拉斯高斯濾波器結(jié)合，用于變轉(zhuǎn)速工況下的特診提取，但濾波器參數(shù)需人為設(shè)定； SAM(Spectral Amplitude Modulation)由MOSHREFZADE等于2019年提出[9]，通過給原始信號賦予不同的權(quán)值重構(gòu)信號，然后計算重構(gòu)信號的平方包絡(luò)譜尋找故障特征頻率，但受噪聲影響較大。因此，本文作者利用麻雀算法自適應(yīng)地選擇MLOG的最優(yōu)參數(shù)，并結(jié)合SAM方法對信號重構(gòu)，從而提取故障特征頻率。

1 理論分析

1.1 拉普拉斯高斯濾波器

高斯濾波可用于消除高斯噪聲，對圖像噪聲處理效果優(yōu)異[10]，得到了廣泛應(yīng)用。由于圖像由不同灰度級的像素組成，不同灰度級即可看作圖像特征信息的突變，在軸承振動信號中，噪聲信號和振動信號頻率、相位、幅值等存在差異，采集同一時刻的故障信號可視為特征信息突變，因此，文中提出利用拉普拉斯高斯濾波器對振動信號去噪。

一階高斯濾波函數(shù)定義式為

(1)

式中：σ為標(biāo)準(zhǔn)差；n為高斯指數(shù)。對式(1)求二階導(dǎo)并將結(jié)果歸一化得：

(2)

對于單位沖擊響應(yīng)濾波器，當(dāng)系數(shù)之和為零時，噪聲處理效果最優(yōu)，此時可很好地提取信號的瞬態(tài)脈沖。式(3)為改進(jìn)的拉普拉斯高斯濾波器數(shù)學(xué)表達(dá)式：

(3)

式中：N為濾波器的階數(shù)。

1.2 麻雀搜索算法

麻雀搜索算法由XUE和SHEN[11]于2020年提出，是一種新的優(yōu)化算法。通過模仿麻雀捕食行為，將種群分為發(fā)現(xiàn)者、加入者和預(yù)警者。發(fā)現(xiàn)者能量儲備較高，負(fù)責(zé)搜尋食物；加入者跟隨發(fā)現(xiàn)者捕食增加能量儲備，當(dāng)發(fā)現(xiàn)有危險靠近時，預(yù)警者預(yù)警，進(jìn)行反哺食行為，此時發(fā)現(xiàn)者將所有麻雀引到安全區(qū)域。生產(chǎn)者為麻雀種群中隨機(jī)選擇能量儲備最高的個體，但數(shù)量保持不變。

設(shè)n個麻雀搜索d維空間，麻雀位置表示為

(4)

麻雀的適應(yīng)度值為

(5)

發(fā)現(xiàn)者位置更新公式為

(6)

式中：t為當(dāng)前迭代次數(shù)；tmax為最大迭代次數(shù);α∈(0,1)；Q服從正態(tài)分布；I是元素都為1的1×d矩陣；R2∈(0,1)，ST∈(0.5,1)；R2

加入者位置更新公式為

(7)

式中：Xp和Xworst分別為當(dāng)前發(fā)現(xiàn)的最優(yōu)位置和全局中當(dāng)前最差位置；A為1×d的矩陣，元素為1或-1隨機(jī)分布，且A+=AT(AAT)-1。

隨機(jī)選擇麻雀總數(shù)的10%～20%作為預(yù)警者，位置更新公式為

(8)

1.3 頻譜調(diào)幅(SAM)

在倒譜預(yù)白化過程中，由于單位幅度的任意設(shè)定而不能更好地表現(xiàn)頻譜的特征，文獻(xiàn)[12]中提出了倒譜預(yù)白化的新方法：

Xm(t,n)=IFT[A(f)nejφ(f)]

(9)

式中：Xm為重構(gòu)信號；n為數(shù)量級。

均值為零的周期性信號x(t)的傅里葉級數(shù)表達(dá)式如下：

(10)

式中：D為傅里葉系數(shù)；q為信號x(t)的長度。

對x(t)做Hilbert變換，有：

(11)

計算原始信號的平方包絡(luò)或范數(shù)，有：

cos[(ωi-ωj)t+φi-φj]

(12)

重構(gòu)信號xm由x(t)的幅值和相位確定，則有：

(13)

計算xm的平方包絡(luò)譜為

SES(fk,n)=(2-δ(k))·

(14)

2 麻雀算法優(yōu)化MLOG數(shù)學(xué)模型

峭度與軸承自身屬性無關(guān)，對沖擊信號極其敏感，正常軸承峭度值為3。研究表明：峭度值會隨軸承的損壞程度而逐漸增大，特別適用于早期故障的診斷問題。丁金玲等[13]以峭度最大值為標(biāo)準(zhǔn)選擇小波包分解后的最優(yōu)頻帶系數(shù)，準(zhǔn)確提取特征頻率。馬晨佩等[14]利用麻雀算法搜索SVM的最優(yōu)參數(shù)，有效解決了軸承的故障診斷問題。因此文中提出基于麻雀算法優(yōu)化拉普拉斯濾波器的方法，以濾波后信號的峭度值為適應(yīng)度函數(shù)，選擇濾波器的最優(yōu)參數(shù)，根據(jù)峭度最大值選擇合適的濾波器參數(shù)，式(15)為峭度計算公式：

(15)

式中：x為振動信號；μ為均值；σ為標(biāo)準(zhǔn)差。

參數(shù)優(yōu)化流程如圖1所示，具體步驟如下：

(1) 輸入故障信號，設(shè)置濾波器尺度參數(shù)和階數(shù)N的上下界，初始化SSA算法的最大迭代次數(shù)和種群數(shù)量；

(2) 濾波器對故障信號濾波，根據(jù)式(15)計算濾波后信號的峭度并保存；

(3)判斷是否達(dá)到終止條件，未達(dá)到則更新位置繼續(xù)迭代。選擇濾波后信號最大峭度值對應(yīng)的N為最優(yōu)參數(shù)；

(4) 用最優(yōu)參數(shù)對原始信號濾波；

(5) 濾波后的信號用SAM處理，分析頻譜特征，識別故障特征頻率。

圖1 故障診斷流程 Fig.1 Flow of fault diagnosis

3 仿真信號分析

為了驗(yàn)證文中方法的有效性，構(gòu)造軸承故障時仿真信號如下：

y(t)=x(t)+n(t)

(16)

(17)

式中：x(t)為周期性沖擊信號；n(t)為隨機(jī)噪聲信號；A為位移常數(shù)；ζ為阻尼系數(shù)；fn為軸承固有頻率。

設(shè)故障沖擊周期為T0=0.01 s，則故障特征頻率為ft=100 Hz，采樣頻率f=20 kHz，采樣點(diǎn)數(shù)N=10 000，令A(yù)=2，ζ=0.2，fn=2 000 Hz，n(t)是人為加入的高斯白噪聲，此處信噪比SNR=-8 dB，圖2所示為y(t)的時域和包絡(luò)譜，從中無法識別任何故障成分。

圖2 仿真信號

設(shè)置最大迭代次數(shù)20，種群數(shù)量10，濾波器階數(shù)N∈[1,10]，標(biāo)準(zhǔn)差σ∈[0.1,5]，收斂曲線如圖3所示，得到濾波器最優(yōu)參數(shù)(N=10，σ=3)。

圖3 仿真信號最大峭度收斂曲線

用選擇的最優(yōu)參數(shù)對信號用拉普拉斯濾波器濾波，將濾波后信號用SAM處理，設(shè)置1≤MO≤2，得SAM的三維圖如圖4所示。選擇最高峰值對應(yīng)的權(quán)重值MO=1.5對信號重構(gòu)，重構(gòu)信號的時域和對數(shù)平方包絡(luò)譜如圖5所示，可明顯看到一倍頻(100 Hz)及其他倍頻成分，可判定軸承故障。

圖4 仿真信號的SAMFig.4 SAM of simulation signal

圖5 重構(gòu)信號

4 實(shí)測數(shù)據(jù)分析

4.1 內(nèi)圈故障數(shù)據(jù)分析

為驗(yàn)證文中方法的實(shí)用性，選取美國凱斯西儲大學(xué)內(nèi)圈故障數(shù)據(jù)進(jìn)行驗(yàn)證，采用6203-2RS JEM SKF深溝球軸承，故障直徑0.177 8 mm，軸承轉(zhuǎn)速n=1 797 r/min，軸承參數(shù)見表1，根據(jù)公式計算得故障特征頻率為fI=147.6 Hz。

表1 6203-2RS JEM SKF軸承參數(shù)Tab.1 Bearing parameters of 6203-2RS JEM SKF

采樣頻率fn=12 000 Hz，采樣點(diǎn)數(shù)N=10 000，振動信號時域和包絡(luò)譜如圖6所示，可看出:轉(zhuǎn)頻fr=37.5 Hz，特征頻率fi=147.9 Hz，與理論計算值基本一致，但其他倍頻成分被噪聲掩蓋。

圖6 內(nèi)圈故障信號

設(shè)置麻雀算法最大迭代次數(shù)20，種群數(shù)量10，濾波器階數(shù)N∈[1,10]，標(biāo)準(zhǔn)差σ∈[0.1,5]，收斂曲線如圖7所示，迭代第2次得最大峭度值K=4.99，得最優(yōu)參數(shù)(N=8，σ=2.24)。

圖7 內(nèi)圈故障信號最大峭度收斂曲線

用選擇的最優(yōu)參數(shù)對內(nèi)圈信號用拉普拉斯濾波器濾波，將濾波后信號用SAM處理，設(shè)置-1≤MO≤0.5，得SAM的三維圖如圖8所示，選擇最高峰值對應(yīng)的權(quán)重值MO=-0.3對信號重構(gòu)，重構(gòu)信號的時域和對數(shù)平方包絡(luò)譜如圖9所示，可明顯看出：一倍頻及其他倍頻成分與計算值接近，說明軸承內(nèi)圈故障。

圖8 內(nèi)圈故障信號SAMFig.8 SAM of inner race fault signal

圖9 重構(gòu)信號

4.2 外圈故障數(shù)據(jù)分析

為再次驗(yàn)證此方法的有效性，選擇西交大全生命周期故障數(shù)據(jù)進(jìn)行驗(yàn)證[15]，軸承型號為LDK UER204，滾動軸承參數(shù)見表2，選擇Bearing 2_5外圈故障數(shù)據(jù)集第120個樣本，采樣頻率fn=25 600 Hz，轉(zhuǎn)速n=2 250 r/min，根據(jù)公式計算得外圈故障特征頻率為fo=115.6 Hz。

表2 LDK UER204軸承參數(shù)Tab.2 Bearing parameters of LDK UER204

外圈故障信號時域如圖10(a)所示，用SSA-MLOG-SAM對信號處理，適應(yīng)度函數(shù)的收斂曲線如圖10(b)所示，得到濾波器最優(yōu)參數(shù)(N=8，σ=1.87)。濾波后信號經(jīng)SAM處理重構(gòu)后的平方包絡(luò)譜如圖10(c)所示，可看到故障特征頻率及其他倍頻成分，證明軸承發(fā)生外圈故障。說明文中方法可準(zhǔn)確識別不同試驗(yàn)臺的不同故障類型。

圖10 外圈故障信號分析

為表明麻雀算法優(yōu)化濾波器參數(shù)的優(yōu)越性，選擇粒子群算法(Particle Swarm Optimization，PSO)和遺傳算法GA對比驗(yàn)證，優(yōu)化迭代曲線如圖11所示。對比結(jié)果顯示：麻雀算法在最少的迭代次數(shù)內(nèi)達(dá)到最大適應(yīng)度值，優(yōu)于PSO和GA。

圖11 不同算法優(yōu)化濾波器參數(shù)迭代曲線

為表明拉普拉斯濾波效果的優(yōu)越性，選擇MCKD-SAM對內(nèi)圈故障數(shù)據(jù)分析，設(shè)置濾波器長度L=160,移位數(shù)M=1，包絡(luò)譜如圖12所示，只能識別一倍頻，其他倍頻被淹沒在噪聲中。不對故障信號預(yù)處理，直接SAM后的包絡(luò)譜如圖13所示：從中只能識別一倍頻和二倍頻。對比結(jié)果顯示：文中提出的方法對軸承早期故障特征提取更有效。

圖12 MCKD-SAM包絡(luò)譜Fig.12 Envelope spectrum of MCKD-SAM

圖13 內(nèi)圈故障信號SAM包絡(luò)譜Fig.13 SAM envelope spectrum of inner race fault signal

5 結(jié)論

文中提出了SSA-MLOG-SAM的滾動軸承早期故障特征提取方法，通過仿真信號的不同試驗(yàn)臺數(shù)據(jù)分析表明：

(1) 利用SSA算法自適應(yīng)選擇拉普拉斯高斯濾波器的最優(yōu)參數(shù)，克服了人為設(shè)定參數(shù)的缺點(diǎn)，且優(yōu)化效率高于傳統(tǒng)粒子群算法和遺傳算法；

(2) 驗(yàn)證了用于圖像去噪的拉普拉斯高斯濾波器也能對振動信號去噪，且效果優(yōu)于傳統(tǒng)濾波方法；

(3)SAM利用非線性濾波分離信號成分，對濾波后信號賦予特定的權(quán)重來重構(gòu)信號，通過平方包絡(luò)譜分析準(zhǔn)確地提取受強(qiáng)噪聲干擾的故障特征頻率，為滾動軸承的故障診斷提供了新的解決方法。