孫 逸，葛 寧，舒 杰

（南京航空航天大學(xué)能源與動(dòng)力學(xué)院，南京210016）

0 引言

在航空發(fā)動(dòng)機(jī)總體設(shè)計(jì)初期，通流方法可以快速得到發(fā)動(dòng)機(jī)性能以及子午面流動(dòng)特征[1-2]，仍是現(xiàn)代航空發(fā)動(dòng)機(jī)設(shè)計(jì)工作的基石和重要的設(shè)計(jì)工具[3-4]。Simon[5]提出基于Navier-Stokes方程的周向平均通流模型，實(shí)現(xiàn)適用于軸流壓氣機(jī)的定常與非定常通流計(jì)算，并對(duì)周向平均Navier-Stokes方程中出現(xiàn)的各種附加項(xiàng)做了詳細(xì)論證；金海良[6]對(duì)某跨聲速單級(jí)風(fēng)扇、某子午加速風(fēng)機(jī)和某跨聲速雙級(jí)風(fēng)扇的研究表明，基于Navier-Stokes方程的周向平均通流模型計(jì)算得到的子午流場(chǎng)和徑向參數(shù)分布與3維計(jì)算結(jié)果十分吻合，在近端壁區(qū)域具有較高的計(jì)算精度；周向平均通流模型的另一重要應(yīng)用是發(fā)動(dòng)機(jī)整機(jī)數(shù)值仿真計(jì)算，可用于發(fā)動(dòng)機(jī)各部件間的流場(chǎng)匹配，研究發(fā)動(dòng)機(jī)各部件間的相互影響[7]；Friederike[8]針對(duì)某大涵道比渦扇發(fā)動(dòng)機(jī)提出1種0D/2D耦合航空發(fā)動(dòng)機(jī)性能預(yù)估方法；Denton等[9]基于40多年研究經(jīng)驗(yàn)開(kāi)發(fā)了多層次葉輪機(jī)械設(shè)計(jì)程序Multall，其中通流計(jì)算作為重要環(huán)節(jié)，可快速分析葉片損失和效率。由于現(xiàn)階段受CFD技術(shù)和計(jì)算機(jī)性能的限制，采用3維數(shù)值模擬方法計(jì)算發(fā)動(dòng)機(jī)整機(jī)需要大量的計(jì)算資源和時(shí)間，而采用通流計(jì)算方法能夠節(jié)約計(jì)算時(shí)間，快速發(fā)現(xiàn)各部件匹配過(guò)程中的問(wèn)題，具有一定的工程應(yīng)用價(jià)值。本文選取某渦扇發(fā)動(dòng)機(jī)低壓部件為研究對(duì)象，采用通流模型耦合計(jì)算，研究各部件相互作用以及過(guò)渡段幾何參數(shù)對(duì)發(fā)動(dòng)機(jī)總體性能的影響。通流模型可用于發(fā)動(dòng)機(jī)總體設(shè)計(jì)初期對(duì)發(fā)動(dòng)機(jī)性能快速預(yù)估，實(shí)現(xiàn)各部件性能匹配。

1 基于Navier-Stokes方程的周向平均通流模型

在相對(duì)柱坐標(biāo)系中，對(duì)Navier-Stokes方程進(jìn)行周向平均[10]，得到通流模型的控制方程為

其中

式中：Q為守恒量；F和H為對(duì)流（無(wú)黏）通量；Fv和Hv為擴(kuò)散（黏性）通量；S為Navier-Stokes方程組在相對(duì)柱坐標(biāo)系下導(dǎo)出的源項(xiàng)；SB為無(wú)黏葉片力項(xiàng)；Sv為黏性葉片力項(xiàng)；fB和fv分別為模化后的無(wú)黏葉片力和黏性葉片力；x、r、θ分別為軸向、徑向和周向坐標(biāo)方向；ρ為密度；v為相對(duì)速度；p為靜壓；e為單位質(zhì)量總能；w為轉(zhuǎn)速；τ為黏性應(yīng)力；q為熱流量；定義葉片堵塞系數(shù)

式中：θs和θp分別為葉片吸力面和壓力面上的角坐標(biāo)；N為葉片數(shù)。

堵塞系數(shù)b在葉片區(qū)小于1，在非葉片區(qū)等于1。基于無(wú)黏葉片力始終垂直于平均流面的假設(shè)，本文將無(wú)黏葉片力的求解分為2步[11-12]：（1）在動(dòng)量方程中不計(jì)入無(wú)黏葉片力；（2）在動(dòng)量方程中僅計(jì)入無(wú)黏葉片力，從而可以修正第1步求得的動(dòng)量，使其與平均流面相切。這就避免了顯式求解無(wú)黏葉片力。黏性葉片力的?；捎弥姆植紦p失模型[13]，引入1種分布的體積力來(lái)計(jì)入流動(dòng)損失的影響。黏性葉片力平行于相對(duì)速度并且反向，所作功僅產(chǎn)生熵增，采用總壓損失系數(shù)計(jì)算葉排進(jìn)口到出口熵增。本文采用時(shí)間推進(jìn)的有限體積法求解上述控制方程，用具有較高間斷分辨率的Roe通量差分分裂方法[14]對(duì)無(wú)黏對(duì)流通量進(jìn)行離散，而黏性通量的離散則是利用網(wǎng)格中心的流場(chǎng)變量將網(wǎng)格面上的黏性通量表示出來(lái)，紊流模型采用B-L模型。

2 通流模型驗(yàn)證

通流耦合計(jì)算中首先要保證通流模型的準(zhǔn)確性，激波位置和激波強(qiáng)度對(duì)壓氣機(jī)流量、效率均有較大影響。因此，選取NASA Rotor67和某3級(jí)風(fēng)扇為驗(yàn)證算例。

2.1 Rotor67算例

Rotor67是1個(gè)低展弦比跨聲速轉(zhuǎn)子，Strazisar[15]于1989年做了細(xì)致的試驗(yàn)。設(shè)計(jì)轉(zhuǎn)速為16043 r/min，設(shè)計(jì)轉(zhuǎn)速下最高效率點(diǎn)流量為34.573 kg/s，壓比為1.642，效率為0.93。最高效率點(diǎn)通流模型計(jì)算結(jié)果與試驗(yàn)值對(duì)比見(jiàn)表1，可見(jiàn)通流計(jì)算結(jié)果與試驗(yàn)值基本吻合。根據(jù)Strazisar的描述，由于在轉(zhuǎn)子出口邊界層處測(cè)點(diǎn)較少，試驗(yàn)測(cè)量的等熵效率可能偏高。

表1 Rotor67通流計(jì)算與試驗(yàn)結(jié)果對(duì)比

Rotor67在最高效率點(diǎn)時(shí)，通流模型計(jì)算得到的葉排出口截面參數(shù)（分別為氣流角、總溫、靜壓、總壓）徑向分布與試驗(yàn)值對(duì)比如圖1所示。從圖中可見(jiàn)，通流計(jì)算得到的靜壓分布與試驗(yàn)值十分吻合，但由于葉尖存在激波，氣流角和總溫在葉尖處與試驗(yàn)值有一定偏差。總壓徑向分布與試驗(yàn)基本吻合，而且能夠預(yù)測(cè)出端壁附面層附近的總壓變化趨勢(shì)。

圖1 Rotor67葉排出口徑向參數(shù)分布

Rotor67在最高效率點(diǎn)的馬赫數(shù)等值線(xiàn)如圖2所示，從圖中可見(jiàn)，葉片區(qū)域存在2道激波，第1道為葉尖進(jìn)口斜激波，第2道為尾緣附近的正激波，波前馬赫數(shù)約為1.4。通流計(jì)算得到的Rotor67在100%轉(zhuǎn)速下的特性線(xiàn)如圖3、4所示。從圖中可見(jiàn)，二者結(jié)果基本吻合，驗(yàn)證了通流模型的準(zhǔn)確性和可靠性。

圖2 Rotor67最高效率點(diǎn)馬赫數(shù)等值線(xiàn)

圖3 流量-效率特性線(xiàn)

圖4 流量-壓比特性線(xiàn)

2.2 某軸流風(fēng)扇算例

選取某帶進(jìn)口導(dǎo)葉的3級(jí)軸流風(fēng)扇進(jìn)行驗(yàn)證，其設(shè)計(jì)狀態(tài)下壓比為3.5，流量為81.0 kg/s。風(fēng)扇子午流場(chǎng)相對(duì)馬赫數(shù)如圖5所示。從圖中可見(jiàn)，第1級(jí)轉(zhuǎn)子葉尖最高波前馬赫數(shù)約為1.5，較強(qiáng)的激波使得葉尖壓比相對(duì)試驗(yàn)值偏高，同時(shí)較強(qiáng)的激波-邊界層相互干擾損失導(dǎo)致轉(zhuǎn)子葉尖效率較低。

圖5 風(fēng)扇子午流場(chǎng)相對(duì)馬赫數(shù)

計(jì)算了該軸流風(fēng)扇在100%轉(zhuǎn)速下的性能特性，并與試驗(yàn)值進(jìn)行了對(duì)比，流量-效率特性線(xiàn)與流量-壓比特性線(xiàn)分別如圖6、7所示。從圖中可見(jiàn)，在風(fēng)扇近設(shè)計(jì)點(diǎn)計(jì)算得到的壓比為3.53，效率為0.84。由于周向平均通流模型捕捉到的激波其物理本質(zhì)與3維計(jì)算有所區(qū)別，因此捕捉到的激波比實(shí)際稍強(qiáng)，導(dǎo)致計(jì)算結(jié)果效率偏低。通流模型計(jì)算得到的特性線(xiàn)與試驗(yàn)值基本吻合，可用于多級(jí)壓氣機(jī)特性計(jì)算。

圖6 風(fēng)扇流量-效率特性

圖7 風(fēng)扇流量-壓比特性

3 通流耦合計(jì)算

3.1 研究對(duì)象

選取某渦扇發(fā)動(dòng)機(jī)除核心機(jī)3大部件外的低壓部件為研究對(duì)象，包括進(jìn)口帶導(dǎo)葉的3級(jí)風(fēng)扇、出口帶支板的2級(jí)低壓渦輪、內(nèi)外涵道、摻混段以及尾噴管。其流道和葉型根據(jù)該發(fā)動(dòng)機(jī)總體及部件性能參數(shù)試驗(yàn)值（見(jiàn)表2），采用流道設(shè)計(jì)程序計(jì)算得到。子午面流道如圖8所示。風(fēng)扇出口分為內(nèi)外涵道2個(gè)計(jì)算域，并通過(guò)外涵道與渦輪后摻混段相連，摻混段下游為尾噴管。內(nèi)、外涵冷熱流體摻混是1個(gè)重要的氣動(dòng)熱力過(guò)程，本文采用流量加權(quán)后定比熱方法計(jì)算摻混后流場(chǎng)。耦合計(jì)算邊界條件均按試驗(yàn)值給定，包括風(fēng)扇進(jìn)口總溫、總壓，內(nèi)涵出口靜壓，低壓渦輪進(jìn)口總溫、總壓，噴管出口給定標(biāo)準(zhǔn)大氣壓力。

表2 通流計(jì)算結(jié)果與試驗(yàn)值對(duì)比

圖8 幾何模型子午面

3.2 計(jì)算結(jié)果及分析

計(jì)算結(jié)果表明，各部件單獨(dú)計(jì)算時(shí)性能均能滿(mǎn)足總體設(shè)計(jì)要求。但多個(gè)部件耦合數(shù)值計(jì)算時(shí)，發(fā)動(dòng)機(jī)總體性能參數(shù)往往會(huì)偏離設(shè)計(jì)狀態(tài)。本文通流耦合計(jì)算的核心就是研究各部件匹配以及過(guò)渡段對(duì)發(fā)動(dòng)機(jī)總體性能的影響。

研究發(fā)現(xiàn)，風(fēng)扇出口內(nèi)、外涵分流板、渦輪出口摻混段整流板是影響部件匹配和發(fā)動(dòng)機(jī)總體性能最重要的2個(gè)過(guò)渡段幾何參數(shù)。分流板不僅分配內(nèi)、外涵流量，控制涵道比，并且對(duì)風(fēng)扇的工作狀態(tài)有較大影響。摻混段整流板對(duì)內(nèi)、外涵摻混有一定影響，而且可以抑制尾椎流動(dòng)分離，對(duì)渦輪和噴管工作狀態(tài)也至關(guān)重要。

通過(guò)調(diào)節(jié)分流板和整流板幾何參數(shù)，使得各部件性能結(jié)算結(jié)果接近設(shè)計(jì)狀態(tài)。通流計(jì)算子午流場(chǎng)相對(duì)馬赫數(shù)如圖9所示。從圖中可見(jiàn)，風(fēng)扇轉(zhuǎn)子葉尖均存在激波，第1級(jí)轉(zhuǎn)子中波前馬赫數(shù)約為1.4，較強(qiáng)的激波損失使得風(fēng)扇效率計(jì)算結(jié)果偏低3.57%。

圖9 子午流場(chǎng)相對(duì)馬赫數(shù)

從表2列出的風(fēng)扇、低壓渦輪、內(nèi)外涵以及尾噴管近設(shè)計(jì)點(diǎn)性能參數(shù)與試驗(yàn)值對(duì)比可見(jiàn)，風(fēng)扇計(jì)算流量偏大3.14%，增壓比偏小0.56%，效率偏低3.57%。低壓渦輪計(jì)算流量偏大4.93%，落壓比偏大0.92%，效率偏低3.41%。通流模型得到的內(nèi)涵流量偏大2.28%，外涵流量偏小0.48%，導(dǎo)致涵道比偏小2.55%。噴管出口計(jì)算流量和排氣速度比試驗(yàn)值偏大，導(dǎo)致發(fā)動(dòng)機(jī)計(jì)算推力比試驗(yàn)值偏大2.67%，耗油率偏低2.52%。從誤差分析來(lái)看，本文通流耦合計(jì)算在近設(shè)計(jì)點(diǎn)能夠滿(mǎn)足發(fā)動(dòng)機(jī)總體性能初步要求，最大誤差不超過(guò)5.0%。

4 結(jié)論

本文通過(guò)對(duì)某渦扇發(fā)動(dòng)機(jī)低壓部件通流耦合計(jì)算，得到如下結(jié)論：

（1）采用周向平均通流模型對(duì)Rotor67和某3級(jí)風(fēng)扇進(jìn)行了計(jì)算并與試驗(yàn)值對(duì)比，結(jié)果證明通流模型能夠正確捕捉激波強(qiáng)度與位置，預(yù)測(cè)的部件特性曲線(xiàn)趨勢(shì)與試驗(yàn)基本吻合。

（2）對(duì)某渦扇發(fā)動(dòng)機(jī)風(fēng)扇、低壓渦輪、內(nèi)外涵道、摻混段和尾噴管等部件進(jìn)行了通流耦合計(jì)算，研究了各部件匹配與過(guò)渡段幾何參數(shù)對(duì)發(fā)動(dòng)機(jī)總體性能的影響，在設(shè)計(jì)狀態(tài)下發(fā)動(dòng)機(jī)沿程各截面計(jì)算特征參數(shù)與試驗(yàn)值最大誤差不超過(guò)5.0%。

（3）對(duì)比分析發(fā)現(xiàn)，風(fēng)扇出口分流板控制涵道比并對(duì)風(fēng)扇性能有一定影響。渦輪出口摻混段整流板對(duì)風(fēng)扇和低壓渦輪性能均有較大影響。通過(guò)調(diào)節(jié)過(guò)渡段幾何參數(shù)，實(shí)現(xiàn)各部件性能匹配，達(dá)到發(fā)動(dòng)機(jī)總體性能要求。

該通流模型能夠快速預(yù)估發(fā)動(dòng)機(jī)總體性能，發(fā)現(xiàn)各部件匹配中存在的問(wèn)題，具有一定的工程應(yīng)用價(jià)值。內(nèi)、外涵摻混段涉及變比熱問(wèn)題，本文采用內(nèi)、外涵流量加權(quán)平均比熱比，噴管出口總溫平均值比試驗(yàn)值偏高3.5%。摻混段計(jì)算有待進(jìn)一步完善，使噴管溫度場(chǎng)與實(shí)際情況更加吻合。

[1]Smith L H.The radial-equilibrium equation of turbomachinery[J].ASME Journal of Engineering for Power，1966，88（1）：1.

[2]Jin H L，Jin D H，Li X J.A time-marching throughflow model and its application in transonic axial compressor [J].Journal of Thermal Sinence，2010，19（6）：519-525.

[3]Denton J D，Horlock J H.A review of some early design practice using computational fluid dynamics and a current perspective[J].Journal of Turbomachinery，2005，127（1）：5-13.

[4]Denton J D，Dawes W N.Computational fluid dynamics for turbomachinery design[J].Archive Proceedings of the Institution of Mechanical Engineers，1998，213（2）：107-124.

[5]Simon J F，Leonard O.Modeling of 3D losses and deviations in a throughflow analysis tool[J].Journal of Thermal Science，2007，16（3）：208-214.

[6]金海良.周向平均方法在多級(jí)軸流風(fēng)扇/壓氣機(jī)設(shè)計(jì)與分析中的應(yīng)用[D].北京：北京航空航天大學(xué)，2011.JIN Hailiang.Application of circumferentional average method in multistage axial fan/compressor design and analysis[D].Beijing：Beijing University of Aeronautics and Astronautics，2011.（in Chinese）

[7]Wan K，Jin H，Jin D.Influence of nonaxisymmetric terms on circumferentially averaged method in fan/compressor[J].Journal of Thermal Science，2013，22（1）：13-22.

[8]Mund F C，Doulgeris G，Pilidis P.Enhanced gas turbine performance simulation using CFD modules in a 2D representation of the low-pressure system for a high-bypass turbofan [J].Gas Turbines Power,2006，129（3）,761-768.

[9]John D Denton.Multall-an open source，CFD based，turbomachinery design system[R].ASME 2017-GT-63993.

[10]Baralon S，Eriksson L E，Hall U.Evaluation of higher-order terms in the throughflow approximation using 3D Navier-Stokes computations of a transonic compressor rotor[R].ASME 1999-GT-074.

[11]Gentry R A，Martin R E，Daly B J.An eulerian differencing method for unsteady compressible flow problems[J].Journal of Computational Physics，1966，1（1）：87-118.

[12]寧方飛.考慮真實(shí)幾何復(fù)雜性的跨音壓氣機(jī)內(nèi)部流動(dòng)的數(shù)值模擬[D].北京：北京航空航天大學(xué)，2002.Ning F F.Numerical investigations of flows in transonic compressors with real geometrical complexities[D].Beijing：Beijing University of Aeronautics and Astonautics，2002.（in Chinese）

[13]Bosman C，Marsh H.An improved method for calculating the flow in turbomachines，including a consistent loss model[J].Archive Journal of Mechanical Engineering Science，1974，31（16）：25-31.

[14]Roe P L.Approximate riemann solvers，parameter vectors，and difference schemes.[J].Journal of Computational Physics，1981，43（2）：357-372.

[15]Strazisar A J，Wood J R，Hathaway M D.Laser anemometer measurements in a transonic axial-flow fan rotor[J].NASA STI/Recon Technical Report N，1981，90（2）：430-437.