朱德軒，楊 旸，余志健，李志鵬，張 冬

(1.中國商飛上海飛機設(shè)計研究院，上海 201210；2.中國科學(xué)院工程熱物理研究所南京未來能源系統(tǒng)研究院，南京 210000；3.中國科學(xué)院工程熱物理研究所先進燃氣輪機實驗室，北京 100190)

民航飛機飛行時燃油箱需要與外界大氣連通以保證油箱結(jié)構(gòu)的安全。同時，飛機惰化系統(tǒng)和環(huán)控系統(tǒng)也需引入外部冷卻空氣以滿足系統(tǒng)運行需求。對于民航飛機來說，目前主要有兩種進氣口設(shè)計，分別為戽斗式進氣口和埋入式進氣口。與戽斗式進氣口相比，埋入式進氣口沒有外部突出物，降低了對飛機氣動阻力的影響。同時無需支撐加固結(jié)構(gòu)，重量較輕。埋入式進氣口主要為(美國)國家航空咨詢委員會(National Advisory Committee for Aeronautics，NACA)埋入式進氣口，該種進氣口已在波音、空客及中國商飛等若干機型上得到廣泛應(yīng)用。

NACA埋入式進氣口，通過側(cè)邊產(chǎn)生一對軸向渦，將壁面氣流卷吸入進氣通道[1]。但由于進氣口通常位于較厚的邊界層內(nèi)，吸入的為邊界層中低能氣流，導(dǎo)致進氣流量減少及出口氣流不均勻度高[2]。同時由于飛機在高馬赫數(shù)(Ma)下飛行，在NACA進氣口前部喉道位置容易產(chǎn)生激波，增加氣動阻力，惡化氣動性能。

因而，需對NACA進氣口進行氣動優(yōu)化設(shè)計。徐尚成等[3]將一體化設(shè)計方法應(yīng)用于進氣道設(shè)計, 采用特征線法設(shè)計了5種不同型面的前體。喬文友等[4]提出一種基于前體激波形狀的一體化設(shè)計方法，以使進氣道捕獲截面和唇口型線與飛行器前體激波匹配。李怡慶等[5]對圓錐流場在不同攻角條件下的氣動特征進行分析, 以流線追蹤技術(shù)為基礎(chǔ)，發(fā)展了一種曲錐前體/三維內(nèi)轉(zhuǎn)進氣道一體化設(shè)計方法, 獲得了3個幾何參數(shù)對一體化方案外形和性能的影響規(guī)律。為設(shè)計低噪聲進氣道，邱昇[6]提出了一套基于聲學(xué)伴隨方法和梯度增強代理模型的高效優(yōu)化設(shè)計框架。李靜等[7]提出了一種三維進氣道沿程結(jié)冰參數(shù)分析方法，并考察了進氣道參數(shù)對結(jié)冰的影響。在實際進氣道工程設(shè)計上，上海飛機設(shè)計研究院薛勇等[8]基于試驗數(shù)據(jù)獲得的進氣口總壓恢復(fù)系數(shù)和截面捕獲流量比關(guān)系，確定飛機極熱天巡航狀態(tài)設(shè)計工況點，而后根據(jù)NACA推薦尺寸比，確定最終尺寸。王赟等[9]提出在滿流量基礎(chǔ)上以燃油代償損失為綜合評價指標(biāo)進行設(shè)計。

國外Soltani等[2]對超音速混合壓縮進氣道進行了實驗研究。以總壓恢復(fù)系數(shù)、質(zhì)量流率和流畸變?yōu)檫M氣性能指標(biāo)。結(jié)果表明，應(yīng)用邊界層抽吸上游的攝入量喉嚨可以大大提高進氣性能的設(shè)計和非設(shè)計工況條件而不影響進氣質(zhì)量流量。Javed等[10]采用計算流體力學(xué)方法計算了超燃沖壓發(fā)動機引擎進氣系統(tǒng)的質(zhì)量捕獲率。

現(xiàn)有進氣道設(shè)計研究主要針對超聲速或發(fā)動機進氣道，考慮減小發(fā)動機前進氣畸變。設(shè)計中兼顧進氣道進氣噪聲及結(jié)冰等特性。對于民用客機公用系統(tǒng)進氣道研究較少。各結(jié)構(gòu)參數(shù)對客機公用系統(tǒng)進氣道性能參數(shù)影響特性研究較少。工程設(shè)計上直接采用基于邊界層的理論公式計算進氣口參數(shù)，各關(guān)鍵幾何尺寸采用NACA早期推薦值，各幾何參數(shù)的具體影響機制和作用機理無法解析。且對于非平行面組成的進氣口，斜坡底部的流動發(fā)散，且外流不平行于斜壁，有渦形成，理論模型精度不足。

為解決上述問題，現(xiàn)采用計算流體力學(xué)(computational fluid dynamics,CFD)方法建立民機公用系統(tǒng)NACA埋入式進氣口數(shù)值模型，對其氣動特性進行分析，并以阻力、總壓恢復(fù)系數(shù)、出口質(zhì)量流量和出口馬赫數(shù)為性能指標(biāo)進行幾何參數(shù)優(yōu)化設(shè)計。首先，理論分析確定幾個關(guān)鍵的進氣口幾何控制參數(shù)，通過試驗設(shè)計(design of experiments，DoE)方法進行參數(shù)優(yōu)化，通過統(tǒng)計回歸找出最優(yōu)的幾何參數(shù)組合結(jié)構(gòu)。并對試驗設(shè)計結(jié)果進行方差分析(analysis of variance, ANOVA)[11]，分析出各控制參數(shù)各自及相互影響規(guī)律。同時基于CFD計算結(jié)果詳細研究各幾何控制參數(shù)對進氣口速度、壓力分布、渦強度及流線等氣動特性影響機理，最終獲得具有最佳氣動特性的進氣口幾何結(jié)構(gòu)參數(shù)。

1 數(shù)值方法

1.1 幾何及邊界條件

所用的NACA進氣道幾何結(jié)構(gòu)如圖1(a)所示，由進氣口、斜坡、轉(zhuǎn)角及尾通道組成。基準工況：出口尺寸68 mm×24 mm，斜坡高30 mm，轉(zhuǎn)角外半徑30 mm，內(nèi)半徑10 mm, 出口高度232 mm。整體計算域如圖1(b)所示，坐標(biāo)原點位于進氣口前緣尖角位置，計算域進口距離前緣尖角1 m，通氣口后1 m結(jié)束，側(cè)面各留0.5 m長，向下擴展0.7 m。氣流經(jīng)過NACA通氣口后拐角轉(zhuǎn)90°向上連接火焰抑制器，計算采用直通口，并預(yù)留一定長度直通道。進氣口附近表面為壁面，主流計算域四周及底面采用壓力遠場邊界，進氣口出口為壓力出口邊界條件。主流參數(shù)如表1所示。

α為斜坡角度；β為斜坡錐角；r為內(nèi)折轉(zhuǎn)半徑；R為折轉(zhuǎn)角半徑

表1 主流參數(shù)

1.2 網(wǎng)格與網(wǎng)格無關(guān)性

主流壁面和NACA進氣口壁面采用3層棱柱層網(wǎng)格，網(wǎng)格第一層高度分別為0.025 mm和0.02 mm，保證無量綱壁面距離y+≈1[12]，由于選用k-omega SST(shear stress transfer)湍流模型，僅研究型面的影響，選擇巡航設(shè)計點作為單一邊界條件，因此對于不同計算工況采用同一套網(wǎng)格分辨率劃分。對整個計算域進行分區(qū)劃分混合網(wǎng)格。對于NACA進氣口及進口向下10 mm單獨分區(qū)采用空間分辨率1.25 mm的結(jié)構(gòu)化多面體網(wǎng)格；進氣口四周一定范圍單獨分區(qū)采用分辨率4 mm的結(jié)構(gòu)化多面體網(wǎng)格；其他區(qū)域采用分辨率15 mm結(jié)構(gòu)化多面體網(wǎng)格。結(jié)構(gòu)化網(wǎng)格以及邊界層網(wǎng)格的交界面處，采用了非結(jié)構(gòu)網(wǎng)格進行填充。由于NACA通氣口為邊界層進氣，非結(jié)構(gòu)填充將極大減小求解過程的數(shù)值耗散問題。整體網(wǎng)格劃分和局部網(wǎng)格細節(jié)如圖2所示。

圖2 整體和局部網(wǎng)格劃分

對網(wǎng)格進行無關(guān)性驗證。在不改變壁面棱柱層網(wǎng)格高度的情況下，對劃分的3個區(qū)域分別進行加密，使得主流、靠近通氣口區(qū)域及通氣口內(nèi)部包含折轉(zhuǎn)部分增加和減少20%的節(jié)點，分別劃分出110萬、197萬和297萬三套網(wǎng)格。

圖3、圖4和圖5為不同網(wǎng)格數(shù)量下進氣口軸向渦發(fā)展截面(x=0.08 m)x方向、y方向及z方向速度分布。x=0.08 m處，邊界層開始有向進氣口內(nèi)部發(fā)展趨勢，兩個軸向渦已經(jīng)發(fā)展完成。將197萬的網(wǎng)格進一步加密到297萬，并不能明顯改變速度分量沿展向的數(shù)值，但加密后計算的速度分布還是有輕微的變化，由于297網(wǎng)格不算太多，因此后續(xù)將采用297萬網(wǎng)格策略。

圖3 不同網(wǎng)格數(shù)量下截面(x=0.08 m)x方向速度分布對比

圖4 不同網(wǎng)格數(shù)量下截面(x=0.08 m)y方向速度分布對比

圖5 不同網(wǎng)格數(shù)量下截面(x=0.08 m)z方向速度分布對比

1.3 數(shù)值方法

穩(wěn)態(tài)RANS(Reynolds-averaged Navier-Stokes)計算采用k-omega SST湍流模型，并采用密度基求解器求解。密度采用理想氣體定律，比熱容采用分段多項式，動力黏度采用Sutherland模型。湍流動能和湍流耗散率采用二階迎風(fēng)格式離散，操作壓力選為零值。壁面均為絕熱條件。

1.4 DoE設(shè)計

根據(jù)三角翼[13]與四面渦發(fā)生器[14]理論，進口坡度對軸向渦產(chǎn)生有重要作用，而對于NACA進氣口需加強軸向渦還是減弱軸向渦，文獻暫未給出結(jié)論?，F(xiàn)有NACA通氣口側(cè)邊根據(jù)文獻[15]設(shè)計，而NACA原始文獻中，該形狀是從三角形出發(fā)，依據(jù)減少進口展向角度修型獲得。對于型面設(shè)計是否影響進氣口氣動性能未有詳細分析。

如圖6所示，NACA進氣口斜坡形狀可以簡化成兩個三角形，并用光滑曲線進行連接[16]。根據(jù)型面曲線結(jié)果，β可保持不變以減少參數(shù)個數(shù)，且民航飛機巡航馬赫數(shù)低于1。

圖6 關(guān)鍵幾何控制參數(shù)

因此主要研究進氣口3個關(guān)鍵幾何參數(shù)，如圖6所示分別為α、θ和折轉(zhuǎn)半徑R。在主流體域及網(wǎng)格策略不變的情況下，固定折轉(zhuǎn)角出口幾何尺寸及其下游幾何域進行試驗設(shè)計。對3個因素進行三水平設(shè)計，確定3個因素對性能的影響是否有非線性；再從中篩選出8個算例，對三因素-兩水平結(jié)果進行方差分析，獲得主要影響參數(shù)。三因素-三水平取值分析如下。

(1)α。該角度為斜坡角度。在保證折轉(zhuǎn)角不變的情況下，增大該角引起進氣口整體縮短。根據(jù)三角翼以及四面體渦發(fā)生器理論，該角度直接影響軸向渦的形態(tài)。主流進氣有3°的攻角，而一般平板邊界層分離擴散角為11°，因此選取11°、13°、和15°共3個幾何變量水平，相應(yīng)的，主流對斜坡的氣動攻角為8°、10°和12°。

(2)θ。該角度為斜坡擴張角，決定斜坡型面拐點位置。θ小則拐點靠近前緣。在NACA報告中沒有對此的描述，而是用光滑樣條曲線連接了出口與前部三角形。為了量化該型面，采用θ來控制拐點位置，取值范圍為36°、40°和44°。

(3)折轉(zhuǎn)角半徑R。在保證出口位置一定時，該弧度會影響進氣口的阻力以及轉(zhuǎn)向時因面積變化導(dǎo)致的流道內(nèi)馬赫數(shù)變化。R的取值為32、40和44 mm。

三因素-三水平的試驗設(shè)計表如表2所示。

表2 水平試驗設(shè)計表

1.5 性能參數(shù)指標(biāo)

選取如下4個性能參數(shù)指標(biāo)作為試驗設(shè)計的因變量。

(1)阻力D。該阻力為作用在NACA進氣口斜坡、斜坡側(cè)面以及所有轉(zhuǎn)角面的總摩擦力與壓力之和，采用絕對值。

(2)總壓恢復(fù)系數(shù)。定義[15]為

(1)

(3)進入進氣口的質(zhì)量流量m，為斜坡出口拐角前截面數(shù)據(jù)。

(4)斜坡出口拐角前截面平均馬赫數(shù)。

2 結(jié)果與分析

2.1 基準工況

圖7為基準工況下進氣口表面極限流線和靜壓分布圖?？梢钥闯觯?/2軸向位置的側(cè)面上有低壓區(qū)，而靠近進氣口折轉(zhuǎn)角位置處靜壓恢復(fù)。

圖7 基準工況進氣口表面極限流線和靜壓分布圖

表面極限流線表明斜坡底部沿流向流線呈現(xiàn)分散趨勢，因此在底部無分離。在側(cè)面靠近通氣口底部的側(cè)邊，通氣口底部流線在此聚集；在側(cè)面靠近主流的一面，側(cè)面流線在此聚集，因此可以推斷，主流在側(cè)邊卷起軸向渦，并向下游發(fā)展。但是這個軸向渦并不貼靠在斜坡底面，否則斜坡底部將會有流線匯聚的現(xiàn)象。

為更形象地分析軸向渦，選取渦量為1.03×108s-1的等值面并用靜壓染色進行分析，如圖8所示。軸向渦從NACA進氣口前緣形成，并且始終靠近斜坡側(cè)面與主流相交的側(cè)邊，而不靠近斜坡的底部，因此推測側(cè)邊的主流與進氣口內(nèi)部壓差，以及側(cè)邊帶來的速度間斷面是產(chǎn)生軸向渦的關(guān)鍵，如果希望增強軸向渦，就需要增大壓差與間斷面，即增大斜坡角并盡量減少側(cè)邊的倒角。側(cè)邊倒角受制造工藝影響，而斜坡角α是可以設(shè)計的。觀察兩支軸向渦，發(fā)現(xiàn)在靠近進氣口處靜壓較低，而在斜坡側(cè)邊拐點處，渦管上靜壓陡升。在側(cè)邊形成的軸向渦繼續(xù)向下游發(fā)展，通過90°轉(zhuǎn)角延伸至出口，而且越向下游，兩支渦管越有相互靠攏的趨勢。在轉(zhuǎn)角附近，速度方向發(fā)生變化，該處可被Q渦識別方法捕獲，并具有較低靜壓。

圖8 渦核區(qū)靜壓染色圖

2.2 DoE結(jié)果

DoE 27個算例求解的性能參數(shù)如表3所示。可以看出，在斜坡角α為15°的所有組合中，進氣口出口馬赫數(shù)均超過1，而在α為11°的所有組合中，出口馬赫數(shù)均小于1，當(dāng)α為13°時，則兩種情況均有。當(dāng)馬赫數(shù)超過1時，跨過激波面總壓損失增大，速度降低，靜壓升高。單就氣動損失而言，是不利的，然而這也是一種短距離內(nèi)快速降低流入進氣口速度的方法。以表3中第2(α=15°)與第27(α=11°)工況為例，在其他參數(shù)不變時，盡管高α角(第2個工況)總壓恢復(fù)系數(shù)較低α角的工況縮小了11.12%，阻力降低了9.62%，但低α角流入通氣口的氣體流量增加了14.57%。

表3 DoE算例計算結(jié)果

對表3數(shù)據(jù)進行三因素-三水平的方差分析，各因素對各性能參數(shù)的主效應(yīng)和交互效應(yīng)如圖9和圖10所示。由圖9可以看出，θ對阻力有非線性的影響，由于整體影響偏小，且上下限相對平均值小于2%，該非線性影響基本可忽略。α對4個因變量影響較其他兩個參數(shù)影響大，隨著α增大，除了出口截面馬赫數(shù)外，其他參數(shù)值均減小。因此減小α對進氣口氣動性能提升有利。θ對4種因變量的影響在3個幾何參數(shù)中最小。除了阻力特性外，其他參數(shù)均隨θ增大而減小。出口轉(zhuǎn)折半徑R對阻力和流量影響大，且均隨著折轉(zhuǎn)半徑R增大而減小。

圖9 各因素對性能參數(shù)主效應(yīng)影響

從圖10交互效應(yīng)可以看出，阻力和質(zhì)量流量上α和θ與折轉(zhuǎn)半徑R相互影響效應(yīng)強；總壓恢復(fù)系數(shù)上，各因素之間相互影響弱；出口馬赫數(shù)上，α和θ之間有相互影響。由主效應(yīng)分析可知，α越小，總壓恢復(fù)系數(shù)和流量越大，進氣口性能越好；而后根據(jù)交互效應(yīng)，在α一定時較小的折轉(zhuǎn)半徑和θ對這兩個參數(shù)有利。表3中工況7和工況27為兩個極端參數(shù)工況，可見工況27性能較佳。

圖10 各因素對性能參數(shù)交互效應(yīng)影響

2.3 工況2 和工況27結(jié)果對比

由圖9可知，α對進氣口氣動特性影響大，且較小α、θ和折轉(zhuǎn)半徑R性能較佳，并進一步比較工況2和工況27性能指標(biāo)，可得α=11°、θ=11°和折轉(zhuǎn)半徑R=36 mm的工況27進氣口綜合性能最好。因而基于工況27(α=11°)和工況2(α=15°)，其他兩參數(shù)保持一致，進一步分析α角的影響機制。

圖11為進氣口渦核區(qū)馬赫數(shù)分布，α越小，進口越大，進氣口前緣氣流擠壓效果越弱，前緣位置馬赫數(shù)越小。

圖11 渦核區(qū)馬赫數(shù)分布

圖12為渦核區(qū)靜壓分布，α越小，進氣口前緣靜壓越大，對減弱壓損有力，斜坡處渦管越長，越有利于卷吸氣流，同時渦管向折轉(zhuǎn)通道下游延伸越長，有利于加強下游氣流摻混，進而利于出口氣流均勻度提升。

圖12 渦核區(qū)靜壓分布

圖13為斜坡各截面馬赫數(shù)分布，可以看出α越小，截面馬赫數(shù)分布越均勻，且數(shù)值越小。圖14 為流入進氣口流線及斜坡各截面速度分布云圖，可以看出α越小，相同位置截面速度分布越均勻。更加均勻的速度分布有利于減小出口畸變，提升進氣口進氣質(zhì)量。從流線上看，α越小，壁面分離效果越弱，有利于減弱壓損和增加出口氣流均勻性。

圖13 斜坡各截面馬赫數(shù)分布

圖14 流入進氣口流線及斜坡各截面速度分布

3 結(jié)論

(1)進氣口軸向側(cè)面上有低壓區(qū)，在斜坡側(cè)邊靠近主流處有一對軸向渦，且越向下游，兩支渦管越靠攏。

(2)α對4個性能參數(shù)影響大，α增大，除出口截面平均馬赫數(shù)外，其他參數(shù)值均減小。θ對性能參數(shù)影響小。出口折轉(zhuǎn)半徑R對平均阻力和平均流量影響大，且均隨著R增大而減小。

(3)交互效應(yīng)可得，阻力和質(zhì)量流量上，α和θ與折轉(zhuǎn)半徑R相互影響效應(yīng)強。α越小，總壓恢復(fù)系數(shù)和流量越大。而后根據(jù)交互效應(yīng)，可確定較小的折轉(zhuǎn)半徑R和θ性能指標(biāo)越好。

(4)α越小，前緣位置馬赫數(shù)越小，靜壓越大，同時渦管向下游延伸越長，同截面馬赫數(shù)和速度分布越均勻，壁面分離越弱，有利于減小壓損，增加卷吸氣流，減小出口畸變，提升進氣質(zhì)量。

綜上分析，較小α、θ和折轉(zhuǎn)半徑R的性能較佳，并進一步比較工況2和工況27性能指標(biāo)，可得α=11°、θ=36°和折轉(zhuǎn)半徑R=36 mm的工況27進氣口綜合性能最好。