航空發動機轉子突加不平衡參數分析及LQR控制技術應用

第一作者岳聰男，博士生，1987年7月生

通信作者任興民男，教授，博士生導師，1960年10月生

航空發動機轉子突加不平衡參數分析及LQR控制技術應用

岳聰1,任興民1， 鄧旺群2， 楊永鋒1

(1. 西北工業大學力學與土木建筑學院，西安710072；2.中國航空動力機械研究所，湖南珠洲412002)

摘要：為分析航空發動機轉子系統升速和正常工作時不平衡突變現象，推導了突變時刻動力參數的瞬時變化以及線性二次型調節器(LQR)在變速轉子上的作用過程。通過有限元法和傳遞矩陣法建立動力渦輪轉子瞬態響應分析模型并驗證其有效性。研究了支承剛度變化和突變量變化對動力參數瞬時響應變化的影響以及基于主動控制技術下突加不平衡瞬時響應的變化過程。研究發現，支承剛度和不平衡量變化對瞬時突變的參數響應具有一定的影響變化規律，LQR能有效降低突變引起的不平衡瞬態響應。

關鍵詞：航空發動機；轉子動力學；突加不平衡；瞬態響應；主動控制

基金項目：國家自然基金(11272257)；航空科學基金(2013ZB08001)；陜西省自然科學基金(2013KJXX-22)

收稿日期：2014-08-04修改稿收到日期：2014-10-11

中圖分類號:TH113.1

文獻標志碼：A

DOI:10.13465/j.cnki.jvs.2015.17.029

Abstract:To analyze aeroengine varying-speed rotor system with sudden mass change in accelerating process and work condition, the instantaneous changing of dynamic parameters and the action process of a linear quadratic regulator(LQR) on the rotor system were derived. Both the finite element method and the transfer matrix method were used to establish the analysis model of the power turbine rotor and its effectiveness was verified. The influences of bearing stiffness, sudden mass change and LQR controller on the transient variation of dynamic parameters were studied. This study showed that bearing stiffness and sudden unbalance change have certain effects on the transient change of dynamic parameters and LQR can effectively reduce the transient response of the rotor system caused by sudden unbalance change.

Parametric analysis due to sudden unbalance change and LQR control application for aeroengine rotor systems

YUECong1,RENXing-min1,DENGWang-qun2,YANGYong-fen1(1. School of Mechanics, Civil Engineering & Architecture, Northwestern Polytechnic University, Xi’an 710072, China;2. China Aviation Dynamical Machinery Research Institute, Zhuzhou 412002, China)

Key words:aeroengine system; rotor dynamics; sudden change unbalance; transient response; active control

轉子是航空發動機的核心部件，由于其工作環境惡劣，受力復雜，在運行過程中受到突加不平衡易發生較大的瞬態振動。據不完全統計，某航空發動機僅因這種故障而造成提前更換發動機的比率占17.3%[1]。某航空發動機轉子系統抵抗突加不平衡響應以及轉子系統是否能再次回復穩定運轉的能力，在航空發動機的實際使用過程中有著十分重要的現實功能。

實際渦輪轉子升速過程是比較短暫的，更多的不平衡量發生突變都是發生在工作轉速時刻，目前對應的高速突加不平衡響應研究還相對較少[2-3]。一些研究學者通常采用階躍脈沖函數或連續函數兩種形式來模擬突加不平衡發生過程時不平衡量的突變過程，并用不同的仿真模型計算了穩定轉速下的突加不平衡響應過程。任興民等[4]研究發現在突加不平衡瞬態響應計算中階躍響應很大，推導了與突加不平衡響應相適應的傳遞矩陣。Braut等[5]采用Jeffcott轉子模型和四階Runge-Kutta法對緊急關機后轉子不平衡突然增加進行了線性和非線性響應研究。

與此同時，針對突加不平衡響應控制技術的研究也有了一定進展。姚國志等[6]設計了用于轉子振動控制用的摩擦板式電流變阻尼器能有效地抑制突加不平衡響應。夏南等[7]分析了雷諾數和系統參數對于轉子——擠壓油膜器系統突加不平衡響應和加速響應特性的影響。Su等[8]以產生突加不平衡的原理為基礎，建立了突加不平衡生成裝置模擬試驗臺并進行了主動平衡試驗，將突加不平衡造成的振幅減小了89.6%。Gunter等[9]采用傳遞矩陣法和改進的有限元程序DyRiBeS對1150MW汽輪機發電機轉子系統進行動態分析，開展了線性和非線性轉子系統突加不平衡響應的研究。

本文根據某渦軸發動機動力渦輪轉子相關數據，在分析瞬時突加不平衡瞬時響應參數獲取和驗證所建立動力渦輪轉子模型準確性的基礎上，研究了不同突變量大小和不同突變方位角對參數瞬時不平衡響應的影響以及LQR主動控制技術對不平衡響應突變的控制效果。

1理論分析

1.1不平衡突變瞬時參數

為了說明動力渦輪轉子系統渦輪盤突變處前后位置的瞬時變化，圖1給出了在某時刻特征圓盤上發生突加不平衡前后穩定旋轉位置的變化示意圖；圖2是以形心重合為標準，將突變前后瞬時穩態位置的形心疊加，觀察突變前后物理關系。

圖1　不平衡位置突變示意圖 Fig.1 Sudden change mass on unbalanced position

圖2　不平衡突變前后位置疊加 Fig.2 The superposition of unbalance position

在圖1～圖2中，XOY是固定坐標系，ξGη是慣性坐標系，G是不平衡截面處的形心，c是質心位置，c′是突變后不平衡量的質心，Gc表示初始不平衡偏心距。GΔc是突變不平衡量偏心距，Gc′為突變后的偏心距。φ(θ)表示初始不平衡量與突變后不平衡量的夾角，θ(t)為瞬時相位角，ψ(t)為瞬時進動角。不平衡突變時轉子系統的瞬時位置和不平衡量發生了改變，表現為ψ(t)和θ(t)的瞬時改變。突變先后各變量的對應關系為：

利用上述關系，得到撓度的瞬時變化值

(1)

(2)

(3)

(4)

通過計算式(4)，就可以得到突變后渦輪盤動力參數隨時間變化的表達式。

1.2LQR控制技術

對于某具有4n個自由度的轉子系統，運動微分方程可以描述為如下形式[11]：

(5)

(6)

式(6)是8n維的狀態方程，其中

式中：Bs為系數矩陣，一般為單位矩陣。u(t)為控制力作用矩陣，其求解的二次型能泛函數[12-13]可表示為

(7)

一般情況下，假定加權矩陣Q(t)，R(t)為定常矩陣，記Q(t)=Q，R(t)=R。上式所求的最優控制就是得到J最小時的控制矩陣u(t)。為了達到最優控制目的，構造如下Hamilton函數：

λ[AZ(t)+Bu(t)+DF(t)]

通過求導方法即可求得最優控制信號u(t)為：

u(t)=-R-1BTP(t)Z(t)

其中：P(t)矩陣即為Riccati方程的解：

P(t)BR-1BTP(t)-Q

(8)

上式一般是多個互相耦合的非線性方程組，除特殊情況外是不存在解析解。通常令tf→∞來求解u(t)的穩態解。于是有：

PA+ATP-PBR-1BTP+Q=0

上式稱作代數Riccati方程，求解形式為

ATXE+ETXA-

(ETXB+S)R-1(BTXE+ST-Q)=0

(9)

式(9)中矩陣X就是代數Riccati方程的解。R為殘留矩陣的Frobenius范數。根據上式的解就可以求得最優控制信號u(t)，將u(t)分解為阻尼矩陣和剛度矩陣分別賦予到式(5)中，即可求得LQR作用下系統響應曲線。

2動力渦輪轉子臨界轉速計算

本文分別通過旋轉動力學專用的有限元分析軟件SAMCEF/ROTOR 和數值仿真算法計算所建立模型的臨界轉速以驗證其準確性。

在SAMCEF/ROTOR平臺下，按照圖紙模型以及相關參數，使用梁單元模型對渦軸發動機動力渦輪轉子結構建立有限元模型。考慮到轉子的復雜性，有如下簡化：①動力渦輪轉子共有4個支承，中間兩支承分別與擠壓油膜阻尼器并聯工作;②由于測扭軸質量較小，且不傳遞動力，不考慮其轉動慣量，將其分作六個集中質量附加到渦輪軸上；③分支處的花鍵當做剛性連接處理。有限元模型見圖3，臨界轉速特性坎貝爾圖見圖4。

為了進一步驗證所建立模型的準確性，利用傳遞矩陣法計算得到模型的前兩階臨界轉速與有限元法得到結果的比較如表1所示。結果表明，兩種算法計算得到的前兩階臨界轉速誤差小于5%，基本保持一致；另一方面，本文計算結果與文獻[14-15]所給出的計算值基本吻合，雖然某些部位的剛性不能與發動機完全一致，但可以用來進行實際的不平衡響應分析。

圖3　T700動力渦輪轉子模型圖 Fig.3 Model of T700 turbine power rotor

圖4　坎貝爾圖 Fig.4 Campbell diagram

計算方法一階臨界轉速二階臨界轉速有限元法7421r/min13452r/min傳遞矩陣法7394r/min13817r/min誤差0.3%2.7%

3突變瞬時不平衡響應分析

為了研究高轉速下動力渦輪轉子突加不平衡瞬態響應，本文分別從支承剛度和突變量的變化對突加不平衡響應的影響進行分析。假設兩個渦輪盤上存在不平衡量，分別為3.75 g·cm∠180°和5.376 g·cm∠36°，根據式(5)知n=2。利用Newmark-β積分法[11]，積分步長為1×10-3s，升速過程中加速度a為600 rad/s2，在此基礎上，模擬在3.5 s時發生突加不平衡，此時轉速為動力渦輪轉子的設計工作轉速(20 000 r/min)。下面分別分析不同支承剛度、不同突變質量大小和方向以及主動控制對工作轉速下瞬時不平衡量突變的影響。

3.1支承剛度影響

四個支承的初始剛度系數分別為：k1=7.26×107，k2=4×107，k3=1×107，k4=1.42×107(單位：N/m)。本節討論了四個支承處剛度系數的變化對工作轉速時突加不平衡產生的瞬時動撓度r突變峰值的影響。在允許的變化范圍內，以支承3剛度系數為例，其他支承剛度以表2所示的組合，觀察隨支承剛度變化下，瞬時不平衡量突變引起Δr變化。兩級渦輪盤上的瞬時響應見圖5～圖6。

圖5　不同剛度下一級渦輪盤突加不平衡瞬時撓度響應 Fig.5 Transient response of deflection under different stiffness combination with sudden change mass on turbine discⅠ

圖6　二級渦輪盤突加不平衡動撓度瞬時響應 Fig.6 Transient response of deflection under different stiffness combination with sudden change mass on turbine discⅡ

由于一、四組只有支承2剛度的不同，可以看到支承3剛度在增大時對渦輪盤突變時刻的動撓度響應影響不是很大，對應圖5～圖6中的響應基本穩定；而二、三組與一、四組相比，突變時刻峰值具有一定的波動周期性。圖中的響應曲線表明，隨著支承3剛度變化，突變時刻瞬時撓度的峰值變化趨勢是基本一致的，而在選取合適的剛度后，瞬時突變不平衡響應帶來的幅值變化較小。

表2　不同支承剛度組合(單位：N/m)

3.2突變量大小變化影響

由于渦軸發動機上兩個渦輪盤之間的距離較近，僅相隔0.06 m，故兩個渦輪盤上的瞬態響應整體變化趨勢總體來說基本一致，只是在數值上有些區別，本文僅以一級渦輪盤上的響應為例進行分析。采用3.1節相同的假設，令二級渦輪盤突變后的不平衡量變為8.96 g·cm∠36°；同一時刻，在一級渦輪盤上施加不平衡量，使不平衡偏心距在0～8×10-5m范圍內變化，取值間隔為5×10-7m，不平衡方位角0～2π內，取值間隔為π/8。仿真計算得到的一級渦輪盤上瞬時動撓度的變化情況分別見圖7。

圖7中，隨著一級渦輪盤上不平衡量在選定范圍內不同，突變時刻的瞬時不平衡響應的幅值大小也會隨之發生改變，并且對應一個二級渦輪盤上的突變不平衡量，總有個對應的一級渦輪盤上的不平衡量能將突變帶來影響降到最低。

3.3突變瞬時對不平衡響應影響

為了反映不平衡突變瞬時響應的識別，圖8中給出了一級渦輪盤在3.2節條件下與無突變下不平衡響應對應時刻瞬時動撓度之差的響應圖。從圖中可以看到，一級渦輪盤取不平衡方位角為5π/8和13π/8時，突變前后兩盤上的瞬態響應差接近0，當不平衡方位角在0和π附近時，突變前后瞬態響應之差最大。

圖7　一級渦輪盤瞬時動撓度變化 Fig.7 Transient deflection change on turbine discⅠ

圖8　一級渦輪盤不平衡突變前后瞬態響應差 Fig.8 Difference between the amplitude response before and after sudden change unbalance on turbineⅠ

3.4LQR控制應用

為了驗證主動控制技術對不同支承剛度下突加不平衡響應的抑制作用，現將LQR主動控制方法應用到不平衡響應突變中去。根據渦輪盤的質量和維數，控制權矩陣Q取值與M相近，R則與E相關：

Q=10×diag(6161616161616161)

R=0.1×E8×8

通過式(7～9)得到每個積分時間下的最優控制信號u(t)后求得瞬時動力渦輪轉子系統響應。圖9給出了動力渦輪轉子加速起動越過前兩階臨界轉速，在1 s發生突加不平衡響應時一級渦輪盤的瞬態響應。從對比的突變響應中可以看到，隨著Q、R取合適的系數后，不僅對轉子系統越過臨界轉速的峰值有效控制，同時對突變時刻的瞬時變化有一定的控制作用。

圖10顯示了表2中第三組所示各支承剛度系數下兩級渦輪盤上有主動控制作用下瞬時響應與無控制下突變響應結果的對比。從圖10中可以看到，在選擇適當Q、R矩陣的參數后，系統對于不同支承剛度下發生突變的瞬時峰值明顯得到了抑制，突變幅值波動幅度也控制在一個較小的范圍內。

圖9　 一級渦輪盤升速過程控制前后 突加不平衡瞬態響應 Fig.9 Unbalance transient response comparison with active control on disc1 in accelerating process

圖10　控制前后突變峰值對比 Fig.10 Instantaneous deflection response with LQR controller

4結論

本文給出了轉子不平衡量突變瞬時的動力特性參數分析和LQR控制的應用方法，同時利用有限元法和數值仿真兩種算法建立了可用于突加不平衡響應分析的渦軸發動機動力渦輪轉子模型。研究發現：

(1)不同組合的支承剛度對不平衡量突變的響應有所差別，當設計有恰當的支承剛度時，能有效降低突加不平衡引起的幅值變化；

(2)隨著不平衡位置處出現不平衡量突變后，其他位置處的不平衡量大小和方位角均會對突變的瞬時響應產生影響。當達到某個值時，可以消除渦輪盤上不平衡突變量帶來的振動影響；

(3)加入LQR主動控制作用后，動力渦輪轉子出現突加不平衡現象時，可以將轉子系統的參數變化幅值控制在在較小范圍內，降低了不平衡突變振動對渦軸發動機動力渦輪轉子正常工作的影響。

參考文獻

[1]唐洪文，陸永忠，廖道訓. 擠壓油膜頭型轉自系統突加質量的響應特性[J]. 華中科技大學學報，2002，30(11): 95-97.

TANG Hong-wen, LU Yong-zhong, LIAO Dao-xun. Characteristics of imbalance response in the sliding bearing-flexible rotor system[J]. Huazhong University of Science and Technology(Nature Science edition), 2002,30(11): 95-97.

[2]Tondl A, Some problems of rotor dynamics[M]. Chapman and Hall,London,1965.

[3]Das A S, Dutt J K, Ray K. Active control of coupled flexural-torsional vibration in a flexible rotor-bearing system using electromagnetic actuator[J]. International Journal of Non-Linear Mechanics, 2011, 46(9):1093-1109.

[4]任興民，顧家柳.航空發動機轉子-支撐系統的突加不平衡響應[J]. 振動工程學報,1991,4(3):57-82.

REN Xing-min, GU Jia-liu.Response of suddenly applied unbalance foraeroengine rotor-bearing systems[J]. Journal of Vibration Engineering, 1991,4(3): 57-82.

[5]Braut S J, Zigulic R, Butkovic M, Numerical and experi-mental analysis of a shaft bow influence on a rotor to stator contact dynamics[J]. Journal of Mechanical Engineering, 2008, 54(10): 693-706.

[6]姚國志，邱陽，方同，等. 多層擠壓型電流變阻尼器在轉子振動控制中的實驗研究[J].航空學報, 1998, 19(5): 570-574.

YAO Guo-zhi,QIU Yang,FANG Tong, et al. Experimental study of a multiple squeeze films electrorh-eological damper for vibration control of rotor system[J]. Acta Aeronautica et Astronautica Sinica, 1998,19(5):570-574.

[7]夏南，孟光，馮心海. 油膜慣性力對雙盤栓子SFD系統突加不平衡和加速響應特性的影響[J]. 航空動力學報, 2000,15(1):71-74.

XIA Nan, MENG Guang, FENG Xing-hai. Inf iuence of fluid inertia force on sudden unbalance and accelerating responses of rotor supported on centralized SFD[J]. Journal of Aerospace Power, 2000, 15(1): 71-74.

[8]Su Y R, He L D. Study on online elimination of sudden unbalance-induced vibration using active balancing technology[J]. High Technology Letters, 2012, 16(2): 210-214.

[9]Gunter E J,Chen W J. Dynamic analysis of an 1150MW turbine-generator[J]. Proceedings of the ASME Power Conference,Chicago,2005,Part A:437-443.

[10]顧家柳等編著. 轉子動力學[M]. 北京：國防工業出版社, 1988: 46-58.

[11]岳聰,任興民,鄧旺群,等. 基于升速響應信息柔性轉子系統的多階多平面瞬態動平衡方法[J]. 航空動力學報, 2013, 28(11): 2953-2599.

YUE Cong, REN Xing-min, DENG Wang-qun, et al. Multi-plane and multi-critical transient dynamic balance method based on rising speed response information of flexible rotor system[J]. Journal of Aerospace Power, 2013, 28 (11): 2953-2599.

[12]Saptarshi D, Indranil P, Kaushik H, et al.LQR based improved discrete PID controller design via optimum selection of weighting matrices using fractional order integral performance index[J]. Applied Mathematical Modeling, 2013, 6(15): 4253-4268.

[13]Zhang G M, Ou L L, Meng Li, et al. LQR control for homogeneous agents with multi-graph topology[J]. Acta Automatica Sinica,2013,39(6): 913-919.

[14]Burgess G, Rio R. T700 power turbine rotor multiplane/multispeed balancing demonstration[R]. MTI TechnicalReport, NASA Contractor Report CR-159586, Prepared for NASA-Lewis Research Center under Contract No. NAS3-18529, 1979.

[15]Wahon J, Lee C, Martm M. High speed balancing applied to the T700 engine[R]. MTI Technical Report No. MTI87TR56, NASA Contractor Report CR-180899, Prepared forNASA-Lewis Research Center under Contract No. NAS3-23929 and NAS3-24633, March 1989.