基于同步提取變換與Vold-Kalman濾波的燃機(jī)階次跟蹤方法研究*

張文海，胡明輝*，江志農(nóng)，馮 坤

(1.北京化工大學(xué) 發(fā)動機(jī)健康監(jiān)控及網(wǎng)絡(luò)化教育部重點(diǎn)實(shí)驗(yàn)室，北京 100029；2.北京化工大學(xué) 高端機(jī)械裝備健康監(jiān)控與自愈化北京市重點(diǎn)實(shí)驗(yàn)室，北京 100029)

0 引 言

目前，燃?xì)廨啓C(jī)已在我國航空、船舶及石油化工等行業(yè)得到了廣泛應(yīng)用，但其在運(yùn)行過程中故障頻發(fā)。當(dāng)燃機(jī)發(fā)生故障時，振動是最敏感的參數(shù)之一。例如轉(zhuǎn)子不平衡、支承不同心及葉片斷裂等故障發(fā)生后，由于對轉(zhuǎn)子狀態(tài)或流場的影響，測得的機(jī)匣振動信號必然發(fā)生變化。

因此，振動分析已經(jīng)被作為燃機(jī)狀態(tài)監(jiān)測的重要手段。然而，目前我國主要通過機(jī)匣的低頻振動速度總值來監(jiān)測燃機(jī)振動狀態(tài)，該參數(shù)可從振動能量角度較好地用于聯(lián)鎖保護(hù)，卻難以實(shí)現(xiàn)振動故障的早期預(yù)警，亦難以向故障診斷提供直接的數(shù)據(jù)支撐。

由于燃?xì)廨啓C(jī)常見的不平衡、不對中、不同心、葉片斷裂等故障，其主要振動特征均為跟轉(zhuǎn)速相關(guān)的諧頻振動即階次成分，因此，燃?xì)廨啓C(jī)振動信號中的階次分量是其振動故障預(yù)警和故障分析的主要依據(jù)。在燃機(jī)的運(yùn)行過程中，包括穩(wěn)定工況和啟停車等非穩(wěn)定工況，有必要對階次振動進(jìn)行跟蹤提取。基于卡爾曼濾波的階次跟蹤(vold-kalman filter based order tracking，VKFOT)是旋轉(zhuǎn)機(jī)械振動階次分量提取方面應(yīng)用較廣的方法之一，其不僅可以用于計算階次分量幅值，還可以得到各階次分量的時域波形，為后續(xù)相位分析等階次分量深入挖掘提供足夠的數(shù)據(jù)支撐。但該方法需要振動信號及與振動同步的轉(zhuǎn)速數(shù)據(jù)同時作為輸入[1]。

然而，由于燃機(jī)轉(zhuǎn)速波動較大，同時很多現(xiàn)場的轉(zhuǎn)速數(shù)據(jù)和振動數(shù)據(jù)往往通過不同設(shè)備進(jìn)行采集，即很難獲得與振動完全同步的精確轉(zhuǎn)速數(shù)據(jù)。在機(jī)匣處測得的振動信號中，工頻成分與直接測量得到的轉(zhuǎn)速頻率存在差異。若將該轉(zhuǎn)速與振動信號作為VKFOT方法的輸入，則必然會因?yàn)檗D(zhuǎn)頻偏差而導(dǎo)致提取到的階次幅值誤差較大，對于高階成分誤差更大。

1995年，Vold和Leuridan首次將經(jīng)典最優(yōu)濾波理論中Kalman濾波器應(yīng)用于階次跟蹤上[2]，但該方法需要振動信號及與振動同步且相同采樣率的轉(zhuǎn)速數(shù)據(jù)同時作為濾波器輸入，否則便會引發(fā)提取到的頻率分量與目標(biāo)不一致問題。針對轉(zhuǎn)速與振動同步性較差或無轉(zhuǎn)速輸入，學(xué)者們首先提出通過設(shè)置合理的濾波器帶寬來補(bǔ)償?shù)姆椒ǎM(jìn)而確保階次提取結(jié)果的準(zhǔn)確性。TUMA[3]通過理論分析，給出了濾波帶寬與共振頻率的關(guān)系；WANG K S等人[4]通過建模仿真得到了共振頻率計算濾波的優(yōu)化帶寬；馮珂等人[5]通過VKF跟蹤的多次循環(huán)計算來進(jìn)行濾波帶寬的優(yōu)化選取。以上方法均可用于補(bǔ)償目標(biāo)階次偏差帶來的VKF跟蹤誤差，但僅適用于結(jié)構(gòu)較為簡單、共振頻率易于獲得的場合，循環(huán)優(yōu)選也不能從根本上解決目標(biāo)階次偏移的問題。另一方面，為獲得同步的轉(zhuǎn)速數(shù)據(jù)，將VKFOT方法應(yīng)用于無轉(zhuǎn)速信號的場合，文獻(xiàn)應(yīng)用了能量重心法頻譜校正、短時傅里葉變換(STFT)及短時chirp-Fourier變換(STCFT)等獲得瞬時轉(zhuǎn)速[6-10]，但由于以上時頻變換方法建立在局部平穩(wěn)性假設(shè)上，不適用于燃機(jī)變負(fù)荷過程轉(zhuǎn)速變化劇烈時的轉(zhuǎn)速提取。于剛[11]在2017年提出了同步提取變換(SET)方法，極大地提高了非平穩(wěn)非線性信號的瞬時頻率分析精度，為快變轉(zhuǎn)速場合轉(zhuǎn)速提取提供了可能。

針對上述問題，為獲取與振動同步的轉(zhuǎn)速數(shù)據(jù)，筆者提出一種基于同步提取變換與Vold-Kalman濾波的燃機(jī)振動階次跟蹤方法。

1 方法基礎(chǔ)理論

1.1 同步提取變換(SET)

同步提取變換是基于短時傅里葉變化(STFT)的時頻分析方法，用于提取瞬時頻率(IF)。

主要實(shí)現(xiàn)過程分如下3步。

(1)STFT時頻譜計算

計算信號x(t)=A·eiω0t的標(biāo)準(zhǔn)STFT時頻譜：

(1)

式中：g(u-t)—可移動的緊支窗函數(shù)。

根據(jù)帕塞瓦爾定理，對式(1)進(jìn)行修正，進(jìn)一步展開可得：

(2)

由式(2)可知，STFT是由一系列與諧波信號頻率ω0一致的時頻系數(shù)組成。在頻率ω=ω0處，時頻系數(shù)幅值最大，距離ω0越遠(yuǎn)時頻系數(shù)幅值越小。

(2)STFT中的瞬時頻率估計

根據(jù)時頻譜的相位信息，可得到瞬時頻率：

(3)

式中：?tGe(t,ω)—Ge(t,ω)對時間t的一階偏導(dǎo)數(shù)。

式(3)表示在二維時頻平面中，只要存在Ge(t,ω)≠0，ω0(t,ω)即為瞬時頻率ω0。

(3)瞬時頻率精確提取

利用同步提取算子將瞬時頻率處的發(fā)散能量去除，得到SET譜：

Te(t,ω)=Ge(t,ω)·δ(ω-ω0(t,ω))

(4)

式中：δ(ω-ω0(t,ω))—同步提取算子。

式(4)通過δ函數(shù)(單位脈沖函數(shù))，僅保留時頻脊線上能量最大(ω=ω0)的時頻系數(shù)，從而得到時頻聚焦良好的時頻分析效果。

原始信號與SET譜的重構(gòu)關(guān)系如下：

(5)

式中：Re{·}—取實(shí)部。

因此，r(t)為SET譜中信號的瞬時頻率軌跡，即實(shí)現(xiàn)了在高分辨率SET譜提取瞬時頻率。

1.2 Vold-Kalman階次跟蹤方法原理

VKFOT的基本思想是使用轉(zhuǎn)速數(shù)據(jù)構(gòu)建跟蹤階次的結(jié)構(gòu)方程作為估計值，以傳感器測量信號形成的數(shù)據(jù)方程作為觀測值，通過令估計誤差與觀測誤差協(xié)方差最小的方法，估計出該狀態(tài)下的跟蹤階次信號。

(1)數(shù)值方程

設(shè)一個振動信號波形x(n)，該信號與其復(fù)包絡(luò)信號x′(n)、幅值A(chǔ)(n)三者的關(guān)系[12]如下：

x(n)=A(n)cos[Θ(n)+φ]=
Re{x′(n)exp(jΘ(n))}

(6)

式中：φ—初相位；Θ(n)—參考軸的角位移。

對于VKFOT數(shù)值方程而言，階次跟蹤就是求解復(fù)包絡(luò)x′(n)的過程，通過式(6)可以將幅值、復(fù)包絡(luò)與波形進(jìn)行相互轉(zhuǎn)換。

VKFOT數(shù)值方程向量形式如下：

y-Cx=η

(7)

式中：y—原振動信號；x—目標(biāo)階次波形，二者均為復(fù)包絡(luò)形式；C(n)—由轉(zhuǎn)速計算求得的角位移；η—誤差項(xiàng)，代表其他諧波分量與隨機(jī)噪聲成分的總和。

由于在數(shù)值方程中引入角位移，該方法被稱為基于角位移的VKF_OT[13]。

(2)結(jié)構(gòu)方程

用包含異構(gòu)項(xiàng)ε(n)的低階多項(xiàng)式表示的結(jié)構(gòu)方程如下(異構(gòu)項(xiàng)ε(n)形容目標(biāo)階次受到噪聲干擾而引起的幅值頻率等變化)：

(8)

結(jié)構(gòu)方程向量形式為：

Ax=ε

(9)

(3)求解過程

VKFOT的求解過程是通過結(jié)構(gòu)方程與數(shù)據(jù)方程聯(lián)立，將其轉(zhuǎn)化為最小均方根問題。

通過最小二乘法建立異構(gòu)項(xiàng)與誤差項(xiàng)的平方和，稱為損失函數(shù)J：

J=r2εTε+ηTη

(10)

式中：r—權(quán)重因子，用于形容平衡誤差項(xiàng)η與異構(gòu)項(xiàng)ε占比。

權(quán)重因子r與跟蹤帶寬Bw的關(guān)系為：

(11)

式中：Δf=Bw/fs，fs—信號采樣頻率。

由于J的展開式中二次項(xiàng)系數(shù)為正值，求J對于xH的一階導(dǎo)數(shù)為0時J取到最小值。求解可得：

x=(r2ATA+E)-1CHy

(12)

求解該方程即可得到跟蹤諧波分量x′(n)。

2 SET與VKFOT結(jié)合的階次跟蹤方法

筆者提出的SET與VKFOT結(jié)合的階次跟蹤方法實(shí)現(xiàn)流程，如圖1所示。

圖1 SET與VKFOT結(jié)合的階次跟蹤方法

具體步驟如下：

(1)由于加速度信號中主要成分為處于高頻段的葉片通過頻率等成分，成分多、波動大，難以直接提取，而工頻成分不明顯，對加速度信號進(jìn)行積分得到工頻成分占主導(dǎo)的速度信號；

(2)對(1)中得到的速度信號進(jìn)行同步提取變換，得到包含多個IF分量的SET時頻圖；

(3)結(jié)合采集所得的同步性較差轉(zhuǎn)速信號，從(2)中SET時頻圖中確定工頻脊線明顯、無噪聲成分的頻率范圍作為工頻提取帶寬，從SET時頻圖中提取工頻信號，即可得到與振動同步的轉(zhuǎn)速信號；

(4)根據(jù)跟蹤目標(biāo)階次的大小，選擇振動加速度信號或速度信號，以及SET提取得到的轉(zhuǎn)速數(shù)據(jù)作為VKF濾波器輸入，分別構(gòu)造數(shù)值方程和結(jié)構(gòu)方程；根據(jù)跟蹤目標(biāo)設(shè)定目標(biāo)階次O，如跟蹤一階分量則取O=1；根據(jù)時頻圖中目標(biāo)階次處無其他噪聲成分最小范圍確定跟蹤帶寬Bw；根據(jù)式(10)計算權(quán)重因子r；

(5)聯(lián)立數(shù)值方程、結(jié)構(gòu)方程、權(quán)重因子，構(gòu)造損失函數(shù)并求其最小值，求解得到目標(biāo)階次時域波形。

3 實(shí)際案例分析

為驗(yàn)證該方法在燃機(jī)振動階次分量提取中的有效性及優(yōu)勢，筆者選取了某型燃?xì)廨啓C(jī)實(shí)際試車數(shù)據(jù)進(jìn)行分析；使用該方法和原始VKFOT方法分別提取目標(biāo)階次時域波形，采用計算階次跟蹤(COT)的方法，分別將二者提取所得的目標(biāo)階次幅值與原信號階次譜中讀出的幅值進(jìn)行對比分析。

(1)數(shù)據(jù)基本情況

某型雙轉(zhuǎn)子燃?xì)廨啓C(jī)啟動及平穩(wěn)運(yùn)行過程級聯(lián)圖如圖2所示。

圖2 啟動過程級聯(lián)圖A—振動加速度；v—轉(zhuǎn)速；f—頻率

振動信號采樣頻率為fs=51 200 Hz，信號長度T=12 s。圖2中標(biāo)出的明顯階次成分頻率8 326 Hz，為低壓轉(zhuǎn)子的57.3階(低壓轉(zhuǎn)子直接測量轉(zhuǎn)速為8 718 r/min)，對應(yīng)低壓壓氣機(jī)某級動葉通過頻率(57階)，但按照測量轉(zhuǎn)速計算得到的階次為57.3，與實(shí)際階次數(shù)57存在較大偏差。下面對這一明顯階次成分進(jìn)行階次提取。

(2)基于振動信號的轉(zhuǎn)速提取

筆者對振動加速度信號進(jìn)行積分得到速度信號，進(jìn)行同步提取變換后得到包含低壓工頻的SET時頻，如圖3所示。

圖3 SET時頻圖T—時間

從圖3可以看出：低壓轉(zhuǎn)子工頻脊線明顯，周圍約6 Hz處無噪聲成分。以直接測量得到的轉(zhuǎn)速為中心頻率，帶寬6 Hz，從SET時頻圖中提取得到低壓轉(zhuǎn)子工頻。

提取得到的轉(zhuǎn)速與直接測量轉(zhuǎn)速對比圖如圖4所示。

圖4 實(shí)測轉(zhuǎn)速與SET提取轉(zhuǎn)速對比

從圖4中可以看出：通過SET譜提取的轉(zhuǎn)速與直接測量所得轉(zhuǎn)速存在一定偏差。當(dāng)跟蹤高階成分時，該誤差會呈倍數(shù)放大。

(3)VKFOT目標(biāo)階次時域波形提取

將振動信號，SET譜提取到的工頻信號作為VKF濾波器輸入，設(shè)定目標(biāo)階次O為57階，由于通過SET譜提取的轉(zhuǎn)頻成分避免了階次偏移的問題，根據(jù)時頻圖中目標(biāo)階次處無其他噪聲成分最小范圍，確定濾波帶寬Bw分別為0.85個階次，最大程度上保留有效階次成分，濾除噪聲干擾。

得到目標(biāo)階次成分時域波形如圖5所示。

圖5 57階目標(biāo)階次時域波形

(4)結(jié)果對比分析

筆者單獨(dú)使用VKFOT進(jìn)行分析，以原始轉(zhuǎn)速及振動信號作為濾波器輸入，設(shè)定相同的目標(biāo)階次與濾波器參數(shù)，得到目標(biāo)階次成分時域波形。

分別對原始振動信號、SET-VKFOT得到的57階時域波形、單獨(dú)VKFOT得到的時域波形以原始轉(zhuǎn)速為參考做階次譜，分別讀取目標(biāo)階次成分(即57階)幅值，得到階次幅值趨勢對比圖，如圖6所示。

圖6 階次幅值趨勢對比

由圖6可知：SET-VKFOT方法得到的階次幅值趨勢與原階次譜幅值趨勢基本吻合，而單獨(dú)使用VKFOT則誤差較大。計算圖中兩種方法與原階次譜幅值趨勢的絕對誤差，SET-VKFOT誤差在15 m/s2以內(nèi)，而VKFOT誤差最大達(dá)到了160 m/s2以上。

由上述對比分析可知：在變負(fù)荷及穩(wěn)定工況下，SET-VKFOT方法提取階次分量幅值誤差較小，相對單獨(dú)使用VKFOT方法有明顯優(yōu)勢。

4 結(jié)束語

針對燃?xì)廨啓C(jī)轉(zhuǎn)速波動較大、測量轉(zhuǎn)速與振動難以同步，進(jìn)而導(dǎo)致難以提取振動信號中目標(biāo)階次幅值變化趨勢的問題，筆者提出了一種同步提取變換與Vold-Kalman結(jié)合的燃?xì)廨啓C(jī)振動階次跟蹤方法；以機(jī)匣振動速度信號為輸入，通過同步提取變換方法獲得與振動同步的轉(zhuǎn)速數(shù)據(jù)；以振動及與其同步的轉(zhuǎn)速數(shù)據(jù)為輸入，使用基于Vold-Kalman濾波的階次跟蹤方法提取振動信號中目標(biāo)階次時域波形。

筆者使用該方法對實(shí)際案例數(shù)據(jù)進(jìn)行了分析，結(jié)果表明，該方法能夠避免轉(zhuǎn)速不精確、不同步帶來的轉(zhuǎn)頻偏差問題，能夠準(zhǔn)確提取目標(biāo)階次時域波形。由此可見，該方法適用于在燃?xì)廨啓C(jī)運(yùn)行過程中跟蹤振動階次分量，可以為振動狀態(tài)監(jiān)測及預(yù)警提供數(shù)據(jù)支撐。