胥永剛,楊苗蕊, 馬朝永
(1.北京工業(yè)大學先進制造技術北京市重點實驗室, 北京 100124;2.北京工業(yè)大學北京市精密測控技術與儀器工程技術研究中心, 北京 100124)
滾動軸承是機械設備中常用的零部件之一,對其進行狀態(tài)監(jiān)測和故障診斷具有重要的理論意義和工程價值[1]?;谡駝有盘柕墓收显\斷方法主要包含特征提取與故障識別2個關鍵環(huán)節(jié),目前常用信號處理作為特征提取的手段。設備振動信號的非線性和非平穩(wěn)性為故障特征提取帶來了困難,傳統(tǒng)信號處理方法往往難以達到理想的效果[2-5]。
基于奇異值分解[6](singular value decomposition, SVD)的信號處理方法是一種有別于傳統(tǒng)分析方法的非線性濾波方法,對于非線性和非平穩(wěn)信號有著較好的效果。利用采集到的振動信號構造軌跡矩陣,然后進行奇異值分解,篩選分解后得到奇異值并重構出新的信號矩陣,可以實現(xiàn)信號的降噪[7],且重構信號沒有相位偏移[8]。結合多分辨奇異值分解包[9]的分解結構和滾動軸承故障信號Hankel 矩陣的奇異值分布特性,黃晨光等[10]提出了延伸奇異值分解包(extended SVD packet, ESVDP)的故障診斷算法。該算法以奇異值分解算法為基礎,實現(xiàn)了振動信號的分解,并以信號能量作為篩選指標,提出了有效分量的篩選準則。該方法具有強魯棒性,極大地改善了奇異值分解包算法的模態(tài)混疊現(xiàn)象。但該方法中的2個關鍵參數(shù)——分解精度和分解層數(shù)是依靠經(jīng)驗值設定的,從而導致了分解結果存在不合理之處。與此同時,分解所得分量中還存在噪聲成分,影響了故障特征的提取效果。
本文以延伸奇異值分解包算法為基礎,利用頻譜趨勢對信號的頻譜進行劃分,以劃分出的頻帶個數(shù)和各頻帶的帶寬來自適應地設定分解精度和分解層數(shù)。同時,結合時域負熵指標,對各頻帶分量進行降噪處理,提高信噪比,最終實現(xiàn)了滾動軸承的故障診斷。
本節(jié)主要對延伸奇異值分解包的理論基礎進行介紹。同時,利用仿真信號來展示原算法的不足。
設有離散信號X=[x1,x2,…,xN],可構造2K行的軌跡矩陣
(1)
式中K為分解精度,K≤2且K為正整數(shù)。
利用奇異值分解算法對軌跡矩陣XA進行分解,即
XA=USVT
(2)
式中:U、V為正交矩陣;S為XA的奇異值對角矩陣。設S中的奇異值為σ1,σ2,…,σ2K,且σ1≥σ2≥…≥σ2K,則S可表示為
S=[diag(σ1,σ2,…,σ2K)]T
(3)
由此,可以得到信號的2K個奇異值。保留排序為2i-1和2i的奇異值,將其余奇異值置零,以此得到新的奇異值矩陣
S′i=[0,diag(σ2i-1,σ2i),0]T
(4)
式中:0表示零矩陣;i∈[1,K]。
再利用新的奇異值矩陣進行重構,得到第1層的第i個分量信號矩陣
(5)
式中:R的上標表示了此分量所在的層數(shù);下標i則表示此分量是該層的第i個分量。
(6)
式中:α=max(1,n-2K+1);β=min(N-2K+1,n);m∈[2,N-1]。
(7)
式中‖·‖L2表示L2范數(shù)。
圖1 延伸奇異值分解包分解結構Fig.1 Decomposition structure of ESVDP
如前所述,在延伸奇異值分解包中存在2個關鍵參數(shù)需要設定——分解精度和分解層數(shù),而在原算法中這2個參數(shù)是依靠經(jīng)驗值進行設定的,從而導致了分解結果不合理。與此同時,分解結果中還存在噪聲成分,影響了后續(xù)的特征提取和故障診斷。
為了直觀地展示延伸奇異值分解包中的不足,構造信號x(t):
x(t)=[4+9sin(2faπt)+8sin(4faπt)]×
sin(2f1πrt)+[5+10sin(2faπt)+8sin(4faπt)]×
sin(2f2π)+[6+9sin(2faπt)+
8sin(4faπt)]×sin(2f3πt)
(8)
式中:fa=15 Hz;f1=300 Hz;f2=500 Hz;f3=700 Hz。
設置采樣點數(shù)為2 048,采樣頻率為2 000 Hz。加入信噪比SNR=1.5的高斯隨機白噪聲后,該信號的波形和頻譜如圖2所示。
圖2 仿真信號x(t)的波形及其傅里葉譜Fig.2 Simulation signal and its Fourier spectrum
利用ESVDP算法處理式(8)所示仿真信號。根據(jù)文獻[8]中提出的關鍵參數(shù)設定范圍,設置分解精度K=3,分解層數(shù)J=5,得到的分解結果如圖3所示。
圖3 延伸奇異值分解包分解結果Fig.3 Decomposition results of the simulation signal by ESVDP
在圖3所示的分解結果中,出現(xiàn)了頻譜混疊的現(xiàn)象。與此同時,以300 Hz和500 Hz為中心的2個頻帶同時出現(xiàn)在了分量1中,且其中的信息變得殘缺不全,這無疑會影響到后續(xù)有效信息的提取。
基于以上幾點不足,本文提出了改進ESVDP算法,利用頻譜趨勢自適應地設置關鍵參數(shù),并結合時域負熵對分解得到的結果進行降噪處理,提高信噪比。
改進ESVDP算法是以延伸奇異值分解包算法為基礎,通過頻譜劃分結果中的頻帶個數(shù)確定分解精度,并利用頻譜劃分結果中的頻帶帶寬確定分解層數(shù),實現(xiàn)信號的自適應分解。同時,通過時域負熵指標提取分解結果中的有效成分,剔除無用的噪聲成分,為后續(xù)的故障診斷奠定基礎。具體流程如圖4所示。
圖4 改進ESVDP算法流程Fig.4 Flow chart of improved ESVDP method
具體步驟為:
1) 參數(shù)設定。預先設置頻譜趨勢估計中的b值和步驟4)中涉及的迭代閾值。
2) 利用關鍵函數(shù)[14-16]對信號進行頻譜趨勢估計。對振動信號x(t)進行快速傅里葉變換(fast Fourier transform, FFT),得到其幅值譜A(f)和相位譜φ(f)。再計算A(f)的傅里葉變換函數(shù),即關鍵函數(shù)
C=FFT(|A(f)|)
(9)
截取關鍵函數(shù)C的前b個點,進行傅里葉逆變換。變換結果近似頻譜的趨勢,稱為趨勢成分Tc(f),即
Tc(f)=iFFT(Cb)
(10)
式中Cb為C的前b個點。
需要指出的是,式(9)中的|A(f)|是一實數(shù)離散序列,對其進行傅里葉變換不會得到x(-t),而是一個新的函數(shù)C。傅里葉逆變換時,b的取值越大,信號趨勢成分中包含的細節(jié)信息越多;反之則細節(jié)信息越少。因此,可以根據(jù)信號的復雜程度選取b值,建議b的取值范圍為[5,15]。
3) 以趨勢成分中的極小值點作為分界點,將頻譜劃分為若干個頻帶。劃分出的頻帶個數(shù)即為分解精度K的值,同時記錄每個頻帶的帶寬和位置作為帶寬范圍ai,i∈[1,K]。
4) 對振動信號進行奇異值分解,同時以能量作為篩選指標篩選最終分解結果。計算每層中第i個有效分量信號的總能量Ei和該信號頻譜中位于帶寬范圍ai之外頻率的信號能量EDi,i∈[1,K]。若EDi與Ei的比值大于1%,則將此分量作為新的待分解信號;反之則將此分量信號保留為最終的分解結果。
5) 為了提高算法效率,在步驟1)中預先設置了迭代閾值。當需要進行下一層分解時,先將當前的分解層數(shù)與迭代閾值進行比較。若當前分解層數(shù)小于迭代閾值,則進行下一層分解;反之則循環(huán)終止,將此分量保留為最終分解結果。分解層數(shù)越多,對計算速度和幅值的影響就越大,經(jīng)多次實驗,筆者認為迭代閾值的設定范圍為[5,15]時的計算效果較好。
(11)
(12)
式中1≤i≤j。
(13)
為了直觀地介紹改進ESVDP算法的實現(xiàn)過程,利用此算法對式(8)構造的仿真信號進行分解。
計算信號的關鍵函數(shù),取其前15個點得到信號趨勢成分,并以趨勢成分中的極小值為分界點劃分頻譜。
將頻譜劃分為了5個頻帶,即分解精度K=5。同時記錄下這5個頻帶的帶寬作為帶寬范圍ai,i=[1,5],如圖5所示。
圖5 頻譜趨勢劃分結果和帶寬范圍Fig.5 Spectrum trend and bandwidth threshold of the simulation signal
利用式(1),將原始信號構造成行數(shù)為10的軌跡矩陣,然后進行奇異值分解,得到第1層的5個有效分量頻譜,如圖6所示。
圖6 第1層有效分量的頻譜Fig.6 Effective components of level 1
由于這5個有效分量的總能量和位于帶寬范圍ai之外頻率的信號能量比值均大于1%,故而將其作為下一層的待分解信號,重復上述分解過程。當分解層數(shù)為7的時候,該層中第1個和第2個有效分量的能量比值小于1%,如圖7所示,將其保留為最終分解結果,其余3個分量繼續(xù)進行分解。
圖7 傅里葉頻譜帶寬對比Fig.7 Comparison of Fourier spectrum bandwidth
圖8中展示了最終的分解結果。這5個分解結果從上到下分別對應了圖5中頻譜趨勢劃分結果的第4個、第2個、第3個、第1個、第5個頻帶。從圖8中可以看出,前3個分解結果包含了有效的信息,后2個分解結果只包含了噪聲干擾。
圖8 分解結果的頻譜Fig.8 Decomposition results of the simulation signal
結合步驟5)中所描述的方法,通過計算分解信號的時域負熵指標最大值,減少了分解結果中的噪聲干擾,結果如圖9所示。
圖9 降噪后分解結果的波形和頻譜Fig.9 Decomposition results after denoising
對比圖8、9可以看出,各個分解結果中的無效成分已經(jīng)基本被去除掉了,突出了有效信息,有利于后續(xù)的特征提取和故障診斷。
為了進一步驗證所提出方法的有效性,采集了滾動軸承故障模擬實驗臺上的振動信號,實驗臺如圖10所示。
圖10 滾動軸承故障實驗臺Fig.10 Rolling bearing test rig
利用該實驗臺采集了6307滾動軸承外圈故障振動信號,電機轉速為1 496 r/min,采樣頻率為15 360 Hz,采樣點數(shù)為8 192。經(jīng)計算,該軸承的外圈故障特征頻率為fo=76.88 Hz。
圖11為采集到的振動信號和根據(jù)其趨勢成分進行頻帶劃分的結果,其中黑色虛線為所劃分的邊界,紅線為頻譜趨勢。如圖所示,信號波形中的周期性沖擊并不明顯,存在較多噪聲成分。
圖11 實驗信號及其頻譜劃分Fig.11 Experimental signal and its spectrum trend
由于噪聲較大,其中的有效信息被湮沒,給后續(xù)的故障診斷帶來了困難。利用提出的改進ESVDP算法對此振動信號進行處理。根據(jù)圖11所示的頻譜劃分結果,可知分解精度K=5,同時獲得了5個帶寬范圍ai,i∈[1,5]?;诖?得到的分解結果如圖12所示。
圖12 實驗信號分解結果Fig.12 Decomposition results of the experimental signal
從圖12所展示的分解圖中不難看出,前2個分量中出現(xiàn)了周期性的沖擊成分,同時其頻譜中包含了明顯的邊頻帶。而第3~5個分量中的有效信息較少,因此在后續(xù)的分析中,只展示了前2個分量的處理結果。
結合第2節(jié)中介紹的算法步驟5),計算分解結果的時域負熵指標,以此來選取需要保留的特征信息,從而實現(xiàn)分解結果的降噪處理。
圖13展示了前2個分量的時域負熵指標,其中紅色標星處為最大值。
圖13 分量1和分量2的時域負熵Fig.13 Time domain negative entropy of decomposition result 1 and 2
降噪處理后的分量1和分量2的波形、頻譜,以及包絡譜(0~800 Hz)如圖14所示。
圖14 降噪處理后的實驗信號分解結果Fig.14 Decomposition results of the experimental signal after denoising
降噪處理后的分量1和分量2的包絡譜中,低頻部分均出現(xiàn)了明顯的峰值,對應著6307軸承外圈故障特征頻率及其倍頻。因此,可以診斷出軸承外圈發(fā)生故障。
1) 通過關鍵函數(shù)的重構得到振動信號的頻譜趨勢成分,實現(xiàn)了頻帶的劃分,從而自適應地設定了分解精度和帶寬范圍,并通過帶范圍確定了分解層數(shù)。
2) 結合時域負熵指標對分解結果進行了降噪處理,有效地去除了干擾成分,提高了信噪比,為后續(xù)的故障診斷打下基礎。
3) 結合仿真和實驗對該方法進行了驗證。結果表明該方法可以從復雜的振動信號中準確地提取出故障特征,最終實現(xiàn)了滾動軸承的故障診斷。