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