離心壓縮機出口動態壓力的最大Lyapunov指數特征研究

劉 雁，高 寬，何 浩，黃 炎，肖 軍

(1.西北工業大學 機電學院,西安 710072；2.合肥通用機械研究院 壓縮機技術國家重點實驗室 ,合肥 230031)

喘振是離心壓縮機固有的一種典型故障，由于存在氣流的劇烈脈動和周期性振蕩，容易加劇壓縮機的磨損、使其機械部件發生變形。這會給離心壓縮機帶來巨大的危害性，甚至引發重大安全事故。因此，在工業現場，實現離心壓縮機喘振的早期檢測非常必要。

目前，對離心壓縮機喘振的研究主要從兩個方面展開：①基于數學模型法，即Greitzer[1]集中參數模型及其衍生的相關模型，主要研究壓縮機結構參數的變化與系統穩定性之間的關聯關系；②基于試驗測試的數據分析法，由于離心壓縮機的工作狀態與壓縮機的出、入口及葉輪處的動態壓力密切相關，且壓縮機是否發生失穩可以通過動態壓力的變化反映出來。因此，通過分析離心壓縮機不同工作狀態下的動態壓力的特征量，分析出當前壓縮機的工作狀態，進而預測壓縮機后期的工作狀態。近年來，研究人員通過分析進口、出口和葉輪的動態壓力研究喘振現象。其中，孫濤等[2]用關聯積分值定量描述了離心壓縮機的失穩和喘振程度。Xue等[3]提取了一個表征離心壓縮機轉速瞬態混沌分岔確定性的特征指標。Sun等[4]演示了通過聲學信號捕捉失穩現象的一種方法，并提供了一種聲學測量技術來檢測壓縮機整個生命周期中的任意時刻的失穩狀態。陳振等[5]提出了一種新的對旋轉失速和喘振的非線性反饋控制策略，即把壓升作為反饋信號，來補償節流閥系數值。Liu等[6-8]分別從分形動力學和熵的角度研究了離心壓縮機出口壓力的非線性特征，并獲得了一些重要的結論。這些研究表明，離心壓縮機的出口動態壓力存在豐富、復雜的非線性特征，因此可以將非線性分析方法應用于離心壓縮機喘振的研究中。

Lyapunov指數作為一類重要的非線性參數，通常用來定量地描述非線性系統的混沌程度，近些年來被廣泛地用于復雜系統的動力學特性研究和數值預測中。如Deng等[9]采用最大Lyapunov指數分析了水輪發電機的隨機全局穩定性和分岔問題。Saikia等[10]研究了Lyapunov指數在帕金森病腦電圖診斷中的意義。Wu等[11]將小波變換、粒子群優化動態灰色模型和Lyapunov指數預測方法相結合，并將其應用于短期風電場的預測。Picano等[12]通過對最大Lyapunov指數的研究，設計了智慧城市的乘客需求預測方法。伍友利等[13-14]研究某型發動機的八級軸流壓氣機級間壓力信號的最大Lyapunov指數，推斷出失速和喘振的起始點具有相同的動力學特性，提出將最大Lyapunov指數的零點作為失速的判據。劉中勝等[15]研究了風力機地震誘導機艙振動非線性特性分析，發現地震誘導振動時程曲線最大Lyapunov指數與地震強度之間存在較大的相關性。劉雁等[16]采用最大Lyapunov指數對局域網絡流量的變化趨勢進行了預測，具有一定的效果。鑒于離心壓縮機出口動態壓力中包含豐富的非線性特征，從非線性角度研究離心壓縮機出口動態壓力的最大Lyapunov指數與出口壓力實時序列之間的關系，具有理論及實踐意義。

1 相關理論

1.1 時間序列的相空間重構

在采用最大Lyapunov 指數對離心壓縮機出口動態壓力進行分析之前，需重構出口動態壓力信號的相空間。由于系統的任意分量的演化都是由系統中與之相關的分量決定的，系統的信息一定會隱含在任意分量的發展過程中，因此通過考察一個分量，對其在一定時間間隔下的測量值進行分析，可以提取系統的非線性特征。20 世紀80年代，Takens[17]提出著名的Takens定理，該定理是相空間重構的理論基礎。根據Takens定理，對于一維時間序列X={x1,x2, …,xN}，其相空間重構的基本過程如下：

步驟1選擇合適的重構參數，即嵌入維數m和時間延遲τ。

步驟2從x1開始取值，每一個值的往后延遲時間為τ，共取m個值，則得到m維相空間中的第1個點

Y1=(x1,x1+τ,x1+2τ,…,x1+(m-1)τ)

步驟3去掉x1，以x2為第一個數，采用步驟2中的方法得到相空間的第2個點

Y2=(x2,x2+τ,x2+2τ,…,x2+(m-1)τ)

步驟4對于長度為N的時間序列，依次可得到由M個相點構成的m維相空間

Y1=(x1，x1+τ，…，x1+(m-1)τ)

Y2=(x2，x2+τ，…，x2+(m-1)τ)

…

YM=(xM，xM+τ，…，xM+(m-1)τ)

(1)

式中，M=N-(m-1)τ，Yi(i=1,2,…,M)為重構m維相空間中的相點。這樣，相點間的連線就可以描述系統在m維相空間中的演化軌跡，該軌跡也稱為重構吸引子，顯然，重構吸引子中總共有M個點。也就是說，經過相空間重構，由一維N個點的時間序列得到了m維的重構吸引子。

實現相空間重構的重要的一個環節在于步驟1，即選擇最優的重構參數嵌入維數m和時間延遲τ。因為，所構建相空間中的吸引子應具有較低的重復性和較強的相關性，因此在選擇重構參數的過程中，需要對重構吸引子幾何上的重復性和相關性做出評估。確定重構參數的方法有多種，主要分為兩類：第一類，認為兩個重構參數的選取是獨立的，如自相關函數法和互量信息法等，這些方法通過評估重構吸引子的自相關性、總體關聯量等信息，分別確定兩個參數；第二類，認為τ和m的選取是相互依賴的，如C-C算法，這種方法是通過統計的手段，同時確定兩個參數[18]。由于C-C算法具有適用于小數據量，計算方便的優點，本文中擬采用C-C算法計算重構參數延遲時間τ和嵌入維數m。

1.2 最大Lyapunov指數的物理含義

Lyapunov指數是定量分析混沌運動的特征參數。混沌運動可能存在于世間一切的復雜運動中，其基本特點就是對初值極為敏感，在相空間中，兩個很靠近的初值所產生的兩條相軌跡，會隨時間的推移按照指數方式分離，Lyapunov指數就是定量描述這一現象的物理量。

設F是Rm→Rm上的m維映射，它可以確定一個m維離散動力系統

xn+1=F(xn)

(2)

假設該系統的初始條件為一個足夠小的m維的小球，經過時間t的演化，小球將變為橢球，將橢球的所有主軸按其長度順序排列，則第i個主軸所對應的Lyapunov指數可以根據第i個主軸的長度Pi(t)的增加速率定義為

(3)

將所有的Lyapunov指數按大小排列為

λ1≥λ2≥λ3≥…≥λm

(4)

式中，λ1為系統的最大Lyapunov指數。

Lyapunov指數可能為正，也可能為負。Lyapunov指數為正表示系統的演化軌跡在對應方向上是發散的。Lyapunov指數為負表示系統的演化軌跡在對應方向上是收縮的。如果所有的Lyapunov指數均為負，則系統將趨于靜止；如果最大Lyapunov指數為零而其余Lyapunov指數為負，則系統做周期運動；如果存在正Lyapunov指數，則系統做混沌運動。因此，識別系統是否做混沌運動只需判斷最大Lyapunov指數是否為正，而不需要計算出所有的Lyapunov指數，這樣就可以極大的減少計算量。并且，最大Lyapunov指數越大，系統中的混沌特性越明顯，混沌程度越高。

1.3 最大Lyapunov指數的計算

對于時間序列重構的相空間，常用的計算最大Lyapunov指數方法有定義法、Jacobian方法、Wolf算法、P-范數法、小數據量法、奇異值分解法等。其中，小數據量法具有小數據組計算可靠，容易操作，計算工作量比較小等優點，本文采用小數據量法計算離心壓縮機出口動態壓力的最大Lyapunov指數λ1。其計算過程如下[19]：

在如式(1)所示的重構相空間中，尋找給定軌跡上點Yj的最近鄰近點Yk，并限制短暫分離，即

(5)

式中，p為時間序列的平均周期，可以通過能量譜平均頻率的倒數估算，也可采用平均延遲時間窗口τw代替p。τw可由式(6)計算

τw=(m-1)τ

(6)

經過i個離散時間步長后，計算距離

dj(j)=‖Yj+i-Yk+i‖i=1,2,…,min(m-j,m-k)

(7)

假設參考點Yj與其近鄰點Yk之間的距離具有指數發散的特性，那么

dj(i)=Cjeλ1(iΔt),Cj=dj(0)

(8)

將式(8)取對數得

lndj(i)=lnCj+λ1(iΔt),j=1,2,…,N

(9)

由式(9)可知，lndj(i)～i線在一定范圍內滿足線性關系，其斜率為λ1Δt。計算所有j對應的lndj(i)，并求其平均值，然后作出曲線

(10)

式中，q為非零dj(i)的個數。最大Lyapunov指數相當于式(10)表示的一簇直線的斜率，通常在y(i)曲線的線性區間內采用最小二乘法進行擬合，該區間內的斜率即為最大Lyapunov指數。

1.4 最大Lyapunov指數與離心壓縮機出口動態壓力

若將離心壓縮機出口管道氣體運動的數學模型作為一個動力系統，則出口動態壓力為這個動力系統中的一個狀態變量。可以認為，該系統是一個復雜系統，其中可能包含混沌運動、周期運動和趨于靜止的運動。最大Lyapunov指數的大小可以反映系統在其動態行為中的狀態。當最大Lyapunov指數較大時，表明所分析時段的出口動態壓力在一定范圍內處于無序波動，此時該系統主要表現為混沌運動，而且最大Lyapunov指數越大，系統表現出的不確定性、不可重復、不可預測越為明顯，其隨機特性越強。而當最大Lyapunov指數接近于零時，出口動態壓力的周期性變強，表現為周而復始的振蕩現象。

2 離心壓縮機出口動態壓力采集及頻率特性

2.1 離心壓縮機出口動態壓力采集系統

數據采集系統的結構方案，如圖1所示。試驗在國內某鼓風機動力集團的離心壓縮機實驗室內完成，該實驗室的壓縮機由一臺800 kW的直流電動機提供動力。當壓縮機工作時，處于室溫狀態的大氣，首先經過進氣過濾器進入工作管路，大量的氣體存儲于風室中等待壓縮，然后進入離心壓縮機壓縮。壓縮后的氣體由兩個直徑分別為250 mm和100 mm的電動防喘振閥控制，經由出口管道及排氣消聲器排入大氣。

圖1 采集系統

系統出口動態壓力變化的采集由動態數據采集系統完成，所采集的壓力為氣體靜壓。動態數據采集系統的傳感器為Kulite的動態壓力傳感器，此傳感器的固有頻率高且具有加速度補償功能，安裝方便。考慮到數據采集的便捷性及可靠性，將傳感器安裝在距離壓縮機出口1 m處的出口管道上。動態數據采集系統為美國晶鉆公司的CoCo-80，此系統可用于振動采集和狀態監測。

為確保喘振試驗的安全性及試驗數據的可靠性，本試驗由壓縮機實驗室的專業技術人員完成試驗操作。技術人員通過控制防喘振閥的開度使離心壓縮機系統由正常工作狀態進入喘振狀態。

由于離心壓縮機在喘振時有噪聲過大、機殼軸承振動、出入口流量變化、轉速波動大等特點[20]，專業技術人員將根據預先設置的喘振線、現場噪聲以及出口動態壓力的變化曲線預估喘振工況。

試驗時壓縮機工作在額定轉速，試驗過程共635.5 s，本文選取150～250 s的數據進行分析，該時段包含穩態，過渡狀態和喘振狀態。以150 s作為起點，在該時段的初始狀態，離心壓縮機處于正常工作工況，從第42～52 s逐漸關閉防喘振閥，使離心壓縮機的工況逐漸進入喘振狀態，第53 s之后，離心壓縮機的工況完全進入喘振狀態。

2.2 出口動態壓力波形及頻率特性

動態壓力采集系統采集到的部分出口動態壓力波形圖，如圖2所示。波形總長為100 s。從圖2可知，出口動態壓力經歷了穩態、過渡過程和喘振三個階段。在20～30 s時段，出口動態壓力處于穩態階段，此時出口動態壓力幅值較小，在-0.5～0.5 kPa 內波動(見圖2(b))；在45～55 s時段，動態壓力逐漸增加，此時壓縮機由過渡狀態進入喘振，且出口動態壓力幅值在52～54 s迅速增大(見圖2(c))。在70～80 s時段，離心壓縮機系統完全進入喘振狀態，出口動態壓力波動值達到-1.7～1.7 kPa，且波形表現出周期振動(見圖2(d))。

圖2 離心壓縮機出口動態壓力

出口動態壓力波動的時域特征可以通過極差、方差、標準差等特征量來衡量。出口動態壓力波動的方差變化圖，如圖3所示。從圖3可知，離心壓縮機由穩態經過喘振時其出口動態壓力的變化。

圖3 出口動態壓力的方差

選取離心壓縮機在穩定運行、過渡過程和喘振三種工作狀態(20～30 s，50～55 s和70～80 s)的動態壓力序列做頻譜分析，其頻譜密度曲線如圖4所示。

結合圖2與圖4的波形可知，隨著失穩狀態的加深，系統進入喘振狀態，出口波形開始表現出規律振動，波動周期變長，此時的頻率分布也越來越向低頻(15～20 Hz)處集中，出口動態壓力逐漸呈現為一個周期運動，其規律性隨著失穩狀態的加深而增強，這一特性與壓縮機喘振的發生機理一致，因為壓縮機在發生喘振時，出口管道的氣體在軸向會出現周期性的低頻晃動。

圖4 出口動態壓力的頻譜圖

3 出口動態壓力的最大Lyapunov指數

3.1 出口動態壓力的最大Lyapunov指數

在計算離心壓縮機出口動態壓力的最大Lyapunov指數之前，需要進行相空間重構。選取出口動態壓力數據中1 s的數據為計算單元，即將100 s的數據分為100段，每段進行一次相空間重構，計算所得的重構參數如圖5所示。

圖5 相空間重構參數

從圖5可知，離心壓縮機正常工作時，延遲時間τ取值范圍小，嵌入維數m的值較小但波動量變化較大；在過渡過程中，延遲時間τ與嵌入維數m明顯的增加，且均出現劇烈的波動；在喘振時，延遲時間τ取值較大，而嵌入維數m的值和波動均比穩態時大。從出口動態壓力的相空間重構參數中延遲時間τ的計算結果可知，離心壓縮機在穩態時的動力特征比喘振時復雜，相空間中吸引子的演化規律相應也較為復雜。

以1 s的時間序列為分析單元，根據表1所示重構相空間的參數，離心壓縮機出口動態壓力在三種工作狀態(20～25 s，51～55 s，70～75 s)的最大Lyapunov指數計算圖,如圖6所示。根據式(10)，對圖6中每一條曲線選擇合適的線性區間，則區間內曲線的斜率即為對應時段的最大Lyapunov指數。由于在不同工作狀態下，離心壓縮機出口動態壓力的最大Lyapunov指數計算圖的線性區間不同，如果線性區間選擇不合適會影響最大Lyapunov指數的準確性和后續分析的合理性。由圖6可知，在穩態時線性區間較小，這里我們選擇0～30即可。而在過渡狀態及喘振時線性區間范圍較寬，選擇0～100更為合適。圖6所示時段的相空間重構參數和最大Lyapunov指數，如表1所示。

表1 動態壓力的混沌參數

圖6 Lyapunov指數計算圖

整個工作時段的最大Lyapunov指數,如圖7所示。由圖7可知，離心壓縮機在正常工作狀態時和發生喘振時的最大Lyapunov指數差別較大。在穩定狀態下，最大Lyapunov指數較大，系統的混沌性較強，可預測性較差，具有較強的蝴蝶效應。而且，在穩態時的最大Lyapunov指數波動明顯，表明系統的復雜度較高，因為系統內在的混沌行為的存在，使系統受到初始條件的影響較大；隨著壓縮機喘振的出現和加深，最大Lyapunov指數較小且變得穩定，并接近于0，系統的復雜度下降，混沌性減弱。可以判定此時出口管道中的氣體的主要波動為近似周期振動。

最大Lyapunov指數變化規律在47～48 s時發生了較大變化，由較大的波動范圍變為較小的波動范圍，由較大的平均值變為較小的平均值(見圖7)。根據圖3及現場專業技術人員的判斷，喘振的出現位于53～54 s。也就是說，在初始喘振前的過渡過程時的最大Lyapunov指數變化很大，針對此次測試處的出口動態壓力信號，我們可以提前約5 s對喘振狀態進行預警，但在實際應用過程中，預測喘振的發生點與閥門關閉速度存在關聯關系，需根據現場情況進行調整。

圖7 最大Lyapunov指數

3.2 離心壓縮機出口動態壓力的最大Lyapunov指數抗干擾性能分析

采用相空間重構和最大Lyapunov指數對離心壓縮機出口動態壓力進行分析時，一維信號會經過一個投影到高維相空間的過程。因此，在分析計算時，微小的干擾信號和底噪信號經過相空間重構，在能量上會進行分散，不會對Lyapunov指數隨喘振變化的規律造成較大影響。因此，該方法應具備抗干擾能力。為驗證這一結論，本文在采集的出口動態壓力的基礎上，再加入一定比例的噪聲信號，對照加入噪聲前后的最大Lyapunov指數計算結果，分析該方法的抗干擾性能。我們對原始出口動態壓力的時間序列中加入壓力峰值的5%的服從均方分布的隨機信號作為干擾信號。同樣以1 s為分析單元計算引入干擾信號的離心壓縮機出口動態壓力的最大Lyapunov指數，如圖8所示。

圖8 加入干擾后的最大Lyapunov指數

和圖7相比較，在加入干擾信號以后，出口動態壓力的最大Lyapunov指數變化趨勢幾乎一致。在離心壓縮機處于穩定時，最大Lyapunov指數的值較高，但波動比原始信號的最大Lyapunov指數要劇烈一些，表明引入干擾信號以后，系統的復雜度有所增加，混沌特征仍十分明顯；而離心壓縮機處于喘振狀態時，引入干擾信號對最大Lyapunov指數影響較小，波形變化不大，接近于零，同樣可以判斷動態壓力在做近似的周期振蕩。因此，系統在不同的工作狀態下最大Lyapunov指數的抗干擾性能存在一定差異，在喘振狀態的抗干擾性能好于穩態。但是，采用最大Lyapunov指數仍可以較為清晰的分析出系統的穩定和失穩狀態，尤其是喘振狀態的識別。

4 結 論

本文使用最大Lyapunov指數研究了離心壓縮機出口動壓力的混沌特征。結果表明，穩定狀態下離心壓縮機出口動態壓力的波動量小，頻譜圖分布較分散，具有較大的最大Lyapunov指數，混沌特征明顯，系統的復雜度較高；當系統進入喘振狀態后，壓力序列表現出周期性，波動量變大，規律性變強，具有較小的最大Lyapunov指數，系統的復雜度降低。最大Lyapunov的值作為判斷離心壓縮機發生喘振的一個特征量，可以清晰的反映系統從穩態到喘振的變化過程，可應用于離心壓縮機的防喘振設計。