趙合陽(yáng) 白廣忱 費(fèi)成巍 邢偉震
(北京航空航天大學(xué) 能源與動(dòng)力工程學(xué)院,北京100191)(中國(guó)人民解放軍空軍93705部隊(duì),遵化064200)
葉尖間隙直接影響著航空發(fā)動(dòng)機(jī)的效率、可靠性、安全性和耗油率[1].其合理設(shè)計(jì)和控制一直是研制高性能、高可靠性航空發(fā)動(dòng)機(jī)亟待解決問(wèn)題之一[2].渦輪盤徑向變形對(duì)葉尖徑向間隙有重要影響,盤的合理設(shè)計(jì)分析對(duì)葉尖間隙設(shè)計(jì)和控制具有重大意義.發(fā)動(dòng)機(jī)一個(gè)工作循環(huán)中,輪盤徑向變形受機(jī)械載荷和熱載荷等方面諸多因素影響,隨工作狀態(tài)的不同而變化,合理選擇徑向變形的分析方法是進(jìn)行渦輪盤設(shè)計(jì)和葉尖徑向運(yùn)行間隙設(shè)計(jì)與控制的基礎(chǔ)[3-4].
目前國(guó)內(nèi)外許多研究機(jī)構(gòu)和學(xué)者對(duì)渦輪盤徑向變形進(jìn)行了數(shù)值分析[1-8].這些研究成果都基于確定性分析,即采用確定的參數(shù)來(lái)確定發(fā)動(dòng)機(jī)工作時(shí)的渦輪盤徑向變形.由于該方法沒(méi)有考慮影響渦輪盤等構(gòu)件變形的各方面因素的隨機(jī)性,只靠留出裕度的方法來(lái)保證葉尖間隙的安全性,因而具有很大的盲目性.事實(shí)上,影響渦輪盤變形的諸多因素都具有很明顯的隨機(jī)性,應(yīng)該考慮各種隨機(jī)因素,從概率的角度研究非線性動(dòng)態(tài)概率分析,才能更客觀描述渦輪盤變形規(guī)律,改善葉尖徑向間隙設(shè)計(jì)和控制的合理性.
因此,基于有限元和響應(yīng)面方法[9-11],考慮材料屬性的非線性、熱載荷和離心載荷的動(dòng)態(tài)性,嘗試對(duì)輪盤徑向變形進(jìn)行非線性動(dòng)態(tài)概率分析.
選取某型航空發(fā)動(dòng)機(jī)一級(jí)高壓渦輪盤為研究對(duì)象,輪緣前端與渦輪軸以圓柱表面定心,用徑向銷釘鏈接,工作時(shí)前后均有空氣對(duì)輪盤側(cè)面進(jìn)行冷卻.渦輪盤的幾何形狀是軸對(duì)稱的,可認(rèn)為渦輪盤承受的外載荷和約束情況是軸對(duì)稱的.簡(jiǎn)化榫槽、銷釘孔等結(jié)構(gòu),建立渦輪盤有限元模型,如圖1所示.
圖1 渦輪盤有限元模型
根據(jù)發(fā)動(dòng)機(jī)冷卻系統(tǒng)特點(diǎn)和發(fā)動(dòng)機(jī)實(shí)測(cè)溫度,輪盤外緣A1、輪盤前側(cè)與封嚴(yán)圈環(huán)結(jié)合處A2、輪盤中心A3等處的溫度邊界可作為第1類邊界條件處理,渦輪盤前、后側(cè)面邊界B1,B2,B3可用熱交換的第3類邊界處理.
渦輪盤非線性動(dòng)態(tài)概率分析的基本思想為:
1)選取某航空發(fā)動(dòng)機(jī)典型飛行剖面參數(shù)[2],將發(fā)動(dòng)機(jī)從地面啟動(dòng)—慢車—起飛—爬升—巡航這一段[0,215 s]作為計(jì)算范圍,取12個(gè)關(guān)鍵點(diǎn)作為計(jì)算點(diǎn),其飛行載荷譜如圖2所示.
2)考慮材料屬性的非線性和載荷的動(dòng)態(tài)性,進(jìn)行渦輪盤變形熱-固耦合分析,計(jì)算出渦輪盤徑向變形規(guī)律,并將渦輪盤最大變形量(即全局最大值)作為概率設(shè)計(jì)的輸出響應(yīng).
3)合理選取隨機(jī)變量,對(duì)隨機(jī)參數(shù)進(jìn)行確定性擬合試驗(yàn),再對(duì)渦輪盤徑向變形輸出響應(yīng)進(jìn)行精確求解,獲得足夠多的樣本值,利用這些樣本值擬合輸出響應(yīng)與隨機(jī)變量之間的響應(yīng)面模型.
4)將響應(yīng)面模型代替有限元模型,基于Monte Carlo法(MC法)對(duì)渦輪盤徑向變形進(jìn)行概率分析和靈敏度分析.
圖2 計(jì)算采用的載荷譜
假設(shè)每一次渦輪盤徑向變形的非線性動(dòng)態(tài)分析的輸出響應(yīng) Y 與輸入?yún)?shù) X=[X1,X2,…,Xr]之間的響應(yīng)面模型為
式中,a0,bi,cij(i=1,2,…,r;j=1,2,…,r)為待定系數(shù).
假設(shè)系統(tǒng)要求渦輪盤最大徑向變形量為δ,則極限狀態(tài)函數(shù)為
由式(2)可知:H<0時(shí)為失效模式;反之為安全模式.若各隨機(jī)變量相互獨(dú)立且均值和方差矩陣為
則
其中,E(·)和D(·)分別表示均值和方差函數(shù).
若 H(X)符合正態(tài)分布,則可靠度[11]為
由式(5)可得各個(gè)隨機(jī)參數(shù)的靈敏度:
本文考慮離心載荷和熱載荷的動(dòng)態(tài)性(如圖2所示)和導(dǎo)熱系數(shù)和膨脹系數(shù)的非線性(如表1所示),利用熱-固耦合的分析方法對(duì)渦輪盤徑向變形進(jìn)行非線性動(dòng)態(tài)概率分析.首先,結(jié)合導(dǎo)熱系數(shù)對(duì)其進(jìn)行非線性動(dòng)態(tài)熱載荷分析;然后,將熱分析換為結(jié)構(gòu)分析,將計(jì)算得到的溫度場(chǎng)作為溫度載荷加入結(jié)構(gòu)分析模型中,在熱載荷和機(jī)械載荷的耦合作用下進(jìn)行徑向變形量計(jì)算.經(jīng)計(jì)算得渦輪盤徑向變形隨時(shí)間變化規(guī)律如圖3所示.
表1 不同溫度下的導(dǎo)熱和膨脹系數(shù)(參考溫度20℃)
圖3 輪盤徑向變形量隨時(shí)間變化規(guī)律
由圖3可以看出:在發(fā)動(dòng)機(jī)啟動(dòng)—慢車—起飛爬升過(guò)程中,輪盤的變形量呈增大趨勢(shì);在加速爬升時(shí)(即 t∈[165 s,190 s]時(shí))輪盤的變形量達(dá)到了最大值;由起飛進(jìn)入巡航狀態(tài)時(shí),渦輪盤變形量又開(kāi)始稍微減小.經(jīng)分析,當(dāng)t=178 s時(shí)(燃?xì)鉁囟茸罡?、轉(zhuǎn)速最大的時(shí)刻)渦輪盤變形量達(dá)到全局最大值:1.236 mm,此時(shí)輪盤中各節(jié)點(diǎn)的最大變形量同時(shí)達(dá)到最大,其中渦輪盤溫度分布和徑向變形云圖如圖4和圖5所示.
圖4 t=178 s時(shí)渦輪盤溫度分布云圖
圖5 t=178 s時(shí)輪盤徑向變形量云圖
由圖4可知:渦輪盤的溫度場(chǎng)不但沿徑向變化,也沿軸向變化.輪緣處的溫度最高,沿徑向渦輪盤溫度逐漸降低,在盤心處取得最小值;由圖5可知:盤的徑向變形在銷釘固定處最小,在盤緣處最大.
選取渦輪盤不同位置的溫度和對(duì)流換熱系數(shù)(根據(jù)各對(duì)象與燃?xì)鉁囟鹊膿Q熱特點(diǎn)[2]計(jì)算),以及轉(zhuǎn)速ω、材料密度ρ、彈性模量E和泊松比υ為隨機(jī)變量,其數(shù)字分布特征如表2所示.假設(shè)所有參數(shù)均服從正態(tài)分布且相互獨(dú)立.
將表2中的隨機(jī)變量分布特征導(dǎo)入有限元模型中,用MC法對(duì)各隨機(jī)參數(shù)進(jìn)行抽樣,得到一系列的樣本點(diǎn),再用這些樣本點(diǎn)對(duì)有限元模型進(jìn)行非線性動(dòng)態(tài)確定性分析,得到281組試驗(yàn)樣本,其中輸出響應(yīng)Y的仿真歷史見(jiàn)圖6.
再利用Box Behnken矩陣法[11]擬合響應(yīng)面函數(shù)式(1),得到響應(yīng)面函數(shù)系數(shù),進(jìn)而建立響應(yīng)面模型,如式(7)所示(忽略系數(shù)小于10-5的項(xiàng)).
表2 輪盤的隨機(jī)變量選取
圖6 有限元模型輸出響應(yīng)Y仿真歷史
將渦輪盤徑向變形的響應(yīng)面模型代替有限元模型,利用MC法進(jìn)行1萬(wàn)次抽樣仿真,其中Y仿真歷史曲線圖和直方圖如圖7和圖8所示.由圖8可知:Y滿足正態(tài)分布,均值為1.2361 mm,方差為 0.030131 mm.
根據(jù)式(2)~式(6)可得渦輪盤徑向變形量在不同設(shè)計(jì)值δ下的可靠度值,如表3所示.
由表3可知:當(dāng)δ=1.33mm時(shí),失效數(shù)為15,可靠度R=0.99852,基本滿足工程需要.因此,在渦輪盤變形設(shè)計(jì)時(shí),建議δ=1.33 mm.
圖7 響應(yīng)面模型輸出響應(yīng)Y仿真歷史
圖8 輪盤徑向變化直方圖
表3 不同設(shè)計(jì)值δ下的可靠度
同時(shí),也可以對(duì)盤的變形進(jìn)行逆概率分析,通過(guò)逆概率分析,其結(jié)果如表4所示.
表4 逆概率分析結(jié)果(部分)
由表4可以看出,通過(guò)逆概率分析,可以得到每個(gè)可靠度下的各個(gè)隨機(jī)參數(shù)的極限值,為渦輪盤設(shè)計(jì)與優(yōu)化提供量化參考.
根據(jù)式(5)和式(6)可計(jì)算出各隨機(jī)變量對(duì)渦輪盤徑向變形的靈敏度和影響概率,如圖9和表5所示.
圖9 隨機(jī)變量對(duì)渦輪盤徑向變形的靈敏度分析
表5 渦輪盤徑向變形靈敏度分析結(jié)果
在圖9中,最重要的隨機(jī)輸入變量(靈敏度最大)在最左邊,其他依次向右排列.在表5中,靈敏度有正負(fù)之分,正表示隨機(jī)輸出變量隨輸入?yún)⒘空兓?,?fù)則表示反變化.由圖9和表5可得以下結(jié)論:
轉(zhuǎn)子轉(zhuǎn)速ω對(duì)輪盤徑向變形的影響最大,靈敏度達(dá)到0.9,影響概率為0.536.符合本文采用的二代發(fā)動(dòng)機(jī)輪盤的工程實(shí)際情況;輪盤溫度Ta1,Tb2,Ta3和 Tb1對(duì)渦輪徑向變形也有非常重要的影響.泊松比、彈性模量、導(dǎo)熱系數(shù)等對(duì)輪盤變形的影響很小.所以,在對(duì)渦輪盤徑向變形設(shè)計(jì)與優(yōu)化中,首要考慮轉(zhuǎn)速、輪緣溫度的影響.同時(shí),也為葉尖徑向運(yùn)行間隙設(shè)計(jì)與控制提供了量化參考.
本文考慮材料屬性的非線性和熱載荷、離心載荷的動(dòng)態(tài)性,對(duì)輪盤徑向變形進(jìn)行了概率分析,得出如下結(jié)論:
1)通過(guò)渦輪盤徑向變形的確定性分析,得到了輪盤徑向變形量及其隨時(shí)間的變化規(guī)律,并且當(dāng)t=178 s時(shí)(燃?xì)鉁囟茸罡?、轉(zhuǎn)速最大時(shí)刻)渦輪盤變形量達(dá)到全局最大值:1.236 mm.
2)通過(guò)概率分析和逆概率分析,建立了輪盤變形概率分析的響應(yīng)面模型;得到了輪盤徑向變形分布特征:均值為 1.236 1 mm,方差為0.030131 mm,且服從正態(tài)分布;還得到不同設(shè)計(jì)值δ下的徑向變形的失效數(shù)和可靠度,建議δ=1.33 mm,此時(shí)可靠度 R=0.998 52,符合設(shè)計(jì)要求.另外還得到不同可靠度下的各隨機(jī)參數(shù)的極限值,為渦輪盤徑向變形設(shè)計(jì)和優(yōu)化,以及葉尖徑向運(yùn)行間隙的設(shè)計(jì)和控制提供了依據(jù).
3)通過(guò)靈敏度分析,得到了各隨機(jī)參數(shù)的靈敏度和影響概率,找到了影響徑向變形變化的主導(dǎo)因素(轉(zhuǎn)速 ω)和主要因素(溫度 Ta1,Tb2,Ta3和Tb1),為渦輪盤徑向變形設(shè)計(jì)和優(yōu)化,以及葉尖徑向運(yùn)行間隙的設(shè)計(jì)和控制提供了量化參考.
References)
[1]Hennecke D K,Trappmann K.Turbine tip clearance control in gas turbine engines[R].NASA N83-29254,1983
[2]Lattime S,Steinetz B.Turbine engine clearance control systems:current practices and future directions[R].AIAA 2002-3790,2002
[3]漆文凱,陳偉.某型航空發(fā)動(dòng)機(jī)高壓渦輪葉尖間隙數(shù)值分析[J].南京航空航天大學(xué)學(xué)報(bào),2003,35(1):63-67 Qi Wenkai,Chen Wei.Tip clearance numerical analysis of an aero-engine HPT[J].Journal of Nanjing University of Aeronautics & Astronautics,2003,35(1):63-67(in Chinese)
[4]Pilidis P,Maccallum N R L.Models for predicting tip clearance changes in gas turbines[R].NASA N83-29258,1983
[5]NSSA Glenn Research Center.HPT clearance control[R].NASA/CR-2005-213970,2005
[6]豈興明,樸英,祝劍虹.某型航空發(fā)動(dòng)機(jī)高壓渦輪葉頂間隙三維數(shù)值分析[J].航空動(dòng)力學(xué)報(bào),2008,3(5):904-908 Qi Xingming,Pu Ying,Zhu Jianhong.3-D numerical analysis of the tip clearance of an aero-engine high pressure turbine[J].Journal of Aerospace Power,2008,3(5):904-908(in Chinese)[7]Eeom Y S,Yoo K S,Park J Y.Reliability-based topology optimization using a standard response surface method for three-dimensional structures[J].Structural and Multidisciplinary Optimization,2011,43(2):287-295
[8]李昌,韓興.基于響應(yīng)面法齒輪嚙合傳動(dòng)可靠性靈敏度分析[J].航空動(dòng)力學(xué)報(bào),2011,26(3):711-715 Li Chang,Han Xing.Analysis of reliability sensitivity for gear engagement based on response surface methods[J].Journal of Aerospace Power,2011,26(3):711-715(in Chinese)
[9]Kypuros J A,Melcher K J.A reduced model for prediction of thermal and rotational effects on turbine tip clearance[R].NASA/TM-2003-212226,2003
[10]Annette E N,Christoph W M,Stehan S.Modeling and validation of the thermal effects on gas turbine transients[J].Journal of Engineering for Gas Turbines and Power,2005,127(3):564-572
[11]Hakan E.Design sensitivity analysis of structures based upon the singular value decomposition[J].Computer Methods in Applied Mechanics and Engineering,2002,191(5):3459 -3476