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

張文海，胡明輝*，江志農，馮 坤

(1.北京化工大學 發(fā)動機健康監(jiān)控及網絡化教育部重點實驗室，北京 100029；2.北京化工大學 高端機械裝備健康監(jiān)控與自愈化北京市重點實驗室，北京 100029)

0 引 言

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

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

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

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

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

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

1 方法基礎理論

1.1 同步提取變換(SET)

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

主要實現過程分如下3步。

(1)STFT時頻譜計算

計算信號x(t)=A·eiω0t的標準STFT時頻譜：

(1)

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

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

(2)

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

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

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

(3)

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

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

(3)瞬時頻率精確提取

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

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

(4)

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

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

原始信號與SET譜的重構關系如下：

(5)

式中：Re{·}—取實部。

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

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

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

(1)數值方程

設一個振動信號波形x(n)，該信號與其復包絡信號x′(n)、幅值A(n)三者的關系[12]如下：

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

(6)

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

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

VKFOT數值方程向量形式如下：

y-Cx=η

(7)

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

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

(2)結構方程

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

(8)

結構方程向量形式為：

Ax=ε

(9)

(3)求解過程

VKFOT的求解過程是通過結構方程與數據方程聯立，將其轉化為最小均方根問題。

通過最小二乘法建立異構項與誤差項的平方和，稱為損失函數J：

J=r2εTε+ηTη

(10)

式中：r—權重因子，用于形容平衡誤差項η與異構項ε占比。

權重因子r與跟蹤帶寬Bw的關系為：

(11)

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

由于J的展開式中二次項系數為正值，求J對于xH的一階導數為0時J取到最小值。求解可得：

x=(r2ATA+E)-1CHy

(12)

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

2 SET與VKFOT結合的階次跟蹤方法

筆者提出的SET與VKFOT結合的階次跟蹤方法實現流程，如圖1所示。

圖1 SET與VKFOT結合的階次跟蹤方法

具體步驟如下：

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

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

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

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

(5)聯立數值方程、結構方程、權重因子，構造損失函數并求其最小值，求解得到目標階次時域波形。

3 實際案例分析

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

(1)數據基本情況

某型雙轉子燃氣輪機啟動及平穩(wěn)運行過程級聯圖如圖2所示。

圖2 啟動過程級聯圖A—振動加速度；v—轉速；f—頻率

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

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

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

圖3 SET時頻圖T—時間

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

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

圖4 實測轉速與SET提取轉速對比

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

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

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

得到目標階次成分時域波形如圖5所示。

圖5 57階目標階次時域波形

(4)結果對比分析

筆者單獨使用VKFOT進行分析，以原始轉速及振動信號作為濾波器輸入，設定相同的目標階次與濾波器參數，得到目標階次成分時域波形。

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

圖6 階次幅值趨勢對比

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

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

4 結束語

針對燃氣輪機轉速波動較大、測量轉速與振動難以同步，進而導致難以提取振動信號中目標階次幅值變化趨勢的問題，筆者提出了一種同步提取變換與Vold-Kalman結合的燃氣輪機振動階次跟蹤方法；以機匣振動速度信號為輸入，通過同步提取變換方法獲得與振動同步的轉速數據；以振動及與其同步的轉速數據為輸入，使用基于Vold-Kalman濾波的階次跟蹤方法提取振動信號中目標階次時域波形。

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