基于時(shí)間推進(jìn)的軸流壓氣機(jī)氣動(dòng)性能計(jì)算模型

陳忠軍，郭 晉，胡 駿

(1.廣州民航職業(yè)技術(shù)學(xué)院，廣州510403；2.南京航空航天大學(xué)能源與動(dòng)力學(xué)院，南京210016)

1 引言

在壓氣機(jī)設(shè)計(jì)階段，準(zhǔn)確預(yù)測(cè)其氣動(dòng)性能至關(guān)重要。近幾十年，隨著計(jì)算流體力學(xué)(CFD)的長(zhǎng)足發(fā)展，使得全三維粘性計(jì)算逐漸應(yīng)用于壓氣機(jī)設(shè)計(jì)體系，為評(píng)估壓氣機(jī)氣動(dòng)性能提供了數(shù)值手段。盡管如此，基于子午流面(S2流面)的壓氣機(jī)氣動(dòng)性能準(zhǔn)三維計(jì)算方法，仍是壓氣機(jī)設(shè)計(jì)階段不可或缺的重要工具。S2流面的概念由吳仲華于1952年提出[1]，目前比較流行的S2流面準(zhǔn)三維計(jì)算方法為流線曲率法[2-3]和流函數(shù)法[4]。這兩種方法理論較為完善，且已應(yīng)用于工程實(shí)際。然而，由于上述方法均基于空間推進(jìn)，當(dāng)壓氣機(jī)內(nèi)部流場(chǎng)軸向馬赫數(shù)局部超聲時(shí)，其將不再有效。為彌補(bǔ)傳統(tǒng)S2流面計(jì)算方法的不足，一種基于時(shí)間推進(jìn)求解帶體積力源項(xiàng)的周向平均的Euler方程或Navier-Stokes方程的計(jì)算方法被提出[5-8]。該方法可同時(shí)適用于亞、跨及超聲速流動(dòng)。

本文基于時(shí)間推進(jìn)這一思路，采用發(fā)展較為完善、考慮因素較為全面的落后角及損失模型，并納入徑向摻混的影響，同時(shí)引入Koch的最大靜壓升系數(shù)法[9]作為判斷壓氣機(jī)穩(wěn)定邊界的手段，建立了一個(gè)基于時(shí)間推進(jìn)的軸流壓氣機(jī)氣動(dòng)性能計(jì)算模型。利用該模型對(duì)兩臺(tái)四級(jí)低速軸流壓氣機(jī)及兩臺(tái)單級(jí)跨聲速壓氣機(jī)的總體特性展開(kāi)計(jì)算，并將計(jì)算結(jié)果與相關(guān)實(shí)驗(yàn)數(shù)據(jù)進(jìn)行了對(duì)比分析。

2 數(shù)值方法

2.1 主控方程

無(wú)葉區(qū)內(nèi)的氣流流動(dòng)采用絕對(duì)圓柱坐標(biāo)系下的軸對(duì)稱(chēng)準(zhǔn)三維流動(dòng)方程求解。同時(shí)，為考慮主流區(qū)內(nèi)的徑向摻混現(xiàn)象，基于Gallimore等[10-11]的思想，在Euler方程組中引入湍流粘性應(yīng)力及湍流熱傳導(dǎo)項(xiàng)以描述該現(xiàn)象。葉片區(qū)內(nèi)，用體積力源項(xiàng)代替葉型作用，并采用葉型堵塞系數(shù)反映葉型厚度對(duì)葉片區(qū)域內(nèi)部流動(dòng)的影響。無(wú)葉區(qū)及葉片區(qū)流動(dòng)方程可統(tǒng)一寫(xiě)為式(1)：

式中：Ev=τθθ=τzθ=τθz=ρ為密度，e為總能量，p為靜壓，v為速度，τ為湍流粘性應(yīng)力，q為湍流熱流密度，T為靜溫，μt為渦粘系數(shù)，kt為湍流熱傳導(dǎo)系數(shù)，Prt為湍流普朗特?cái)?shù)，??v→表示軸對(duì)稱(chēng)假設(shè)下的速度散度。Ev和Gv中的下標(biāo)1～3分別表示軸向(z)、周向(θ)及徑向(r)。葉型堵塞系數(shù)b定義為，N為葉片數(shù)，θp、θs分別表示相鄰葉型壓力面及吸力面的周向坐標(biāo)；無(wú)葉區(qū)中，b指定為1，Sb及SF中各項(xiàng)均為0。ω為壓氣機(jī)轉(zhuǎn)速。F表示葉片體積力源項(xiàng)，代表葉型對(duì)氣流產(chǎn)生的作用。

采用有限體積法對(duì)方程(1)進(jìn)行離散，時(shí)間項(xiàng)采用五步龍格-庫(kù)塔法進(jìn)行顯式時(shí)間推進(jìn)，并采用當(dāng)?shù)貢r(shí)間步長(zhǎng)加速收斂過(guò)程。無(wú)粘通量采用LDFSS格式求解[12]。渦粘系數(shù)的計(jì)算參考Gallimore等[10-11]的方法，Prt給定為0.9。粘性應(yīng)力計(jì)算時(shí)忽略軸向偏導(dǎo)數(shù)的影響。為節(jié)約計(jì)算時(shí)間，模型將壁面邊界處理為無(wú)粘滑移邊界。為考慮端壁附面層對(duì)流動(dòng)造成的堵塞影響，對(duì)輪轂及輪緣徑向坐標(biāo)做經(jīng)驗(yàn)性修正。

2.2 體積力計(jì)算

體積力計(jì)算方法遵循了Marble的思想[13]，將體積力分解為垂直于相對(duì)速度的無(wú)粘力及平行反向于相對(duì)速度的有粘力兩個(gè)部分。有粘力由葉排進(jìn)出口熵增求得，而熵增與基元葉型損失系數(shù)?之間存在一一對(duì)應(yīng)關(guān)系，因此可建立損失系數(shù)與粘性力之間的關(guān)系。無(wú)粘力則通過(guò)先確定葉片區(qū)內(nèi)的周向速度分布，再根據(jù)穩(wěn)態(tài)周向動(dòng)量方程獲得。周向速度的確定方法為：基于Taddei等[8]提出的方法，在葉片前緣區(qū)域依據(jù)葉片進(jìn)口周向速度給定線性分布，在剩余葉片區(qū)域內(nèi)依據(jù)葉片出口落后角δ給定氣流角分布，再根據(jù)當(dāng)前時(shí)間步子午面速度計(jì)算獲得絕對(duì)周向速度。通過(guò)上述方式可建立落后角及損失系數(shù)與體積力之間的關(guān)系。本文所采用的落后角及損失模型主要基于文獻(xiàn)[14]、[15]的研究工作，將落后角及損失系數(shù)分解為對(duì)應(yīng)基準(zhǔn)攻角的基準(zhǔn)部分和偏離基準(zhǔn)攻角的非基準(zhǔn)部分；同時(shí)對(duì)落后角及損失系數(shù)進(jìn)行三維修正，并考慮了雷諾數(shù)對(duì)其的影響。

2.3 穩(wěn)定邊界判定

采用Koch提出的最大靜壓升系數(shù)法進(jìn)行壓氣機(jī)穩(wěn)定邊界的判定[9]。Koch基于大量的實(shí)驗(yàn)數(shù)據(jù)總結(jié)了一組預(yù)測(cè)壓氣機(jī)最大失速壓升能力的曲線，將壓氣機(jī)的失速靜壓升系數(shù)與基元葉柵的尺寸、葉尖和軸向間隙以及雷諾數(shù)等關(guān)聯(lián)起來(lái)。Koch定義的靜壓升系數(shù)見(jiàn)公式(2)。此外，CFD計(jì)算中也常采用計(jì)算發(fā)散作為判穩(wěn)準(zhǔn)則。由于模型為基于時(shí)間推進(jìn)的CFD技術(shù)，且Koch所總結(jié)的預(yù)測(cè)壓升能力的曲線較依賴(lài)于工程經(jīng)驗(yàn)，其有效適用范圍有限。因此，若計(jì)算提前發(fā)散，則認(rèn)定為壓氣機(jī)失穩(wěn)。

3 算例驗(yàn)證

為展現(xiàn)所建立模型的工程應(yīng)用能力，對(duì)南京航空航天大學(xué)兩臺(tái)四級(jí)低速軸流壓氣機(jī)LSC4-1、LSC4-2，以及NASA兩臺(tái)單級(jí)跨聲速壓氣機(jī)NASA Stage 35、NASA Stage 37展開(kāi)計(jì)算，并與相關(guān)實(shí)驗(yàn)數(shù)據(jù)進(jìn)行對(duì)比分析。

3.1 四級(jí)低速軸流壓氣機(jī)對(duì)比分析

LSC4-1、LSC4-2是為現(xiàn)代高壓壓氣機(jī)后面級(jí)的低速模擬研究而設(shè)計(jì)，第三級(jí)為其模擬級(jí)。設(shè)計(jì)轉(zhuǎn)速為900 r/min，外徑為1.5 m，輪轂比為0.88；導(dǎo)向器數(shù)目為60，轉(zhuǎn)子、靜子葉片數(shù)分別為72和120。兩臺(tái)壓氣機(jī)的區(qū)別在于第三級(jí)的轉(zhuǎn)子和靜子造型不同。文獻(xiàn)[16]給出了壓氣機(jī)設(shè)計(jì)結(jié)構(gòu)參數(shù)和詳細(xì)的實(shí)驗(yàn)測(cè)量結(jié)果。

圖1、圖2分別為模型計(jì)算獲得的不同轉(zhuǎn)速下LSC4-1和LSC4-2的總體特性與實(shí)驗(yàn)數(shù)據(jù)的對(duì)比，可見(jiàn)計(jì)算獲得的總體特性及穩(wěn)定邊界與實(shí)驗(yàn)結(jié)果吻合良好。圖3為L(zhǎng)SC4-1不同轉(zhuǎn)速下各級(jí)有效靜壓升系數(shù)隨流量的變化關(guān)系。可見(jiàn)，隨著流量的增加，各級(jí)有效靜壓升系數(shù)不斷接近并最終達(dá)到級(jí)失速靜壓升系數(shù)。

圖4示出了設(shè)計(jì)點(diǎn)模型獲得的兩臺(tái)壓氣機(jī)第三級(jí)靜子進(jìn)口攻角及損失系數(shù)沿葉高分布與實(shí)驗(yàn)數(shù)據(jù)的對(duì)比，可見(jiàn)計(jì)算結(jié)果與實(shí)驗(yàn)結(jié)果的分布形態(tài)較為一致。

3.2 單級(jí)跨聲速壓氣機(jī)對(duì)比分析

NASA Stage 35和NASA Stage 37是NASA于20世紀(jì)70年代設(shè)計(jì)的典型單級(jí)跨聲速壓氣機(jī)[17-18]。Stage 35的設(shè)計(jì)轉(zhuǎn)速為17 188.7 r/min，質(zhì)量流量為20.19 kg/s，總壓比為1.82，轉(zhuǎn)子、靜子葉片數(shù)分別為36和46。Stage 37的設(shè)計(jì)轉(zhuǎn)速和葉片數(shù)均與Stage 35的相同，質(zhì)量流量為20.20 kg/s，總壓比為2.05。

圖5、圖6分別為Stage 35和Stage 37模型計(jì)算獲得的總體特性與實(shí)驗(yàn)數(shù)據(jù)的對(duì)比。可發(fā)現(xiàn)，程序捕獲的兩臺(tái)壓氣機(jī)高轉(zhuǎn)速下的堵點(diǎn)流量與實(shí)驗(yàn)結(jié)果基本相當(dāng)；堵點(diǎn)流量均由計(jì)算發(fā)散確定，與實(shí)驗(yàn)數(shù)據(jù)吻合較好。對(duì)于Stage 35，文獻(xiàn)[17]中記錄的100%轉(zhuǎn)速下實(shí)驗(yàn)測(cè)得的失速流量為18.20 kg/s，但文獻(xiàn)[19]中指出在該點(diǎn)壓氣機(jī)可能已經(jīng)進(jìn)入旋轉(zhuǎn)失速狀態(tài)，因此該點(diǎn)實(shí)驗(yàn)工況在圖5中并未給出。總體上看，計(jì)算獲得的特性曲線與實(shí)驗(yàn)數(shù)據(jù)吻合良好，初步說(shuō)明了該模型對(duì)跨聲壓氣機(jī)的適用性。

4 結(jié)論

利用基于時(shí)間推進(jìn)的技術(shù)，通過(guò)求解周向平均帶體積力源項(xiàng)的流動(dòng)控制方程，融合發(fā)展較為完善的落后角及損失模型，考慮徑向摻混影響，同時(shí)引入Koch穩(wěn)定邊界模型，建立了一個(gè)壓氣機(jī)性能計(jì)算模型。該性能計(jì)算模型能較為準(zhǔn)確地預(yù)測(cè)多級(jí)壓氣機(jī)的總體特性，且具有較好的工程應(yīng)用前景。后續(xù)將繼續(xù)開(kāi)展對(duì)多級(jí)跨聲速壓氣機(jī)的計(jì)算分析工作，以進(jìn)一步提高模型的預(yù)測(cè)精度，擴(kuò)大模型的應(yīng)用范圍。

參考文獻(xiàn)：

[1]Wu C H.A general theory of three dimensional flow in subsonic and supersonic turbomachines of axial,radial,and mixed flow types[R].NACA TN 2604，1952.

[2]Novak R A.Streamline curvature computing procedures for fluid-flow problems[J].ASME Journal of Engineering for Gas Turbines and Power，1967，89(4)：478—490.

[3]胡 駿.均勻與非均勻進(jìn)氣條件下多級(jí)軸流壓氣機(jī)性能計(jì)算-均勻進(jìn)氣及徑向畸變的影響[J].航空動(dòng)力學(xué)報(bào)，2000，15(3)：225—228.

[4]Marsh H.A digital computer program for the through-flow fluid mechanics in an arbitrary turbomachine using a ma?trix method[R].London：Aeronautical Research Council R and M，1968.

[5]Baralon S，Hall U，Eriksson L E.Viscous throughflow mod?elling of transonic compressors using a time-marching fi?nite volume solver[R].ISABE 97-7071，1997.

[6]Simon J F，Léonard O.A throughflow analysis tool based on the Navier-Stokes equations[C]//.Proceedings of 6th European Turbomachinery Conference.France Lille，2005.

[7]PaccianiR，RubechiniF，MarconciniM，etal.A CFD-based throughflow method with an explicit body force model and an adaptive formulation for the S2 stream?surface[J].Journal of Power and Energy，2016，230(1)：16—28.

[8]Taddei S R，Larocca F.CFD-based analysis of multistage throughflow surfaces with incidence[J].Mechanics Re?search Communications，2013，47：6—10.

[9]Koch C C.Stalling pressure rise capability of axial flow compressor stages[J].ASME Journal of Engineering for Gas Turbines and Power，1981，103(4)：645—656.

[10]Gallimore S J，Cumpsty N A.Spanwise mixing in multi?stage axial flow compressors:part I:experimental investi?gation[R].ASME 86-GT-20，1986.

[11]Gallimore S J，Cumpsty N A.Spanwise mixing in multi?stage axial flow compressors:part II:throughflow calcula?tions including mixing[R].ASME 86-GT-21，1986.

[12]Edwards J R.A low-diffusion flux-splitting scheme for Navier-Stokes calculations[J].Computers&Fluids,1997，26(6)：635—659.

[13]Marble F E.Three-dimensional flow in turbomachines,aerodynamics of turbines and compressors[M].Princeton：Princeton University Press，1964.

[14]趙 勇.風(fēng)扇/壓氣機(jī)非設(shè)計(jì)點(diǎn)性能計(jì)算和進(jìn)氣畸變影響預(yù)測(cè)方法研究[D].南京：南京航空航天大學(xué)，2008.

[15]Banjac M，Petrovic M V，Wiedermann A.Secondary flows,endwall effects,and stall detection in axial compressor de?sign[R].ASME GT2014-25115，2014.

[16]張晨凱.多級(jí)軸流壓氣機(jī)內(nèi)復(fù)雜流動(dòng)結(jié)構(gòu)的實(shí)驗(yàn)和數(shù)值研究[D].南京：南京航空航天大學(xué)，2015.

[17]Reid L，Moore R D.Performance of a single stage axi?al-flow transonic compressor with rotor and stator aspect ratio of 1.19 and 1.26,respectively,andwith design pressure ratio of 1.82[R].NASA TP-1338，1978.

[18]Moore R D，Reid L.Performance of single-stage axi?al-flow transonic compressor with rotor and stator aspect ratios of 1.19 and 1.26,respectively,and with design pres?sure ratio of 2.05[R].NASA TP-1659，1980.

[19]Chima R V.SWIFT code assessment for two similar tran?sonic compressors[R].AIAA 2009-1058，2009.