一階線性系統(tǒng)隨機共振在轉(zhuǎn)子軸故障診斷中的應(yīng)用研究

冷永剛，田祥友

（1．天津大學(xué) 機械工程學(xué)院，天津 300072；2．天津大學(xué) 機構(gòu)理論與裝備設(shè)計教育部重點實驗室，天津 300072）

大型汽輪機、水輪機、燃氣輪機等旋轉(zhuǎn)動力機械常工作在高溫、高壓和高轉(zhuǎn)速等較惡劣的工作環(huán)境下，其正常運行對于保證工廠的安全生產(chǎn)意義重大。因此需要對旋轉(zhuǎn)機械進行狀態(tài)監(jiān)測，提取振動信號并進行分析以判別是否存在早期微弱故障。但在工程實際中，信號測點往往離故障源較遠，使測得的信號容易被強噪聲污染，加上早期故障信號本身就比較微弱，這些因素使轉(zhuǎn)子軸的早期微弱故障診斷較為困難。

目前常用的轉(zhuǎn)子軸類故障診斷方法有：基于諧波窗、EEMD濾波的轉(zhuǎn)子軸心軌跡提純方法［1－2］，利用小波或小波包對故障信號進行分解和重構(gòu)［3］方法，以及基于非線性系統(tǒng)隨機共振檢測微弱信號的振動分析等［4－6］。這些檢測方法各有優(yōu)點，在一定程度上能反映出轉(zhuǎn)子軸的故障特征，但同時也存在一些局限性。如轉(zhuǎn)子軸心軌跡提純法步驟較繁瑣，測試數(shù)據(jù)時對硬件設(shè)備要求較高；小波或小波包降噪時如何選擇閾值和進行閾值優(yōu)化比較困難；利用非線性系統(tǒng)的隨機共振原理降噪時存在系統(tǒng)參數(shù)調(diào)節(jié)等復(fù)雜問題。本文以一階線性系統(tǒng)隨機共振方法為基礎(chǔ)，通過模型簡化，提出一種基于線性系統(tǒng)調(diào)參廣義隨機共振的轉(zhuǎn)子軸類早期微弱故障診斷方法。

1 一階線性系統(tǒng)的調(diào)參廣義隨機共振現(xiàn)象

很多學(xué)者［7－10］關(guān)于線性系統(tǒng)隨機共振現(xiàn)象的研究都是通過復(fù)雜的數(shù)學(xué)模型推導(dǎo)與仿真得出系統(tǒng)輸出存在傳統(tǒng)的或是廣義的隨機共振現(xiàn)象，但這樣的系統(tǒng)中信號與噪聲的作用方式比較復(fù)雜，很難應(yīng)用于實際故障信號的處理分析中。本文僅針對周期信號加高斯白噪聲作用的一階線性系統(tǒng)中出現(xiàn)的廣義隨機共振現(xiàn)象［11］能否用于轉(zhuǎn)子軸類早期微弱故障信號的檢測進行討論。

一階線性系統(tǒng)動力學(xué)模型

式中：a是系統(tǒng)參數(shù)，輸入信號 sn（t）＝A sin（2πf0t）＋ξ（t），A sin（2πf0t）是頻率為 f0、幅值為 A、相位為 0的正弦信號，ξ（t）是功率譜密度為2D的噪聲，其中ξ（t）是均值為0，方差為1的高斯白噪聲。

求解系統(tǒng)模型（1），得系統(tǒng)輸出功率譜密度函數(shù)為

為比較系統(tǒng)輸出端信號與噪聲的大小關(guān)系，引入輸出信噪比。這里信噪比的定義采用信號處理及通信等工程實際中常用的信噪比定義，即系統(tǒng)輸出端信號總功率與噪聲總功率之比［12－13］

式中：S（f）為輸出端信號的功率譜密度函數(shù)，SN（f）為輸出端噪聲功率譜密度函數(shù)。

由輸出信噪比表達式（4）可以看出，系統(tǒng)參數(shù)a不變時，SNRout隨噪聲強度D單調(diào)遞減變化，表明周期信號與加性噪聲驅(qū)動的一階線性系統(tǒng)中不會出現(xiàn)傳統(tǒng)意義上的隨機共振現(xiàn)象。但仔細觀察式（4）可知，分母中時等號成立），故在噪聲強度D不變的情況下，系統(tǒng)參數(shù)a取2πf0時，系統(tǒng)輸出信噪比取得最大值SNRmax＝當(dāng)輸入信號的頻率與系統(tǒng)參數(shù)取得某種協(xié)同匹配關(guān)系時，系統(tǒng)輸出能夠達到共振狀態(tài)。這種輸出信噪比隨系統(tǒng)參數(shù)非單調(diào)變化的現(xiàn)象稱之為調(diào)參廣義隨機共振［11］。

下面通過數(shù)值模擬驗證系統(tǒng)的上述特性，設(shè)定一組仿真參數(shù)：A＝0．1，f0＝0．01 Hz，采樣頻率 fs＝5 Hz，計算數(shù)據(jù)點數(shù)N＝4 000。用四階Runge-Kutta方法（以下同）對方程（1）進行數(shù)值計算，觀察輸出信噪比分別隨參數(shù)D、a的變化趨勢，計算結(jié)果如圖1。算法中的輸出信噪比定義為

式（5）是給定計算數(shù)據(jù)點數(shù)條件下的信噪比，其中PS表示輸出功率譜中信號總功率，PN表示輸出功率譜中噪聲總功率。

圖1 輸出與輸入信噪比分別隨參數(shù)D、a的變化趨勢Fig．1 The variation tendency of SNRout and SNRin by the parameter D＆a．Parameter

同樣由圖1仿真曲線可以看出，線性系統(tǒng)的輸出信噪比隨輸入噪聲強度D單調(diào)遞減變化，隨系統(tǒng)本身參數(shù)a非單調(diào)變化，存在a值（a＝0．064，約等于2πf0）使輸出信噪比達到極值。圖中虛線表示相應(yīng)情況下的輸入信噪比，其定義為在給定數(shù)值分析長度條件下，信號總功率與噪聲總功率之比。

2 一階線性系統(tǒng)調(diào)參廣義隨機共振的微弱信號檢測

由前述分析知，在給定輸入信號與噪聲時，系統(tǒng)取參數(shù)a＝2πf0輸出信噪比達到極大值。為了滿足實際工程應(yīng)用，文獻［11］分析了系統(tǒng)參數(shù)a＝2πf0固定不變情況下，采樣頻率與信號頻率比值對輸出頻譜的影響，并給出了檢測特征信號可辨識性好的比值范圍。這表明，實際應(yīng)用中，使系統(tǒng)輸出信噪比最大的參數(shù)a值并不一定就能使輸出特征信號達到理想的可辨識性，因為信噪比大，只說明系統(tǒng)輸出的信號能量與噪聲能量比值大，而特征信號的可辨識性與信號幅值、噪聲大小以及噪聲的分布有關(guān)。因此需要進一步分析，實際應(yīng)用時應(yīng)如何調(diào)節(jié)a值才能使特征信號的辨別性最優(yōu)。下面仿真說明高斯白噪聲條件下，輸出頻譜上特征信號可辨識性與參數(shù)a的變化關(guān)系。

根據(jù)文獻［11］比較輸出頻譜圖上特征信號可辨識程度的方法，即提取輸出譜圖上信號頻率處的峰值h，將其與整個譜圖上除去信號頻率處譜峰后其余譜值中最大的一個值hr作比較，比值越大表明特征信號的可辨識性越好，并將這個比值定義為辨別率r＝h／hr．，一般只有r＞1時特征信號才容易識別出來。在采樣頻率與信號頻率比值不變的情況下，改變參數(shù)a值進行仿真，計算不同a值下輸出譜圖上特征信號的辨別率r。仿真參數(shù)同上：A＝0．1，f0＝0．01 Hz，D＝0．8，采樣頻率fs＝5 Hz，計算數(shù)據(jù)點數(shù) N＝4 000，令參數(shù) a與2πf0的比值在0．1～10之間變化，結(jié)果如圖2所示。

由圖2可以看出，保持采樣頻率與信號頻率比值fs／f0＝500不變的情況下，輸出譜圖特征信號的辨別率r與參數(shù)a的取值有很大關(guān)系，隨a／2πf0比值的變化趨勢是先增大至某一極大值，并幾乎保持一段水平值，然后再減小，當(dāng)a／2πf0＝4～9時存在一個水平極大值。由此可見選取一個合適的參數(shù)a值對特征信號頻率的判別有著直接的影響，實際信號檢測時應(yīng)盡量選取使辨別率達到最佳的a值。圖3比較了參數(shù)a取兩個不同值時輸出譜圖中特征信號的可辨識性，其余參數(shù)同圖2。

圖2 系統(tǒng)輸出譜圖上特征信號辨別率r隨 a／2πf0變化曲線Fig．2 The characteristic signal recognition r changes with a／2πf0 in system output spectrum

由圖3可見，圖（b）比圖（a）的特征信號辨別性明顯要好。雖然圖（a）中f0＝0．01 Hz處的峰值已經(jīng)是整個譜圖上的最大值，但與其附近峰值高度相接近，容易造成誤判，而圖（b）在f0＝0．01 Hz處的譜值明顯高于附近其他峰值，很容易確定這就是所尋找的特征信號頻率成分。

圖3 不同參數(shù)a值下系統(tǒng)輸出頻譜Fig．3 The output spectrum of different parameter a

這里需要指出的是，圖2得到的辨別率r隨比值a／2πf0變化關(guān)系是在采樣頻率fs與信號頻率f0比值保持恒定的情況下得到的，如果改變采樣頻率與信號頻率的比值，那么r隨比值a／2πf0的變化關(guān)系曲線也會發(fā)生變化。為綜合分析r隨比值a／2πf0和fs／f0的變化關(guān)系，圖4給出了三者變化關(guān)系的三維圖。從圖4可見，對應(yīng)于每一組較優(yōu)的fs／f0值，均可以找到合適的a／2πf0值，使特征信號辨別率r達到最大。另外，應(yīng)用此方法處理含噪微弱信號時，可不必局限在小頻率范圍，只要采樣頻率合適，可以用于檢測任意頻率大小信號［11］。

圖4 特征信號辨別率r隨比值a／2πf0和 fs／f0變化三維圖Fig．4 The 3D graph of characteristic signal recognition r changes with a／2πf0 and fs／f0

根據(jù)以上分析，應(yīng)用線性系統(tǒng)調(diào)參廣義隨機共振檢測實際微弱信號時的步驟可總結(jié)為：① 根據(jù)已知條件估算故障信號特征頻率f0的可能范圍，參考特征信號辨別率隨比值 a／2πf0和 fs／f0變化的三維圖，確定合適的采樣頻率fs和參數(shù)a使信號可識別性最佳。②若輸出譜圖上信號特征頻率不明顯，可通過微調(diào)參數(shù)a的值進行修正或重新確定一組采樣頻率fs和參數(shù)a采樣，再次識別計算。

3 轉(zhuǎn)子軸早期微弱故障識別分析

在旋轉(zhuǎn)機械中軸彎曲是一種常見故障，轉(zhuǎn)子軸彎曲是指轉(zhuǎn)軸各截面幾何中心連線與旋轉(zhuǎn)軸線不重合，從而使轉(zhuǎn)子產(chǎn)生偏心質(zhì)量，進而引起不平衡振動。實驗在滑動軸承故障試驗臺上進行，轉(zhuǎn)子軸中心偏離旋轉(zhuǎn)軸線0．38 mm，在試驗臺遠離軸承座向外0．5 m處布置一個加速度傳感器，一方面模擬工程實際中由于空間限制無法在軸承座附近布置傳感器采集信號的情況，另一方面使故障信號強度進一步衰減，模擬早期微弱故障。

為獲得這一早期微弱故障特征信號，利用一階線性系統(tǒng)調(diào)參共振方法進行處理時，首先根據(jù)軸的轉(zhuǎn)頻和故障類型，判定彎曲故障信號頻率f0值，實驗中轉(zhuǎn)子軸的轉(zhuǎn)速約為1 700 r／min，則故障信號頻率 f0約在28 Hz左右。利用NI數(shù)據(jù)采集儀和Lab VIEW采集軟件進行數(shù)據(jù)采集，依據(jù)特征信號可辨識性最好的采樣頻率取值原則［11］，設(shè)置采樣頻率為 fs＝1 kHz，采樣點數(shù)N＝5 000，實測信號的時域波形及低頻頻譜如圖5所示?？芍?，信號譜干擾成分較多，很難通過直接觀察的方法從頻譜圖上準(zhǔn)確判斷出故障特征。參考圖4取線性系統(tǒng)參數(shù) a＝1．1×2πf0，將實測信號輸入式（1）系統(tǒng)進行處理，得到圖6時頻域結(jié)果。

圖5 實測轉(zhuǎn)子軸彎曲故障信號的時域圖及頻域圖Fig．5 The time domain and frequency domain of observed rotor shaft bending fault signal

圖6 轉(zhuǎn)子軸彎曲故障信號經(jīng)過線性系統(tǒng)處理后的時域及頻域圖Fig．6 The time domain and frequency domain of rotor shaft bending fault signal after being processed in linear system

圖7 轉(zhuǎn)子軸故障信號經(jīng)過非線性雙穩(wěn)系統(tǒng)處理后的時域及頻域圖（參數(shù) a＝1，b＝10）Fig．7 The time domain and frequency domain of rotor shaft fault signal after being processed in non-linear bistable system．（parameter a＝1，b＝10）

由圖6可以看到，經(jīng)過處理后輸出頻譜圖上27．8 Hz處有一明顯譜峰，且遠高于周圍其他譜峰，由此識別出滑動軸承轉(zhuǎn)子存在“早期微弱”彎曲不平衡故障。

作為比較，以下采用非線性雙穩(wěn)隨機共振方法處理該微弱故障信號［5－6］，非線性雙穩(wěn)隨機共振檢測模型是

式中：a，b為系統(tǒng)參數(shù)，sn（t）為待測含噪信號。將上述實驗采集數(shù)據(jù)代入模型（6）進行處理。由于圖5顯示信號幅值很低，能量不足以引起粒子在兩勢阱間的躍遷，所以首先對采集數(shù)據(jù)進行幅值放大處理，將原始含噪信號幅值同比增大 200倍，即 sn′（t）＝200×sn（t）。其次，根據(jù)產(chǎn)生隨機共振的小參數(shù)特性，按照文獻［5－6］變尺度方法對原始數(shù)據(jù)進行二次采樣，取二次采樣頻率 fsr＝5 Hz。令參數(shù) a＝1，b＝10，將經(jīng)過上述處理的數(shù)據(jù)代入式（6），用四階Runge-Kutta法計算，并進行尺度還原得到系統(tǒng)輸出時域圖和低頻頻譜如圖7所示。

由圖7可以看出，在相應(yīng)參數(shù)條件下系統(tǒng)為欠共振狀態(tài)，頻譜圖上不易看出明顯的特征信號。嘗試改變系統(tǒng)參數(shù)a＝1，b＝50，重新計算得到系統(tǒng)輸出時域圖和低頻頻譜如圖8所示。經(jīng)過參數(shù)的再次調(diào)節(jié)，圖8的頻譜可以分辨出故障特征信號成分。

圖8 轉(zhuǎn)子軸故障信號經(jīng)過非線性雙穩(wěn)系統(tǒng)處理后的時域及頻域圖（參數(shù)a＝1，b＝50）Fig．8 The time domain and frequency domain of rotor shaft fault signal after being processed in non-linear bistable system．（parameter a＝1，b＝50）

從兩種檢測方法的對比中可以看出，利用非線性系統(tǒng)隨機共振法檢測微弱信號時，多參數(shù)的協(xié)調(diào)選擇是一個難點，這些參數(shù)調(diào)節(jié)一般需要根據(jù)經(jīng)驗選取，多次嘗試后才能得到較理想的結(jié)果。相比較用一階線性系統(tǒng)的調(diào)參隨機共振法檢測，采樣頻率選定后只需要調(diào)節(jié)系統(tǒng)參數(shù)a，就可以得到有效結(jié)果，因此該方法在實際應(yīng)用上更簡便有效。

4 結(jié) 論

本文研究了一階線性系統(tǒng)調(diào)參廣義隨機共振機理，及其用于檢測微弱故障信號的處理方法。周期信號與高斯白噪聲作用下的一階線性系統(tǒng)，通過調(diào)節(jié)系統(tǒng)參數(shù)，在輸出端信噪比能取到極大值，但這并不能保證輸出譜圖上特征信號能達到理想的識別性。為在輸出譜圖上得到更清晰的特征信號，文中討論了信號特征頻率的可辨識性與系統(tǒng)參數(shù)、信號頻率、采樣頻率等參數(shù)之間的變化關(guān)系，給出實際應(yīng)用時各參數(shù)取值的參考圖。隨后滑動軸承試驗臺上轉(zhuǎn)子軸彎曲故障的實例分析，驗證了一階線性系統(tǒng)調(diào)參共振檢測早期微弱故障信號的有效性與簡便性。

［1］張文斌，周曉軍，楊先勇，等．基于諧波窗方法的轉(zhuǎn)子軸心軌跡提純［J］．振動與沖擊，2009，28（8）：74－77．ZHANG Wen-bin，ZHOU Xiao-jun，YANG Xian-yong，et al．Harmonic window method for purification of axis trace［J］．Journal of Vibration and Shock，2009，28（8）：74－77．

［2］陳仁祥，湯寶平，呂中亮．EEMD濾波的轉(zhuǎn)子軸心軌跡提純方法［J］．重慶大學(xué)學(xué)報，2012，35（11）：15－20．CHEN Ren-xiang， TANG Bao-ping， L Zhong-liang． A mthod of rotor orbit purification based on ensemble empirical mode decomposition filter［J］． Journal of Chongqing University，2012，35（11）：15－20．

［3］楊文志，馬文生，任學(xué)平．小波包降噪方法在滑動軸承故障診斷中的應(yīng)用研究［J］．噪聲與振動控制，2009，29（4）：50－53．YANG Wen-zhi， MA Wen-sheng， REN Xue-ping．Application of wavelet packet denoising method in fault diagnosis of bearings［J］．Noise and Vibration Control，2009，29（4）：50－53．

［4］石鵬，冷永剛，范勝波，等．雙穩(wěn)系統(tǒng)處理微弱沖擊信號的研究［J］．振動與沖擊，2012，31（6）：150－154．SHI Peng，LENG Yong-gang，F(xiàn)AN Sheng-bo，et al．A bistable system for detecting a weak pulse signal［J］．Journal of Vibration and Shock，2012，31（6）：150－154．

［5］Leng Y G，Wang T Y，Guo Y，et al．Engineering signal processing based on bistable stochastic resonance［J］．Mechanical Systems and Signal Processing，2007，21：138－150．

［6］冷永剛．雙穩(wěn)調(diào)參高頻共振機理［J］．物理學(xué)報，2011，60（2）：020503＿1－7．LENG Yong-gang．Mechanism of high frequency resonance of parameter-adjusted bistable system ［J］． Acta Physica Sinica，2011，60（2）：020503＿1－7．

［7］Berdichevsky V，Gitterman M．Stochastic resonance in linear systems subject to multiplicative and additive noise［J］．Phys．Rev．E，1999，60（2）：1494－1499．

［8］Gitterman M．Harmonic oscillator with multiplicative noise：Nonmonotonic dependence on the strengthand the rate of dichotomous noise［J］．Phys．Rev．E，2003，67：057103．

［9］張莉，劉立，曹力．過阻尼諧振子的隨機共振［J］．物理學(xué)報，2010，59（3）：1494－1498．ZHANG Li，LIU Li，CAO Li．Stochastic resonance in an overdamped harmonic oscillator［J］．Acta Physica Sinica，2010，59（3）：1494－1498．

［10］陸志新，曹力．輸入方波信號的過阻尼諧振子的隨機共振［J］．物理學(xué)報，2011，60（11）：110501．LU Zhi-xin，CAO Li．Stochastic resonance of square wave signal in an overdamped harmonic oscillator［J］．Acta Physica Sinica，2011，60（11）：110501．

［11］田祥友，冷永剛，范勝波．一階線性系統(tǒng)的調(diào)參隨機共振研究［J］．物理學(xué)報，2013，62（2）：020505．TIAN Xiang-you， LENG Yong-gang， FAN Sheng-bo．Parameter-adjusted stochastic resonance of first-order linear system［J］．Acta Physica Sinica，2013，62（2）：020505．

［12］Mitaim S，Kosko B．Adaptive stochastic resonance［J］．Proceedings of the IEEE，1998，86：2152－2183．

［13］ Gingl Z，Vajtai R，Kiss L B．Signal-to-noise ratio gain by stochastic resonance in a bistable system［J］．Chaos，Solitons＆Fractals，2000，11：1929－1932．