周翰瑋，陳勇，*，譚兆光，司江濤，李杰，李棟

1. 中國(guó)商用飛機(jī)有限責(zé)任公司 上海飛機(jī)設(shè)計(jì)研究院，上海 201210 2. 西北工業(yè)大學(xué) 航空學(xué)院，西安 710072

翼身融合(BWB)布局飛機(jī)是一種新型布局的運(yùn)輸機(jī)。由于其高度的集成特性，其浸潤(rùn)面積也大大小于同量級(jí)的傳統(tǒng)布局飛機(jī)，因此，與同量級(jí)傳統(tǒng)布局的飛機(jī)相比，BWB布局飛機(jī)具有較輕的重量、更高的升阻比和更小的燃油消耗量。除上述優(yōu)勢(shì)外，BWB布局飛機(jī)還具有其他更好的特性，比如更好的舒適性等[1]。歐美等航空工業(yè)發(fā)達(dá)地區(qū)相繼投入大量資源進(jìn)行BWB研究[2-3]，美國(guó)以波音公司為首攜同Langley研究中心以及斯坦福大學(xué)等科研機(jī)構(gòu)對(duì)BWB進(jìn)行了廣泛而深入的研究，涉及氣動(dòng)設(shè)計(jì)[4-5]、結(jié)構(gòu)分析與客艙設(shè)計(jì)[6-7]、風(fēng)洞試驗(yàn)[8]、操穩(wěn)特性評(píng)估[9]等。中國(guó)相關(guān)高校和科研機(jī)構(gòu)也針對(duì)BWB布局飛機(jī)進(jìn)行了探索和研究[10-13]，西北工業(yè)大學(xué)張彬乾等已針對(duì)BWB布局總體方案、氣動(dòng)設(shè)計(jì)、飛行控制等多個(gè)領(lǐng)域開展了研究工作[14-15]。

然而由于總體布局、適航、噪聲等方面的限制，BWB布局飛機(jī)不能像傳統(tǒng)布局飛機(jī)采用翼吊或者在機(jī)身上安裝發(fā)動(dòng)機(jī)，因此在背部安裝發(fā)動(dòng)機(jī)成為首選，將發(fā)動(dòng)機(jī)布置在BWB布局飛機(jī)背部，能夠有效減小發(fā)動(dòng)機(jī)短艙對(duì)飛機(jī)翼面流場(chǎng)的干擾和破壞，顯著提升飛機(jī)的氣動(dòng)性能。但背部安裝發(fā)動(dòng)機(jī)容易帶來一系列的空氣動(dòng)力干擾問題，其發(fā)動(dòng)機(jī)與機(jī)體之間的干擾和影響機(jī)制與常規(guī)布局飛機(jī)不同。發(fā)動(dòng)機(jī)與機(jī)身一體化的氣動(dòng)設(shè)計(jì)已成為BWB發(fā)展的關(guān)鍵領(lǐng)域。本文針對(duì)BWB布局飛機(jī)的機(jī)體和發(fā)動(dòng)機(jī)之間的氣動(dòng)干擾進(jìn)行了數(shù)值研究，主要包括噴流短艙及通流短艙的影響，以及短艙高度、展向位置、弦向位置等位置參數(shù)影響。

1 計(jì)算方法概述

1.1 控制方程與離散格式

本文控制方程為三維積分形式的雷諾平均Navier-Stokes方程[16],其表達(dá)式為

(1)

式中：v為控制體體積；s為控制體表面；Q為守恒量；fiv為對(duì)流項(xiàng)；fvs為黏性通量；n為單元表面的外法向單位矢量?？刂品匠滩捎糜邢摅w方法進(jìn)行離散，時(shí)間推進(jìn)采用LU-SGS(Lower Upper-Symmetric Gauss Seidel)方法，空間離散采用Roe’s FDS格式，高階重構(gòu)采用三階MUSCL格式，為了抑制激波振蕩采用van Albada限制器，湍流計(jì)算采用了兩方程k-ωSST(Shear Stress Transport)湍流模式[17]。

1.2 渦輪風(fēng)扇發(fā)動(dòng)機(jī)計(jì)算模型

渦輪風(fēng)扇發(fā)動(dòng)機(jī)流動(dòng)特性非常復(fù)雜。但是由于發(fā)動(dòng)機(jī)對(duì)飛機(jī)部件的干擾主要是體現(xiàn)在進(jìn)排氣效應(yīng)上，而發(fā)動(dòng)機(jī)內(nèi)部流場(chǎng)并不是主要關(guān)注點(diǎn),因此利用特定的計(jì)算邊界條件對(duì)其進(jìn)排氣流場(chǎng)進(jìn)行模擬就能夠?qū)崿F(xiàn)整個(gè)干擾流場(chǎng)的模擬。對(duì)于一個(gè)典型的渦輪風(fēng)扇發(fā)動(dòng)機(jī)，其計(jì)算模型可簡(jiǎn)化為如圖1 所示[18-19]。

圖1 渦扇發(fā)動(dòng)機(jī)簡(jiǎn)化模型[18-19]Fig.1 Simplified model for turbo fan engine[18-19]

1.3 推阻力的定義

對(duì)推阻力進(jìn)行準(zhǔn)確定義是計(jì)算飛機(jī)動(dòng)力對(duì)氣動(dòng)特性影響的重要前提，本節(jié)簡(jiǎn)要介紹采用的推力和阻力分析方法[20]。

流向亞聲速飛機(jī)動(dòng)力裝置的氣流，其流管形狀主要由發(fā)動(dòng)機(jī)工作狀態(tài)決定。通過對(duì)流過動(dòng)力裝置的內(nèi)流捕獲流管和發(fā)動(dòng)機(jī)控制剖面取控制體，可以獲得發(fā)動(dòng)機(jī)凈推力參數(shù)。分析發(fā)動(dòng)機(jī)的安裝推力，還要考慮發(fā)動(dòng)機(jī)推力的安裝損失，這主要包含流入發(fā)動(dòng)機(jī)流管表面上的阻力(即附加阻力)、短艙外表面上的壓差阻力、摩擦阻力以及掛架阻力等因素。

如圖2所示，對(duì)于民用飛機(jī)發(fā)動(dòng)機(jī)的典型工作狀態(tài)，從遠(yuǎn)場(chǎng)(位置0)開始，外部繞過動(dòng)力裝置的氣流逐漸減速到進(jìn)氣道唇口的駐點(diǎn)位置(位置i)，接著氣流從駐點(diǎn)加速直到發(fā)動(dòng)機(jī)短艙直徑的最大位置(位置M)。和機(jī)翼繞流類似，從位置i到位置M 的外罩上作用的氣動(dòng)力合力在軸向的分量成為外罩吸力。由于黏性、分離以及跨聲速工作時(shí)局部激波的影響，吸力和附加阻力不會(huì)全部抵消，其差值稱為溢流阻力。在不考慮發(fā)動(dòng)機(jī)安裝掛架阻力以及干擾阻力的條件下，短艙阻力可以表示為

(2)

式中:p為靜壓；p0為短艙入口處壓力；A為投影面積；τw為短艙壁面摩阻系數(shù)。

發(fā)動(dòng)機(jī)推力可以表示為

T=?A[ρu(V·n)+(p-p0)nx-τx·n]dA-

?A[ρu(V·n)]dA-?A(-τx·n)dA

(3)

圖2 短艙推阻力分析示意圖Fig.2 Schematic diagram of analysis of thrust and drag of nacelle

式中：ρ為密度；V為速度；u為速度在流向上的分量；nx為n在流向上的分量；τx為流向上的剪應(yīng)力。

2 計(jì)算方法驗(yàn)證

為了驗(yàn)證本文所采用的動(dòng)力模型簡(jiǎn)化方法、網(wǎng)格生成策略以及流場(chǎng)計(jì)算方法的正確性，對(duì)單獨(dú)的帶動(dòng)力發(fā)動(dòng)機(jī)模型流場(chǎng)進(jìn)行了數(shù)值模擬，并與試驗(yàn)結(jié)果進(jìn)行比較。

本文中采用的單獨(dú)發(fā)動(dòng)機(jī)為一軸對(duì)稱渦扇發(fā)動(dòng)機(jī)模型，取自日本航空宇宙技術(shù)研究所“NAL-AERO-02-01” TPS (Turbine Powered Simulation)風(fēng)洞試驗(yàn)?zāi)Ｐ蚚21]。文獻(xiàn)[21]給出了該模型的二維半模輪廓線數(shù)據(jù)，利用建模軟件，將該輪廓線繞軸線旋轉(zhuǎn)360°即得到了三維模型，圖3為生成的該模型表面網(wǎng)格[21-22]。

發(fā)動(dòng)機(jī)的動(dòng)力效應(yīng)，通過進(jìn)排氣邊界條件確定。本文在風(fēng)扇進(jìn)氣邊界給定質(zhì)量流量，在內(nèi)涵和外涵出口邊界給定總壓、總溫和流動(dòng)方向。

如表1所示,對(duì)單獨(dú)發(fā)動(dòng)機(jī)高速巡航狀態(tài)下的進(jìn)、排氣流場(chǎng)進(jìn)行了數(shù)值模擬和分析(雷諾數(shù)Re均為1×106，基于發(fā)動(dòng)機(jī)最大直徑)。表1中，Ma∞為來流馬赫數(shù)，α為迎角，MFR為質(zhì)量流量比，F(xiàn)PR為風(fēng)扇出口總壓比，F(xiàn)TR為風(fēng)扇出口總溫比，CPR為內(nèi)涵出口總壓比，CTR為內(nèi)涵出口總溫比。

圖4為兩種狀態(tài)下風(fēng)扇整流罩和渦輪整流罩表面壓力分布計(jì)算值與試驗(yàn)值的對(duì)比，圖中Cp為壓力系數(shù)，x為沿流向的坐標(biāo)。兩種狀態(tài)下的計(jì)算值與試驗(yàn)值吻合得很好，從而可以驗(yàn)證本文對(duì)于單獨(dú)發(fā)動(dòng)機(jī)數(shù)值模擬所采用的網(wǎng)格分塊策略和進(jìn)、排氣邊界處理方法是合適的。

圖3 “NAL-AERO-02-01”模型表面網(wǎng)格[21-22]Fig.3 Surface mesh of “NAL-AERO-02-01” model[21-22]

表1 狀態(tài)參數(shù)Table 1 Condition parameters

圖4 兩種狀態(tài)下表面壓力計(jì)算值與試驗(yàn)值對(duì)比Fig.4 Comparison of surface pressure between computational and experimental results under two conditions

3 BWB全機(jī)流場(chǎng)數(shù)值模擬

3.1 計(jì)算模型和網(wǎng)格

本文針對(duì)一種BWB布局的飛機(jī)進(jìn)行了數(shù)值計(jì)算研究。分別對(duì)該飛機(jī)的兩種不同構(gòu)型(帶動(dòng)力短艙(Power Nacelle,PN)和帶通流短艙(Through Flow Nacelle,TFN))的流場(chǎng)進(jìn)行計(jì)算分析,如圖5所示。

圖5 計(jì)算采用的兩種不同構(gòu)型Fig.5 Two different configurations for computation

本文中的計(jì)算采用多塊結(jié)構(gòu)化網(wǎng)格,如圖6所示。各塊網(wǎng)格可根據(jù)各自區(qū)域和流場(chǎng)的特點(diǎn)，靈活選擇合適的拓?fù)浣Y(jié)構(gòu)，并安排合理的疏密分布，提高了結(jié)構(gòu)網(wǎng)格處理復(fù)雜外形的能力。在近物面區(qū)域使用“O”型網(wǎng)格，可以很好地保證近物面網(wǎng)格的正交性，非常適合于模擬飛行器附面層流動(dòng)。在流場(chǎng)的其他區(qū)域使用相對(duì)簡(jiǎn)單的“H”型網(wǎng)格。

文中對(duì)飛機(jī)-發(fā)動(dòng)機(jī)一體化分析，采用統(tǒng)一的模式通過設(shè)定進(jìn)氣、排氣口邊界條件來模擬發(fā)動(dòng)機(jī)的動(dòng)力影響效應(yīng)。對(duì)風(fēng)扇進(jìn)氣口邊界進(jìn)行簡(jiǎn)化，通過給定質(zhì)量流量、壓力或者速度等，根據(jù)式(4) 可以求得進(jìn)氣口邊界的質(zhì)量流量比MFR，進(jìn)而求得進(jìn)氣口邊界的速度、密度、溫度、靜壓等參數(shù)。

(4)

式中：AHL為來流入口邊界面積；Afan為風(fēng)扇入口面積；ρfan為風(fēng)扇入口處密度；ρ∞為來流密度；qfan為風(fēng)扇入口處動(dòng)壓；q∞為來流動(dòng)壓。

尾噴口邊界分為外涵道噴口和內(nèi)涵道噴口兩部分，內(nèi)外涵噴口處的馬赫數(shù)可采用與風(fēng)扇進(jìn)氣口邊界處相同的方法求出。

圖6 兩種不同構(gòu)型的網(wǎng)格Fig.6 Mesh of two different configurations

3.2 兩種構(gòu)型計(jì)算結(jié)果

針對(duì)上述兩種構(gòu)型在Ma∞=0.8，Re=4.8×107狀態(tài)的全機(jī)氣動(dòng)特性進(jìn)行計(jì)算。圖7為兩種構(gòu)型在Ma∞=0.8，α=2°時(shí)，表面壓力系數(shù)的對(duì)比。圖8為兩種構(gòu)型在Ma∞=0.8時(shí)，升阻極曲線的對(duì)比(圖中CL和CD分別升力系數(shù)和阻力系數(shù))。從對(duì)比結(jié)果看，帶動(dòng)力短艙的全機(jī)阻力系數(shù)在相同升力系數(shù)下比帶通流短艙的要大，在巡航迎角(約2°)附近要大約0.002，這主要是由動(dòng)力短艙和通流短艙在尾噴口的表面壓力的差別造成的，這也說明在進(jìn)行飛機(jī)的阻力CFD預(yù)測(cè)時(shí)要對(duì)這部分阻力格外注意。

圖7 兩種構(gòu)型表面壓力系數(shù)對(duì)比(Ma∞=0.8，α=2°)Fig.7 Comparison of surface pressure coefficient for two configurations (Ma∞=0.8，α=2°)

圖8 兩種構(gòu)型升阻極曲線對(duì)比(Ma∞=0.8)Fig.8 Comparison of lift to drag polar curves for two configurations (Ma∞=0.8)

4 發(fā)動(dòng)機(jī)位置參數(shù)的影響

本文對(duì)發(fā)動(dòng)機(jī)短艙不同高度位置、展向位置、流向位置對(duì)飛機(jī)性能及發(fā)動(dòng)機(jī)推力性能的影響進(jìn)行了研究。此研究是為了發(fā)現(xiàn)相關(guān)規(guī)律性的問題，因此并沒有考慮飛機(jī)的結(jié)構(gòu)、布置等的可實(shí)現(xiàn)性，有些位置變化會(huì)非常大。同時(shí)，為了單純研究發(fā)動(dòng)機(jī)位置變化的影響，將發(fā)動(dòng)機(jī)掛架去掉，并采用通流短艙構(gòu)型。本節(jié)針對(duì)Ma∞=0.8，Re=4.8×107，α=2°狀態(tài)下的全機(jī)氣動(dòng)性能及發(fā)動(dòng)機(jī)推力性能進(jìn)行計(jì)算。

4.1 短艙高度位置對(duì)飛機(jī)性能的影響

不同發(fā)動(dòng)機(jī)高度位置示意如圖9所示。選擇了發(fā)動(dòng)機(jī)軸線距離機(jī)身上表面由低到高的8個(gè)發(fā)動(dòng)機(jī)位置進(jìn)行研究，升力系數(shù)及升阻比隨發(fā)動(dòng)機(jī)高度(h)增加而變化的曲線如圖10所示。從圖中可以看出，隨著發(fā)動(dòng)機(jī)高度的增加，升力系數(shù)顯著上升，升阻比上升。

為了分析發(fā)動(dòng)機(jī)高度變化導(dǎo)致升阻力變化的原因，做出發(fā)動(dòng)機(jī)中間截面壓力分布，如圖11所示。從圖中可以看出，短艙位置較低時(shí)，短艙下部和機(jī)身上表面會(huì)產(chǎn)生相互干擾，出現(xiàn)一個(gè)較大的低壓區(qū)，對(duì)阻力產(chǎn)生影響。隨著短艙位置的升高，低壓區(qū)范圍減小。

圖12是發(fā)動(dòng)機(jī)推力隨短艙高度的變化，圖中橫線是CFD計(jì)算得到的單獨(dú)發(fā)動(dòng)機(jī)的推力值。發(fā)動(dòng)機(jī)推力隨著離機(jī)身高度的增加是先增大再減小。高度較小(0.5 m以內(nèi))時(shí)，由于受氣流分離的影響，推力相對(duì)單獨(dú)發(fā)動(dòng)機(jī)推力值減小，而隨著分離的消失，機(jī)身反而會(huì)產(chǎn)生有利的干擾，使發(fā)動(dòng)機(jī)的推力增大，但高度到一定程度后，有利干擾越來越小，推力就又慢慢減小。這個(gè)現(xiàn)象與常規(guī)布局民機(jī)的現(xiàn)象不太一樣。這一研究也為后續(xù)BWB布局飛發(fā)一體化設(shè)計(jì)提供了一些參考。

圖9 不同短艙高度位置Fig.9 Different height locations of nacelle

圖10 不同短艙高度對(duì)升力系數(shù)、升阻比的影響Fig.10 Influence of different nacelle heights on lift coefficient and lift-to-drag ratio

圖11 不同高度短艙中間截面壓力分布Fig.11 Distribution of pressure on middle slice for different height nacelles

圖12 推力隨短艙高度的變化Fig.12 Variation of thrust with nacelle height

4.2 短艙展向位置對(duì)飛機(jī)性能的影響

選擇發(fā)動(dòng)機(jī)軸線距離飛機(jī)對(duì)稱面由內(nèi)向外3個(gè)位置進(jìn)行研究，如圖13所示。

升力系數(shù)及升阻比隨發(fā)動(dòng)機(jī)展向位置變化的曲線如圖14所示。由圖可知，發(fā)動(dòng)機(jī)展向位置越靠?jī)?nèi)，升力越大，升阻比越?。辉娇客?，升力越小，升阻比越大。

圖13 不同短艙展向位置Fig.13 Different spanwise locations of nacelle

圖14 不同短艙展向位置對(duì)升力系數(shù)、升阻比的影響Fig.14 Influence of different nacelle spanwise locations on lift coefficient and lift-to-drag ratio

4.3 短艙流向位置對(duì)飛機(jī)性能的影響

選擇發(fā)動(dòng)機(jī)入口截面距離機(jī)頭由前到后3個(gè)位置進(jìn)行研究，升力系數(shù)及升阻比隨發(fā)動(dòng)機(jī)展向位置變化的曲線如圖15所示。由圖可知，發(fā)動(dòng)機(jī)展向位置越靠后，升力越大，升阻比也越大；越靠前，升力越小，升阻比也越小。

圖15 不同短艙流向位置對(duì)升力系數(shù)、升阻比的影響Fig.15 Influence of different nacelle streamwise locations on lift coefficient and lift-to-drag ratio

5 結(jié) 論

1) 所采用發(fā)動(dòng)機(jī)簡(jiǎn)化模型和推阻力定義方法能夠準(zhǔn)確地模擬發(fā)動(dòng)機(jī)的動(dòng)力效應(yīng)。

2) 帶動(dòng)力短艙后全機(jī)的阻力較通流短艙會(huì)有明顯增加。

3) 隨著發(fā)動(dòng)機(jī)離機(jī)身高度增大，發(fā)動(dòng)機(jī)推力先增大后減小，飛機(jī)升力增大，升阻比增大。

4) 隨著發(fā)動(dòng)機(jī)展向位置向外移動(dòng)，飛機(jī)升力減小，升阻比增大。

5) 隨著發(fā)動(dòng)機(jī)流向位置向后移動(dòng)，飛機(jī)升力增大，升阻比增大。

致 謝

感謝西北工業(yè)大學(xué)張彬乾教授課題組的支持與幫助，感謝上海飛機(jī)設(shè)計(jì)研究院同事們的指導(dǎo)與建議。