汽輪機(jī)轉(zhuǎn)子系統(tǒng)中共振組織的全局拓?fù)湟?guī)律分析

徐 璐， 褚衍東， 楊 瓊， 李險(xiǎn)峰

(1.蘭州交通大學(xué) 機(jī)電工程學(xué)院，蘭州 730070； 2.蘭州交通大學(xué) 數(shù)理學(xué)院，蘭州 730070)

《世界能源展望》中提到可再生能源的全球份額在發(fā)電領(lǐng)域迅猛增長,現(xiàn)有常規(guī)發(fā)電廠的目標(biāo)之一是補(bǔ)償剩余負(fù)荷，即使在將來的能源情景中，汽輪機(jī)依舊起著重要作用。因此，獲取更多汽輪機(jī)轉(zhuǎn)子系統(tǒng)動(dòng)力學(xué)信息，對(duì)其機(jī)理進(jìn)一步地認(rèn)識(shí)，最終進(jìn)行非線性動(dòng)力學(xué)設(shè)計(jì)，具有重要的現(xiàn)實(shí)意義和廣泛的工程背景。

研究大型轉(zhuǎn)子系統(tǒng)的穩(wěn)定性與分岔理論具有重要意義。描述系統(tǒng)動(dòng)力學(xué)問題的微分方程具有非線性特性，在一定條件下必然發(fā)生分岔和混沌現(xiàn)象[1-5]。劉小峰等[6]利用分岔圖、Poincaré截面圖等，揭示了交叉剛度對(duì)碰摩轉(zhuǎn)子系統(tǒng)運(yùn)行狀態(tài)的影響。徐文標(biāo)等[7]基于Paris裂紋擴(kuò)展規(guī)律，在參數(shù)平面內(nèi)分析了裂紋擴(kuò)展對(duì)含拉桿裂紋的轉(zhuǎn)子系統(tǒng)動(dòng)力學(xué)特性的影響，發(fā)現(xiàn)了多種分岔現(xiàn)象。張金劍等[8]利用單參數(shù)數(shù)值工具探討了質(zhì)量偏心、勵(lì)磁電流和氣隙動(dòng)靜偏心對(duì)水電機(jī)組碰摩轉(zhuǎn)子-軸承系統(tǒng)的影響。我國對(duì)于汽流激振的研究起步較晚，是從1980年開始的。柴山等[9-11]對(duì)汽輪機(jī)的汽流激振力做了較詳細(xì)的研究，包括葉片優(yōu)化設(shè)計(jì)、汽流激振力計(jì)算公式的推導(dǎo)及計(jì)算程序的設(shè)計(jì)等。接著，文獻(xiàn)[12-17]充分利用分岔圖、Poincaré截面圖、軸心軌跡圖、頻譜圖等單參數(shù)數(shù)值工具討論了高壓或低壓汽輪機(jī)轉(zhuǎn)子系統(tǒng)的分岔和混沌行為。最近，曹麗華等[18-19]以耦合了非線性汽流激振力和油膜力的超超臨界汽輪機(jī)轉(zhuǎn)子系統(tǒng)為研究對(duì)象，基于分岔圖、最大Lyapunov指數(shù)圖等，揭示了非線性汽流激振力對(duì)系統(tǒng)的影響。翁雷等[20]對(duì)汽輪機(jī)轉(zhuǎn)子系統(tǒng)做了較全面的分析，分別探討了葉尖汽流激振對(duì)汽輪機(jī)轉(zhuǎn)子、基礎(chǔ)松動(dòng)轉(zhuǎn)子[21]、裂紋轉(zhuǎn)子[22]、碰摩轉(zhuǎn)子[23]及裂紋-碰摩[24]轉(zhuǎn)子系統(tǒng)的影響。羅躍剛等[25]分析了汽流激振和偏心量對(duì)密封-轉(zhuǎn)子系統(tǒng)的影響。

綜上所述，盡管關(guān)于轉(zhuǎn)子系統(tǒng)非線性動(dòng)力學(xué)行為的研究有豐碩的成果，但是都是通過單參數(shù)分岔分析實(shí)現(xiàn)，導(dǎo)致動(dòng)力學(xué)信息有限且不完整。這表明我們對(duì)大型旋轉(zhuǎn)機(jī)械的雙參數(shù)結(jié)構(gòu)的理解還很遠(yuǎn)，更不用說更高維的參數(shù)空間。因此，針對(duì)大型旋轉(zhuǎn)機(jī)械系統(tǒng)的非線性動(dòng)力學(xué)行為的研究，真正的挑戰(zhàn)是發(fā)現(xiàn)并總結(jié)二維參數(shù)平面內(nèi)分岔行為的拓?fù)淦毡樾浴?/p>

基于以上原因，本文以油膜支撐，在汽流激振力、非線性油膜力和不平衡離心力共同作用下的汽輪機(jī)轉(zhuǎn)子系統(tǒng)為研究對(duì)象，利用高清的雙參數(shù)平面周期穩(wěn)定相圖、分岔圖、Lyapunov指數(shù)圖和時(shí)間序列圖，分析了轉(zhuǎn)子系統(tǒng)內(nèi)周期吸引子的全局拓?fù)湟?guī)律，目的是補(bǔ)充和完善已得到的研究成果，更全面地掌握汽輪機(jī)轉(zhuǎn)子系統(tǒng)的非線性動(dòng)力學(xué)特性，揭示更多隱藏在系統(tǒng)中的動(dòng)力學(xué)信息，為大型旋轉(zhuǎn)機(jī)械實(shí)際設(shè)計(jì)及故障預(yù)測提供有價(jià)值的依據(jù)。

1 系統(tǒng)力學(xué)模型及運(yùn)動(dòng)方程

1.1 非線性間隙汽流激振力計(jì)算公式

圖1為系統(tǒng)受力簡圖，其中fax，fay是汽輪機(jī)非線性汽流激振力Fa在x和y方向的分力。圖2為沿環(huán)向的兩個(gè)橫剖面展開圖。

圖1 系統(tǒng)受力圖Fig.1 Loading diagram of the system

圖2 汽流在靜、動(dòng)葉片間流動(dòng)示意圖Fig.2 The schematic diagram of air-exciting force moving between stator and rotor blades

則非線性汽流激振力的無量綱模型為：

Fa=A1δE+A3δ3E3

(1)

其中

1.2 油膜力計(jì)算公式

忽略扭轉(zhuǎn)振動(dòng)和陀螺力矩，在圖3中：O1,O2分別為軸承內(nèi)瓦和轉(zhuǎn)子幾何中心；O3為轉(zhuǎn)子質(zhì)心；轉(zhuǎn)子兩端由半徑為R，長度為L的滑動(dòng)軸承支撐；m1,c1和m2,c2分別為轉(zhuǎn)子在軸承處和圓盤處的等效集中質(zhì)量和結(jié)構(gòu)阻尼；Fx,Fy為非線性油膜分量；h為平均油膜厚度。

圖3 油膜支撐轉(zhuǎn)子系統(tǒng)示意圖Fig.3 A schematic diagram of the rotor system supported by oil film

圖3是油膜支撐轉(zhuǎn)子示意圖，利用短軸承模型[26]，則設(shè)x，y為軸承位移，則無量綱油膜力在x和y方向的分力為

(2)

其中

(3)

(4)

(5)

(6)

1.3 系統(tǒng)動(dòng)力學(xué)方程

設(shè)x1，y1為轉(zhuǎn)子左端軸承處的徑向位移；x2，y2為轉(zhuǎn)盤處的徑向位移，則系統(tǒng)的無量綱運(yùn)動(dòng)方程為

(7)

式中，τ=ωt。式(7)中的參數(shù)意義及取值，如表1所示。

表1 汽輪機(jī)轉(zhuǎn)子系統(tǒng)參數(shù)意義及取值Tab.1 The parameter meanings and values of the system

2 計(jì)算方法

在本文，雙參數(shù)平面內(nèi)的數(shù)值結(jié)果主要使用isoperiodic圖[27]呈現(xiàn)。Isoperiodic圖不僅顯示了一個(gè)周期內(nèi)尖峰數(shù)量的積累，即波形的復(fù)雜性，以及在參數(shù)平面內(nèi)尖峰數(shù)量積累的地方以及積累的快慢，而且也可以有效地用于處理試驗(yàn)數(shù)據(jù)。對(duì)于大型旋轉(zhuǎn)機(jī)械系統(tǒng)，制作isoperiodic圖是一項(xiàng)非常耗時(shí)的任務(wù)。因此，采用GPU并行計(jì)算技術(shù)。

為了生成isoperiodic穩(wěn)定相圖，借助GPU并行計(jì)算的優(yōu)勢，將感興趣的參數(shù)窗口劃分為N×N的等距點(diǎn)網(wǎng)格，通?？蛇_(dá)2 000×2 000=4×106。對(duì)于每個(gè)點(diǎn)，通過對(duì)方程進(jìn)行數(shù)值積分來獲得時(shí)間演化。對(duì)方程組式(7)使用標(biāo)準(zhǔn)四階Runge-Kutta算法。為了進(jìn)行積分，從任意選擇的初始條件開始，從左到右水平掃描參數(shù)，通過“跟隨吸引子”向右進(jìn)行。在所有情況下，都會(huì)放棄起初2×105步以避免瞬態(tài)響應(yīng)。在對(duì)時(shí)間步長進(jìn)行后續(xù)80×105積分后，可以獲得isoperiodic穩(wěn)定相圖。每一張isoperiodic圖會(huì)使用19種顏色的調(diào)色板顯示每個(gè)周期的尖峰數(shù)量。通過循環(huán)19個(gè)基本顏色“模19”來繪制超過19個(gè)尖峰的振蕩，即通過為它們分配一個(gè)顏色索引。黑色用于表示“混沌或擬周期”，即缺乏數(shù)字可檢測的周期性。

同時(shí)，本文采用一種新的方法，即根據(jù)運(yùn)動(dòng)方程初始小的擾動(dòng)副本(“克隆”)的時(shí)間演化，計(jì)算模型的李雅普諾夫指數(shù)譜，其中雅可比矩陣的復(fù)雜計(jì)算替換為線性化版本的運(yùn)動(dòng)方程的時(shí)間演化?！翱寺》ā崩钛牌罩Z夫指數(shù)的簡要并行計(jì)算代碼請(qǐng)參考文獻(xiàn)[28-29]。

3 具體數(shù)值結(jié)果分析

在圖4中偏心量e和旋轉(zhuǎn)速度ω同時(shí)變化，進(jìn)汽速度v保持200不變，剛度k=2.5×107。通過右側(cè)的顏色棒很容易區(qū)分周期窗口的周期數(shù)、大小及形狀。將圖4(a)的周期序列劃分為兩組，這兩組周期結(jié)構(gòu)在圖4(b)和圖4(c)中進(jìn)一步放大詳述。這些周期區(qū)域具有自相似性。

圖4 汽輪機(jī)轉(zhuǎn)子系統(tǒng)在(ω,e)參數(shù)平面內(nèi)的穩(wěn)定周期相圖Fig.4 The isoperiodic diagrams of the steam turbine rotor system in (ω,e) parameter plane

在圖4(b)中，13條周期數(shù)不同的“細(xì)絲”嵌入在大片的黑色區(qū)域內(nèi)。為了更加深入地理解這組周期序列的演化細(xì)節(jié)，沿著藍(lán)色直線做單參數(shù)分岔圖。在圖5(a)中，水平掃描e，隨著e的增大且旋轉(zhuǎn)速度ω固定在2 000，13個(gè)周期窗口依次鎖定在非周期區(qū)域內(nèi)：P-3→P-17→P-14→P-11→P-8→P-13→P-18→P-5→P-17→P-12→P-7→P-16→P-9。圖5(b)呈現(xiàn)了通過“克隆”法求得的最大Lyapunov指數(shù)圖，可見堆積在其中的周期區(qū)域?qū)M周期區(qū)域瓜分，這是一個(gè)周期解與擬周期解頻繁轉(zhuǎn)換的過程。擬周期解隨參數(shù)變化在某些參數(shù)域鎖相到周期解，因此這些周期區(qū)域稱為共振組織或鎖模結(jié)構(gòu)。

圖5(a)中A點(diǎn)處(e=0.060 8×10-3)的狀態(tài)傳遞矩陣的特征值為：λ=[-0.760 7±0.131 8i,0.359 6±0.658 1i, 1.012 1, 0.356 8]，主導(dǎo)Floquet乘子|λ|max=1.012 1>1，即|λ|max通過正實(shí)數(shù)穿過單位圓，因此在A點(diǎn)處發(fā)生了鞍結(jié)分岔(SN)。在B點(diǎn)處(e=0.070 5×10-3)的狀態(tài)傳遞矩陣的特征值為：λ=[-0.542 7±0.086 0i, 0.446 1±0.702 0i, 1.052 7, 0.211 8]，主導(dǎo)Floquet乘子|λ|max=1.052 7>1，即|λ|max通過正實(shí)數(shù)穿過單位圓，則在B點(diǎn)處發(fā)生了鞍結(jié)分岔(SN)。

由此可得：周期解通過鞍結(jié)分岔進(jìn)入擬周期解，即兩根“細(xì)絲”之間的邊界為鞍結(jié)分岔曲線。

旋轉(zhuǎn)數(shù)不斷變化，結(jié)果形成了“魔鬼樓梯”[32](如圖5(c))。隨著e的不斷增加，r/s不是線性增長，而是一系列的跳躍和不同大小的臺(tái)階，不同的臺(tái)階代表著不同的鎖模頻率，鎖模的大小則由臺(tái)階的寬度來反映。同時(shí)，這些臺(tái)階具有分形的自相似性。沿著臺(tái)階，整個(gè)“魔鬼樓梯”平穩(wěn)地單調(diào)遞增，表明周期運(yùn)動(dòng)被鎖定在關(guān)于e的一個(gè)不斷增加且有限范圍內(nèi)的某個(gè)方向上。

圖5 沿著圖4(b)中藍(lán)色直線的分岔圖、指數(shù)圖及“魔鬼樓梯”Fig.5 The bifurcation diagram,Lyapunov exponent diagram and “Devil’s staircase” along the blue line of Fig.4(b)

圖4(c)放大了第二組共振組織。正如圖6所示，有7個(gè)周期窗口嵌入在擬周期區(qū)域內(nèi)：P-11→P-18→P-7→P-17→P-10→P-13→P-16。通過計(jì)算在C點(diǎn)處(ω=1 812)的狀態(tài)傳遞矩陣的特征值：λ=[0.530 3±0.721 0i,-0.382 4±0.508 4i,0.130 5±0.495 2i]，由3對(duì)共軛復(fù)數(shù)組成，主導(dǎo)Floquet乘子|λ|max=1.002 5>1，即|λ|max通過其中一對(duì)共軛復(fù)數(shù)穿過單位圓，因此在C點(diǎn)處發(fā)生了Hopf分岔。由于Hopf分岔的存在導(dǎo)致了共振組織的出現(xiàn)。

圖6 沿著圖4(c)中白色直線的分岔圖及“魔鬼樓梯”Fig.6 The bifurcation diagram and “Devil’s staircase” along the white line of Fig4.(c)

這組共振組織的旋轉(zhuǎn)數(shù)的分布也滿足法里樹序列的拓?fù)湟?guī)律，且形成了單調(diào)遞增的“魔鬼樓梯”。

圖8(c)～圖(i)顯示了圖7(a)中被淡藍(lán)色的A,B,...,G點(diǎn)所標(biāo)記的參數(shù)組合下的強(qiáng)度振蕩的時(shí)間演化過程，并通過波形變形說明了新尖峰的產(chǎn)生，揭示了一種按有規(guī)律的方式所組織的周期振蕩。綠色垂直線標(biāo)記振蕩周期T1。圖8(c)～圖8(i)中的時(shí)間演化表明了在轉(zhuǎn)子模式的穩(wěn)定復(fù)雜性中隱藏的規(guī)律：每一族似乎都是一些準(zhǔn)相同模式的固定組合的串聯(lián)，尖峰越多，波形就越多。

圖7 k×v參數(shù)平面內(nèi)穩(wěn)定振蕩和混沌振蕩的曲折分布Fig.7 The distribution of periodic and chaotic oscillations in k×v parameter plane

圖8 沿著圖7(a)中黃色直線的分岔圖、指數(shù)圖及時(shí)間歷程圖Fig.8 The bifurcation diagram,Lyapunov exponent diagram and time series diagrams along the yellow line of Fig7.(a)

在圖7(b)中，隨著v的增大，順著水平方向掃描圖像，有8個(gè)周期窗口堆積在黑色的擬周期區(qū)域內(nèi)。通過僅僅調(diào)整一個(gè)參數(shù)而得到的單參數(shù)分岔圖9(a)，相應(yīng)地，圖9(b)中有8個(gè)負(fù)指數(shù)被鎖在虛直線以下及圖9(c)中單調(diào)遞減的“魔鬼樓梯”。

圖9 沿著圖7(b)中藍(lán)色直線的分岔圖、指數(shù)圖及“魔鬼樓梯”Fig.9 The bifurcation diagram,Lyapunov exponent diagram and “Devil’s staircase” along the blue line of Fig7.(b)

圖10呈現(xiàn)了參數(shù)空間k×ω，k×e，v×e內(nèi)共振組織的全局拓?fù)湟?guī)律。很容易識(shí)別出這些不同的雙參數(shù)平面具有相同的拓?fù)湟?guī)律，即法里樹序列。因此，滿足法里和的共振組織分布是該系統(tǒng)的共性。在旋轉(zhuǎn)機(jī)械中，轉(zhuǎn)子的質(zhì)量偏心會(huì)產(chǎn)生振動(dòng)頻率為ω的周期性激振力。若轉(zhuǎn)子彈性力有非線性就會(huì)出現(xiàn)許多非線性共振現(xiàn)象。

圖10 k×ω，k×e，v×e平面內(nèi)共振組織的全局拓?fù)錂C(jī)制Fig.10 The global topological mechanism of periodic oscillations of k×ω，k×e，v×e parameter planes

4 結(jié) 論

針對(duì)大型旋轉(zhuǎn)機(jī)械系統(tǒng)，對(duì)其中的動(dòng)力學(xué)現(xiàn)象做詳細(xì)的分類是一個(gè)費(fèi)時(shí)費(fèi)力的工作，由于GPU并行計(jì)算方法的卓越性能和算法的不斷改進(jìn)，本文從二維參數(shù)角度揭示了汽輪機(jī)轉(zhuǎn)子系統(tǒng)中共振組織排列的拓?fù)湟?guī)律，可以得到：滿足法里和序列的共振組織分布是該系統(tǒng)的共性。從非線性動(dòng)力學(xué)角度來看，二維參數(shù)平面相圖顯示了波形的復(fù)雜性，以及周期吸引子在控制參數(shù)空間中積累的速度和位置。因此，這一共性有利于預(yù)測在更寬的參數(shù)空間動(dòng)力學(xué)行為的演變；從參數(shù)匹配的角度來看，周期振蕩的展開會(huì)表明系統(tǒng)對(duì)某些參數(shù)的敏感程度，因此實(shí)際試驗(yàn)也可以從以上的圖表中受益，并提供需要什么樣的精度才能區(qū)分相圖中呈現(xiàn)的周期和混沌窗口的復(fù)雜交替。從實(shí)際應(yīng)用角度來看，拓?fù)湟?guī)律的提取是將復(fù)雜的機(jī)械現(xiàn)象轉(zhuǎn)化為簡單的數(shù)學(xué)問題，為優(yōu)化設(shè)計(jì)、參數(shù)匹配及故障預(yù)測提供了可視化、可比較的圖表，有利于對(duì)轉(zhuǎn)子系統(tǒng)非線性動(dòng)力學(xué)特性更進(jìn)一步的理解。