劉庭申，馮立好，王晉軍

北京航空航天大學(xué) 航空科學(xué)與工程學(xué)院，北京 100083

進(jìn)氣道的錘激波載荷由發(fā)動(dòng)機(jī)強(qiáng)喘振引起。發(fā)動(dòng)機(jī)喘振主要由壓氣機(jī)葉片的部分或全部失速所產(chǎn)生，可導(dǎo)致進(jìn)入發(fā)動(dòng)機(jī)的流量突然減少，同時(shí)伴隨著發(fā)動(dòng)機(jī)進(jìn)口壓力的急劇升高。這種壓力的突變會(huì)產(chǎn)生向上游傳播的壓力波，一般稱(chēng)為錘激波。錘激波的峰值壓力遠(yuǎn)遠(yuǎn)大于定常狀態(tài)的壓力值，通常能達(dá)到定常值的兩倍以上。因此盡管錘激波現(xiàn)象并不常見(jiàn)，錘激波載荷卻是飛機(jī)進(jìn)氣道結(jié)構(gòu)強(qiáng)度設(shè)計(jì)的決定性因素之一[1]。如果進(jìn)氣道的結(jié)構(gòu)強(qiáng)度無(wú)法承受錘激波載荷的沖擊，發(fā)動(dòng)機(jī)出現(xiàn)喘振時(shí)可能會(huì)導(dǎo)致進(jìn)氣道的破壞，進(jìn)而威脅飛機(jī)的飛行安全。如果進(jìn)氣道結(jié)構(gòu)強(qiáng)度設(shè)計(jì)過(guò)于保守，則會(huì)使得進(jìn)氣道重量太大影響到飛機(jī)的性能。所以，對(duì)進(jìn)氣道錘激波現(xiàn)象進(jìn)行研究，獲得較為準(zhǔn)確的進(jìn)氣道錘激波載荷具有非常重要的價(jià)值。

國(guó)外已經(jīng)進(jìn)行了一些關(guān)于進(jìn)氣道錘激波現(xiàn)象的研究。一些研究者采用試驗(yàn)方法探究了飛機(jī)進(jìn)氣道的錘激波現(xiàn)象[1-6]，通過(guò)各種方法誘使發(fā)動(dòng)機(jī)喘振，進(jìn)而進(jìn)行數(shù)據(jù)的測(cè)量。Marshall[2]依據(jù)大量的試驗(yàn)測(cè)量數(shù)據(jù)總結(jié)了發(fā)動(dòng)機(jī)喘振時(shí)壓力變化過(guò)程以及喘振超壓的估算方法，可以在只做較少次數(shù)的發(fā)動(dòng)機(jī)喘振試驗(yàn)的情況下估算多種不同狀態(tài)的發(fā)動(dòng)機(jī)喘振壓升情況，為進(jìn)氣道的錘激波研究提供了基礎(chǔ)。Nugent和Holzman[5]對(duì)F-111A飛機(jī)進(jìn)行飛行試驗(yàn)，獲得了多種飛行狀態(tài)下發(fā)生喘振時(shí)發(fā)動(dòng)機(jī)和進(jìn)氣道壓力隨時(shí)間變化情況以及錘激波的傳播速度等數(shù)據(jù)。在數(shù)值計(jì)算方面，早期的研究主要針對(duì)二維管道開(kāi)展[7-8]。Hindash等[7]通過(guò)對(duì)二維管道的錘激波模擬分析了進(jìn)氣道的受載情況并與試驗(yàn)數(shù)據(jù)對(duì)比驗(yàn)證了計(jì)算的準(zhǔn)確性，同時(shí)文章進(jìn)一步分析了進(jìn)氣道出口不均勻壓升產(chǎn)生的錘激波的傳播特點(diǎn)以及不同進(jìn)氣道幾何形狀對(duì)錘激波傳播過(guò)程的影響。

此后一些學(xué)者陸續(xù)進(jìn)行了三維數(shù)值模擬研究，并且應(yīng)用到了多種戰(zhàn)斗機(jī)進(jìn)氣道錘激波現(xiàn)象的研究中[9-12]。Goble等[9]進(jìn)行了F-22飛機(jī)進(jìn)氣道的錘激波三維數(shù)值模擬研究，他們采用在進(jìn)氣道出口設(shè)定隨時(shí)間變化的壓力作為邊界條件來(lái)產(chǎn)生錘激波，計(jì)算了多種飛行條件下的錘激波情況，總結(jié)了進(jìn)氣道所受的壓力分布規(guī)律，并發(fā)現(xiàn)了進(jìn)氣道的支路和抽吸孔排氣可以降低進(jìn)氣道所受的錘激波峰值壓力。類(lèi)似的錘激波模擬方法應(yīng)用到了瑞典JAS-39飛機(jī)的進(jìn)氣道錘激波研究中[10]，作者討論了動(dòng)態(tài)計(jì)算時(shí)間步長(zhǎng)的選取，以及加快計(jì)算速度的方法，探討了不同計(jì)算模型在計(jì)算速度和準(zhǔn)確性方面的優(yōu)劣。根據(jù)計(jì)算結(jié)果獲得了進(jìn)氣道壁面的錘激波載荷數(shù)據(jù)，并為飛機(jī)的進(jìn)氣道強(qiáng)度設(shè)計(jì)提供了依據(jù)。數(shù)值方法的應(yīng)用可以較為準(zhǔn)確地獲得進(jìn)氣道所受錘激波壓力的分布情況，為進(jìn)氣道結(jié)構(gòu)強(qiáng)度設(shè)計(jì)提供更豐富的信息。Gridley等[13]進(jìn)一步提出了將不同飛行狀態(tài)的錘激波載荷和各狀態(tài)出現(xiàn)的概率相結(jié)合的方法，允許在小概率情況下出現(xiàn)錘激波載荷超出進(jìn)氣道結(jié)構(gòu)設(shè)計(jì)強(qiáng)度，從而進(jìn)一步降低進(jìn)氣道的結(jié)構(gòu)重量。

國(guó)內(nèi)關(guān)于進(jìn)氣道錘激波載荷的相關(guān)文獻(xiàn)較少。主要由朱宇和沈天榮[14]歸納總結(jié)了國(guó)外錘激波載荷的研究方法和成果，而關(guān)于試驗(yàn)和數(shù)值計(jì)算的研究成果很少。目前，國(guó)內(nèi)主要采用工程估算的方法評(píng)估發(fā)動(dòng)機(jī)喘振超壓比，以此作為進(jìn)氣道錘激波載荷。與此同時(shí)，如果能輔以數(shù)值模擬計(jì)算錘激波沿管道運(yùn)動(dòng)過(guò)程中載荷變化情況，可以更加有效地為進(jìn)氣道的結(jié)構(gòu)設(shè)計(jì)提供幫助。

雖然國(guó)外已經(jīng)對(duì)錘激波現(xiàn)象進(jìn)行了一些研究，但是大部分只提供了進(jìn)氣道壁面壓力分布，較少關(guān)注錘激波過(guò)程中進(jìn)氣道內(nèi)部的流動(dòng)結(jié)構(gòu)變化，尤其是內(nèi)部復(fù)雜的三維流動(dòng)特性。本文對(duì)進(jìn)氣道錘激波進(jìn)行研究，探討適用于進(jìn)氣道錘激波計(jì)算的三維非定常數(shù)值模擬方法，獲得錘激波傳播過(guò)程中進(jìn)氣道壁面壓力變化情況，并且對(duì)進(jìn)氣道內(nèi)部的三維流動(dòng)結(jié)構(gòu)和動(dòng)態(tài)演化過(guò)程進(jìn)行分析總結(jié)。

1 模型和計(jì)算方法

本文采用的計(jì)算模型為M2129型S形進(jìn)氣道。這是一種圓形截面進(jìn)氣道。它由AGARD(Advisory Group for Aerospace Research and Development)設(shè)計(jì)，并進(jìn)行了風(fēng)洞試驗(yàn)[15]。試驗(yàn)在來(lái)流馬赫數(shù)為0.21下進(jìn)行，包含了兩種不同的流量條件，兩種流量下喉道的馬赫數(shù)分別為0.412和0.794。試驗(yàn)所測(cè)得的進(jìn)氣道壁面上的壓力分布數(shù)據(jù)可以獲得，經(jīng)常被用于數(shù)值模擬結(jié)果的驗(yàn)證對(duì)比[16-18]。

M2129型S形進(jìn)氣道的剖面幾何形狀見(jiàn)圖1。該進(jìn)氣道模型包括外部壁面、筆直的入口段和出口段，以及中間的S段。為了接近常見(jiàn)飛行器的進(jìn)氣道尺寸，本文中的進(jìn)氣道模型尺寸相比于原試驗(yàn)中的模型擴(kuò)大了4.65倍。其入口段的直徑為0.6 m，出口段直徑約為0.7 m，入口段長(zhǎng)度為0.6 m，中間S段的中軸線(xiàn)長(zhǎng)度為2.184 m。由于模型是對(duì)稱(chēng)的，僅選取了一半的結(jié)構(gòu)進(jìn)行計(jì)算。為了記錄進(jìn)氣道壁面壓力隨時(shí)間的動(dòng)態(tài)變化，選取了上下壁面各11個(gè)點(diǎn)進(jìn)行記錄，即P1～P11和S1～S11。22個(gè)點(diǎn)的編號(hào)和位置也在圖1中標(biāo)出，第1個(gè)點(diǎn)同進(jìn)氣道入口的距離和后面每相鄰兩個(gè)點(diǎn)沿x方向的距離均為0.3 m。

圖1 M2129型S形進(jìn)氣道幾何模型Fig.1 M2129 S-shaped inlet geometry model

圖2 進(jìn)氣道壁面網(wǎng)格劃分Fig.2 Grid distribution on inlet wall

整個(gè)計(jì)算域使用結(jié)構(gòu)網(wǎng)格進(jìn)行劃分，進(jìn)氣道壁面上網(wǎng)格如圖2所示。為了精確模擬進(jìn)氣道內(nèi)部的流動(dòng)，進(jìn)氣道內(nèi)部流場(chǎng)采用了較密的網(wǎng)格。近壁面的網(wǎng)格也進(jìn)行了加密來(lái)模擬壁面流動(dòng)。網(wǎng)格數(shù)量一共約200萬(wàn)。

本文采用非定常雷諾平均Navier-Stokes方程進(jìn)行數(shù)值計(jì)算，采用耦合求解器求解。外流場(chǎng)設(shè)定為壓力遠(yuǎn)場(chǎng)邊界條件，給定來(lái)流的馬赫數(shù)和靜壓，在進(jìn)氣道出口給定反壓。首先采用定常計(jì)算的方法獲得初值，再在進(jìn)氣道出口給定壓力隨時(shí)間變化的函數(shù)來(lái)模擬發(fā)動(dòng)機(jī)喘振的壓力變化，進(jìn)而進(jìn)行錘激波的動(dòng)態(tài)計(jì)算。這種錘激波的數(shù)值模擬方法將在后文中詳細(xì)介紹。動(dòng)態(tài)計(jì)算的時(shí)間步長(zhǎng)為2.5×10-5s。計(jì)算所選用的湍流模型為k-ε模型。這種湍流模型已經(jīng)由Nichols[16]在M2129型進(jìn)氣道模型上進(jìn)行了靜態(tài)的計(jì)算，并與試驗(yàn)結(jié)果對(duì)比驗(yàn)證了正確性。

本文為了驗(yàn)證計(jì)算的準(zhǔn)確性，將風(fēng)洞試驗(yàn)數(shù)據(jù)和計(jì)算所獲得的結(jié)果進(jìn)行了對(duì)比。試驗(yàn)數(shù)據(jù)來(lái)自文獻(xiàn)[15]。因?yàn)樵囼?yàn)的低流量狀態(tài)和后面用于錘激波計(jì)算的流量狀況比較接近，選取該狀態(tài)的試驗(yàn)數(shù)據(jù)進(jìn)行對(duì)比驗(yàn)證。

試驗(yàn)條件是來(lái)流馬赫數(shù)為0.21，進(jìn)氣道喉道處的馬赫數(shù)為0.412，雷諾數(shù)為1.158×106。圖3為試驗(yàn)和計(jì)算獲得的進(jìn)氣道上下壁面的壓力分布結(jié)果，圖中橫坐標(biāo)為各點(diǎn)x坐標(biāo)與該方向總長(zhǎng)度L的比值，縱坐標(biāo)為靜壓p與來(lái)流總壓pt的比值。上壁面壓力的計(jì)算結(jié)果與試驗(yàn)數(shù)據(jù)非常吻合，下壁面計(jì)算結(jié)果數(shù)值略高于試驗(yàn)數(shù)據(jù)，但是誤差在3%以?xún)?nèi)。表明本文所用方法可以對(duì)M2129型進(jìn)氣道內(nèi)部的流動(dòng)情況進(jìn)行準(zhǔn)確的模擬。

圖3 進(jìn)氣道上下壁面壓力計(jì)算與試驗(yàn)結(jié)果對(duì)比Fig.3 Comparison between numerical and experimental results of inlet upper and lower walls

為進(jìn)一步驗(yàn)證非定常計(jì)算的準(zhǔn)確性，對(duì)激波管問(wèn)題進(jìn)行了數(shù)值模擬，并與理論解進(jìn)行了比較，計(jì)算所用激波管長(zhǎng)度為1 m，均分為等長(zhǎng)的高壓段和低壓段，之間由薄片隔開(kāi)。高壓段壓力為1 Pa，密度為1 kg/m3。低壓段壓力為0.1 Pa，密度為0.125 kg/m3，計(jì)算從兩段氣體之間的薄片突然撤離開(kāi)始。圖4為時(shí)間t=0.2 s時(shí)激波管內(nèi)的壓力p、密度ρ、速度v分布情況的計(jì)算結(jié)果和理論值的對(duì)比。結(jié)果表明本文所用數(shù)值方法可以準(zhǔn)確地模擬激波傳播的非定常過(guò)程。

圖4 激波管內(nèi)壓力、密度和速度的數(shù)值與理論結(jié)果對(duì)比Fig.4 Comparison of pressure,density and velocity between numerical and theoretical results of shock tube

2 錘激波波形設(shè)置

通常所采用的計(jì)算方法是直接給定喘振時(shí)的壓力變化波形作為邊界條件[7-11]，進(jìn)一步進(jìn)行錘激波的計(jì)算，本文采用同樣的方法。壓力波形的確定依賴(lài)于發(fā)動(dòng)機(jī)喘振的實(shí)際測(cè)量數(shù)據(jù)，不同發(fā)動(dòng)機(jī)發(fā)生喘振時(shí)可能會(huì)有不同的壓力變化波形。波形的差異會(huì)影響錘激波的計(jì)算結(jié)果，因此應(yīng)當(dāng)盡量準(zhǔn)確地模擬實(shí)際的壓力變化情況。

通過(guò)對(duì)喘振時(shí)發(fā)動(dòng)機(jī)壓力波形的試驗(yàn)測(cè)量發(fā)現(xiàn)不同的發(fā)動(dòng)機(jī)型號(hào)以及不同的喘振狀況產(chǎn)生的壓力變化波形會(huì)有所差異。但是其隨時(shí)間變化規(guī)律通常為迅速達(dá)到峰值之后緩慢地下降[7]。因此本文選用圖5所示三角形波形進(jìn)行錘激波的數(shù)值模擬，其超壓比為2.1，即峰值壓力為定常值的2.1倍。出口壓力一開(kāi)始保持恒定，在t=0.01 s時(shí)瞬間達(dá)到峰值。之后壓力呈線(xiàn)性下降，經(jīng)過(guò)0.1 s，即在t=0.11 s時(shí)壓力恢復(fù)到初始值。在此之后，壓力保持不變，直到計(jì)算結(jié)束。這種波形的壓力曲線(xiàn)如圖5所示，圖中縱坐標(biāo)為各時(shí)刻壓力同定常時(shí)壓力pss的比值。

圖5 錘激波壓力波形Fig.5 Pressure waveform of hammershock

3 壓力變化與流動(dòng)結(jié)構(gòu)

本節(jié)進(jìn)行錘激波計(jì)算的狀態(tài)為來(lái)流馬赫數(shù)Ma∞=0.65，來(lái)流靜壓p∞=101 325 Pa，靜溫T∞=288.15 K，進(jìn)氣道定常時(shí)的換算流量Wc=60 kg/s(即進(jìn)氣道出口馬赫數(shù)Mae=0.44，入口馬赫數(shù)Ma0=0.63)。進(jìn)氣道上下壁面各點(diǎn)靜壓隨時(shí)間變化記錄在圖6中。各點(diǎn)所對(duì)應(yīng)的位置已經(jīng)在圖1中顯示，圖中將各點(diǎn)的靜壓值除以來(lái)流靜壓p∞從而進(jìn)行無(wú)量綱處理。圖6(a)為上壁面各點(diǎn)的結(jié)果，圖6(b)為下壁面各點(diǎn)的結(jié)果。從圖中可以明顯地看出沿著進(jìn)氣道入口到出口各點(diǎn)壓力曲線(xiàn)的變化規(guī)律。

圖6 進(jìn)氣道上下壁面壓力隨時(shí)間變化曲線(xiàn)Fig.6 Time history of pressure on inlet upper and lower walls

圖6(a)中沿實(shí)線(xiàn)箭頭方向依編號(hào)次序逐個(gè)經(jīng)過(guò)P1～P11的曲線(xiàn)，說(shuō)明越靠近進(jìn)氣道出口的點(diǎn)壓力開(kāi)始升高的時(shí)間越早，反映了錘激波沿著進(jìn)氣道逐漸向上游傳播的過(guò)程。觀察圖6(a)中壓力曲線(xiàn)在達(dá)到峰值之后的下降段，從P1～P11各點(diǎn)的曲線(xiàn)也沿著虛線(xiàn)箭頭方向排布，表明越靠近進(jìn)氣道出口的點(diǎn)，壓力達(dá)到峰值之后下降得越慢。由于該種變化規(guī)律，使得曲線(xiàn)的形態(tài)呈現(xiàn)出越靠近進(jìn)氣道入口的點(diǎn)的曲線(xiàn)顯得更尖，而越接近進(jìn)氣道出口的點(diǎn)的曲線(xiàn)更寬。關(guān)于下壁面的點(diǎn)S1～S11，圖中箭頭也標(biāo)示出了相似的變化規(guī)律。根據(jù)最靠近進(jìn)氣道入口的點(diǎn)P1和S1的壓力達(dá)到峰值的時(shí)間，可計(jì)算出錘激波在進(jìn)氣道內(nèi)移動(dòng)的平均速度為306 m/s，馬赫數(shù)約為0.9。而錘激波相對(duì)于上游來(lái)流的平均傳播速度為481 m/s，馬赫數(shù)達(dá)1.42。

圖7統(tǒng)計(jì)了上下壁面各11個(gè)點(diǎn)的峰值壓力和定常時(shí)的壓力值。可以看到總體上靠近入口位置的載荷略大于出口附近載荷，最大峰值壓力位于上壁面靠近入口的彎折處。上下壁面各點(diǎn)的峰值壓力均達(dá)到來(lái)流靜壓的2.3～2.65倍左右，說(shuō)明錘激波能夠誘導(dǎo)產(chǎn)生比定常狀態(tài)更高的載荷。

圖7 進(jìn)氣道上下壁面峰值壓力和定常狀態(tài)壓力Fig.7 Peak pressure and steady pressure on the inlet upper and lower walls

圖8顯示了錘激波經(jīng)過(guò)進(jìn)氣道時(shí)壁面上的壓力同來(lái)流靜壓比值(p/p∞)的分布云圖；圖9為上下壁面沿程壓力分布曲線(xiàn)，從圖中可以看到錘激波沿進(jìn)氣道向上游傳播時(shí)壓力變化的過(guò)程。圖8(a)為t=0.01 s時(shí)的壓力云圖，此時(shí)進(jìn)氣道出口壓力剛產(chǎn)生突變，在進(jìn)氣道壁面上還未產(chǎn)生明顯的影響。在t=0.015 s時(shí)，即圖8(b)所示，錘激波剛經(jīng)過(guò)進(jìn)氣道的第1個(gè)彎折處，此時(shí)錘激波位置的下壁面受到的壓力略高于上壁面。在t=0.017 5 s和t=0.02 s兩個(gè)時(shí)刻，即圖8(c)和圖8(d)所示，錘激波正經(jīng)過(guò)進(jìn)氣道第2個(gè)彎道，此時(shí)可以明顯看出錘激波處上壁面的壓力較高。由圖8所示的壓力云圖可以發(fā)現(xiàn)，錘激波向上游傳播的過(guò)程中進(jìn)氣道受到的最大壓力的區(qū)域位于進(jìn)氣道入口附近彎折處的上壁面，這一區(qū)域所受峰值壓力達(dá)到來(lái)流靜壓的2.65倍左右。在錘激波傳出進(jìn)氣道之后，進(jìn)氣道出現(xiàn)了中部壓力較高，入口和出口附近壓力較低的特點(diǎn)，如圖8(e)所示。在這一中部高壓區(qū)逐漸后移并消散之后，進(jìn)氣道便呈現(xiàn)出越靠近入口壓力越低的特點(diǎn)。

圖8 進(jìn)氣道壁面壓力云圖Fig.8 Pressure contours on inlet wall

圖9 進(jìn)氣道上下壁面沿程壓力分布Fig.9 Pressure distribution along inlet upper and lower walls

對(duì)于圖9所示壓力沿程分布曲線(xiàn)。在t=0.01 s時(shí)剛剛產(chǎn)生錘激波，因而可以看到相應(yīng)時(shí)刻曲線(xiàn)在出口處(x=3.4 m)壓力突增，而其他部位的壓力分布與穩(wěn)態(tài)時(shí)相同。在t=0.02 s時(shí)，錘激波傳至進(jìn)氣道入口附近?？梢钥吹剑诩げ嫦掠危舷卤诿鎵毫ρ剡M(jìn)氣道軸向以相近的梯度逐漸下降，與所設(shè)定的錘激波波形的形態(tài)類(lèi)似。錘激波傳出進(jìn)氣道之后，壁面壓力整體上開(kāi)始逐漸下降。在t=0.03 s時(shí)，進(jìn)氣道壁面壓力呈現(xiàn)入口處低、出口處高的特點(diǎn)，并且上下壁面沿軸向壓力以近似的梯度增高。此時(shí)入口處壓力已經(jīng)降至接近當(dāng)?shù)胤€(wěn)態(tài)時(shí)壓力。而在之后時(shí)刻，入口附近的壓力會(huì)進(jìn)一步下降直至低于當(dāng)?shù)胤€(wěn)態(tài)時(shí)的壓力，典型的如圖中展示的t=0.08 s時(shí)刻的曲線(xiàn)。在該時(shí)刻，進(jìn)氣道前段所受壓力低于其穩(wěn)態(tài)時(shí)壓力，而進(jìn)氣道后段所受的壓力仍高于穩(wěn)態(tài)壓力。值得注意的是此時(shí)進(jìn)氣道在靠近入口彎道處的下壁面區(qū)域。由于進(jìn)氣道壁面在此處突然開(kāi)始彎曲，在穩(wěn)態(tài)時(shí)，該處便是一個(gè)明顯的低壓區(qū)域。但在t=0.08 s時(shí)刻，該區(qū)域的壓力進(jìn)一步降低，同時(shí)其壓力梯度也更大。這會(huì)使得進(jìn)氣道在該處受到較大的吸力載荷以及較大的剪切應(yīng)力。

不同時(shí)刻進(jìn)氣道對(duì)稱(chēng)面上壓力云圖以及流線(xiàn)圖顯示在圖10中，圖中壓力均除以來(lái)流靜壓以進(jìn)行無(wú)量綱處理。圖10(a)～圖10(d)顯示了錘激波在進(jìn)氣道內(nèi)的傳播過(guò)程，激波面的法向基本與進(jìn)氣道中軸線(xiàn)方向一致。位于錘激波上游的壓力和流線(xiàn)未受到影響，依舊與錘激波發(fā)生前一樣，但是在錘激波下游，則可以觀察到回流的產(chǎn)生。值得注意的是，圖10(c)錘激波經(jīng)過(guò)進(jìn)氣道靠近入口的彎折處時(shí)，由于進(jìn)氣道的彎曲，導(dǎo)致錘激波下游回流的流線(xiàn)向上壁面傾斜，錘激波處上壁面出現(xiàn)明顯的壓力集中，使得此處成為進(jìn)氣道內(nèi)部峰值壓力最高的部位。在錘激波傳播過(guò)程中進(jìn)氣道下游產(chǎn)生回流，而上游氣流仍從進(jìn)氣道入口流入，因而在進(jìn)氣道內(nèi)部堆積大量流體，形成高壓區(qū)域。圖10(d)中錘激波即將從進(jìn)氣道傳出，激波面后方出現(xiàn)的一段高壓區(qū)域，使得流動(dòng)出現(xiàn)分離，一部分繼續(xù)向上游流動(dòng)，另一部分則流向下游，并在進(jìn)氣道中部誘導(dǎo)產(chǎn)生復(fù)雜的螺旋結(jié)構(gòu)。在錘激波離開(kāi)進(jìn)氣道后，聚集在進(jìn)氣道中部的高壓區(qū)導(dǎo)致氣流同時(shí)向進(jìn)氣道入口和出口流出。由于在進(jìn)氣道出口所設(shè)定的波形持續(xù)時(shí)間較長(zhǎng)，因而入口附近壓力下降更快，因此這一高壓區(qū)域逐漸向后移動(dòng)。圖10(e)中高壓區(qū)位于進(jìn)氣道中部，流動(dòng)在中部分離分別向入口和出口流出，而圖10(f)由于高壓區(qū)的后移，流動(dòng)分離區(qū)域也移動(dòng)到了進(jìn)氣道的后段。中部高壓區(qū)持續(xù)了較短時(shí)間，在圖10(f)之后不久便完全消散，此時(shí)距離錘激波的產(chǎn)生時(shí)間僅過(guò)了0.02 s，而此時(shí)進(jìn)氣道出口的壓力還處于較高的狀態(tài)，進(jìn)而導(dǎo)致進(jìn)氣道再次出現(xiàn)回流，回流從進(jìn)氣道入口流出，在進(jìn)氣道外部產(chǎn)生旋渦，如圖10(g)所示。此時(shí)呈現(xiàn)出越靠近進(jìn)氣道入口，壓力下降越快的特點(diǎn)?；亓鞒掷m(xù)了較長(zhǎng)的時(shí)間，并隨著進(jìn)氣道出口壓力的下降而逐漸減弱。在進(jìn)氣道出口壓力恢復(fù)到初始?jí)毫r(shí)，進(jìn)氣道內(nèi)流動(dòng)尚未恢復(fù)正常，而是直到0.195 s以后，即距錘激波發(fā)生0.185 s之后，流動(dòng)恢復(fù)到初始方向并逐漸穩(wěn)定，如圖10(h)所示。

圖10 不同時(shí)刻進(jìn)氣道對(duì)稱(chēng)面流線(xiàn)圖和壓力云圖Fig.10 Streamlines and pressure contours on symmetry plane of inlet at different time

圖11 不同時(shí)刻進(jìn)氣道內(nèi)部三維流線(xiàn)圖 (左：側(cè)視；右：俯視)Fig.11 3D streamlines in inlet duct at different time (left：side view; right：top view)

為了進(jìn)一步分析進(jìn)氣道內(nèi)部的三維流動(dòng)特征，圖11給出了不同時(shí)刻進(jìn)氣道內(nèi)部的三維流線(xiàn)圖，圖中每一時(shí)刻均由兩種不同視角流線(xiàn)圖組成。圖中流線(xiàn)顏色表示速度的大小，單位為m/s。t=0.01 s時(shí)的流線(xiàn)圖如圖11(a)所示，此時(shí)進(jìn)氣道出口剛產(chǎn)生壓力突變，內(nèi)部流動(dòng)還未受到影響。在t=0.015 s時(shí)(圖11(b))，錘激波傳播至進(jìn)氣道中部，可以看到錘激波下游的回流速度遠(yuǎn)小于錘激波上游氣流的流速。從俯視圖中可以看到此時(shí)進(jìn)氣道內(nèi)部流動(dòng)沒(méi)有出現(xiàn)明顯的側(cè)向運(yùn)動(dòng)，流線(xiàn)基本與進(jìn)氣道中軸線(xiàn)方向一致。t=0.02 s時(shí)(圖11(c))，位于進(jìn)氣道中心的旋渦使得流線(xiàn)旋轉(zhuǎn)并互相纏繞著向進(jìn)氣道兩側(cè)移動(dòng)。同時(shí)值得注意的是，在錘激波下游的回流中，處于進(jìn)氣道中心旋渦處的氣流流速非常低(小于5 m/s)，而旋渦周邊靠近壁面的氣流流速相對(duì)較高，尤其是下壁面附近，流速為錘激波下游回流中最高。t=0.025 s時(shí)(圖11(d))，位于進(jìn)氣道中部的旋渦可以在圖10(e)中觀察到其在對(duì)稱(chēng)面的二維投影。而在進(jìn)氣道后段的上壁面附近還可以觀察到另外兩個(gè)對(duì)稱(chēng)的旋渦。t=0.03 s時(shí)(圖11(e))，氣流從高壓區(qū)同時(shí)向進(jìn)氣道入口和出口流出，在進(jìn)氣道內(nèi)越遠(yuǎn)離高壓區(qū)的地方流速越快。圖11(f)為t=0.08 s中部高壓區(qū)消失后出現(xiàn)的回流，此時(shí)回流的流速遠(yuǎn)高于錘激波經(jīng)過(guò)進(jìn)氣道時(shí)其下游回流的流速，流線(xiàn)基本與進(jìn)氣道中軸線(xiàn)平行，沒(méi)有出現(xiàn)明顯的三維流動(dòng)特性。

4 進(jìn)氣道換算流量的影響

進(jìn)一步選用3種不同的流量條件進(jìn)行錘激波過(guò)程的數(shù)值模擬。3種計(jì)算條件下穩(wěn)態(tài)時(shí)的換算流量分別為30、45、60 kg/s(對(duì)應(yīng)的進(jìn)氣道出口馬赫數(shù)分別為Mae=0.21，0.31，0.44)。通過(guò)改變進(jìn)氣道出口的反壓來(lái)達(dá)到預(yù)計(jì)的流量條件，而來(lái)流條件以及非定常計(jì)算時(shí)錘激波波形均與第3節(jié)所述相同。

為了對(duì)比換算流量變化的影響，圖12詳細(xì)列出了上下壁面各點(diǎn)在3種不同流量條件下的峰值壓力?？梢钥闯?，隨著換算流量(或進(jìn)氣道出口馬赫數(shù))的減小，進(jìn)氣道各處所受的峰值壓力逐漸增大。

圖12 不同流量條件時(shí)上下壁面各點(diǎn)的峰值壓力Fig.12 Peak pressure of each point on upper and lower walls for different flow rate conditions

圖13 不同流量條件時(shí)點(diǎn)P1處壓力隨時(shí)間變化曲線(xiàn)Fig.13 Time history of pressure at probe point P1 for different flow rate conditions

圖13給出了點(diǎn)P1在3種流量條件下壓力隨時(shí)間的變化過(guò)程，其中點(diǎn)P1的位置已經(jīng)在圖1中標(biāo)出?？梢钥闯鲈诓煌瑩Q算流量時(shí)，除了壓力的數(shù)值呈現(xiàn)差異，壓力開(kāi)始變化的時(shí)間也有所不同，說(shuō)明錘激波的傳播速度不同。依據(jù)進(jìn)氣道入口附近點(diǎn)達(dá)到峰值的時(shí)間，在Wc=30 kg/s(出口馬赫數(shù)Mae=0.21)時(shí)，錘激波傳播速度達(dá)391 m/s(馬赫數(shù)Ma=1.15)，而在Wc=60 kg/s(出口馬赫數(shù)Mae=0.44)時(shí)，錘激波相對(duì)進(jìn)氣道的平均傳播速度為306 m/s(Ma=0.9)。

5 超壓比的影響

錘激波波形的超壓比(Over Pressure Ratio, OPR)表示波形壓力變化峰值與初始?jí)毫Φ谋戎?。為分析超壓比的差異?duì)進(jìn)氣道所受壓力的影響，將錘激波波形變更為峰值壓力達(dá)3倍初始?jí)毫Γ鴫毫ν蛔優(yōu)榉逯档臅r(shí)間以及壓力恢復(fù)到初始值的時(shí)間均與圖5所示的原波形相同。除了波形的改變外，計(jì)算時(shí)其他條件均保持相同。

兩種超壓比條件下上下壁面峰值壓力如圖14所示。表明增大波形的超壓比顯著增大了進(jìn)氣道壁面所受的峰值壓力。超壓比由2.1增大到3使得設(shè)定的波形的峰值壓力增大了約1.43倍。而進(jìn)氣道各點(diǎn)的峰值壓力的增大比例為1.34～1.42倍，略低于波形中峰值壓力的增大比例。

圖15給出了進(jìn)氣道上壁面點(diǎn)P1的壓力變化曲線(xiàn)，可以看出在超壓比為3的條件下，點(diǎn)P1達(dá)到峰值的時(shí)間更早，說(shuō)明錘激波在進(jìn)氣道內(nèi)傳播速度更快。在超壓比為3時(shí)，錘激波在進(jìn)氣道內(nèi)傳播的平均速度為379 m/s(Ma=1.11)。

圖14 不同超壓比時(shí)上下壁面各點(diǎn)的峰值壓力Fig.14 Peak pressure of each point on upper and lower walls for different Over Pressure Ratios(OPRs)

圖15 不同超壓比時(shí)點(diǎn)P1處壓力隨時(shí)間變化曲線(xiàn)Fig.15 Time history of pressure at probe point P1 for different over pressure ratios

6 中心線(xiàn)幾何形狀的影響

進(jìn)氣道的幾何形狀會(huì)對(duì)錘激波壓力的分布情況產(chǎn)生影響。對(duì)于S形進(jìn)氣道幾何造型方面的研究較多，其中Lee和Boedicker[19]提出的關(guān)于中心線(xiàn)形狀的3種曲線(xiàn)方程得到了非常廣泛的應(yīng)用。本文選取其中兩條曲線(xiàn)與原曲線(xiàn)進(jìn)行對(duì)比，探討其錘激波壓力分布的差異。3條曲線(xiàn)如圖16所示，曲線(xiàn)1為原M2129型進(jìn)氣道中心線(xiàn)，曲線(xiàn)2和3來(lái)自文獻(xiàn)[19]，圖中橫縱軸分別表示中心線(xiàn)各點(diǎn)x坐標(biāo)和z坐標(biāo)與對(duì)應(yīng)方向總長(zhǎng)度的比值。對(duì)于前半段(Δx<0.5)，相比于曲線(xiàn)1，曲線(xiàn)2初始時(shí)曲率較小，方向變化緩慢，之后曲率迅速增大，而曲線(xiàn)3則一開(kāi)始曲率最大，方向變化最快，之后變得最為平緩。在后半段(Δx>0.5)，曲線(xiàn)2和3變化規(guī)律相反。計(jì)算時(shí)除了中心線(xiàn)形狀外，其他參數(shù)的設(shè)定完全相同。

圖16 進(jìn)氣道中心線(xiàn)幾何形狀Fig.16 Geometry shape of inlet centerlines

圖17 不同中心線(xiàn)時(shí)上下壁面各點(diǎn)的峰值壓力Fig.17 Peak pressure of each point on upper and lower walls with different centerlines

圖17給出了3種不同中心線(xiàn)形狀的進(jìn)氣道上下壁面各點(diǎn)在錘激波演化過(guò)程中的峰值壓力。對(duì)比可以看到峰值壓力的變化在進(jìn)氣道前半段的上壁面最為明顯，最大峰值壓力的位置和大小均有所不同。曲線(xiàn)2開(kāi)始最平緩之后曲率迅速變大，相應(yīng)的進(jìn)氣道靠近入口位置峰值壓力較小，最大峰值所處位置更靠后。曲線(xiàn)3開(kāi)始曲率很大，之后平緩，相應(yīng)的進(jìn)氣道模型受到最大峰值壓力的地方更靠前。而進(jìn)氣道后半段各點(diǎn)的峰值壓力受幾何形狀的影響較小。通過(guò)對(duì)比可以看出進(jìn)氣道中心線(xiàn)形狀的變化對(duì)其錘激波壓力分布具有重要的影響。但是峰值壓力還受到其他因素的影響，如所處進(jìn)氣道的位置、進(jìn)氣道截面積變化等。因而各處峰值壓力的變化并不與曲率分布完全對(duì)應(yīng)。

7 結(jié) 論

本文選取M2129型S形進(jìn)氣道模型進(jìn)行了錘激波動(dòng)態(tài)過(guò)程的數(shù)值模擬。利用在進(jìn)氣道出口設(shè)定壓力波形的方法模擬喘振時(shí)發(fā)動(dòng)機(jī)的壓力變化。詳細(xì)分析了來(lái)流馬赫數(shù)為0.65，進(jìn)氣道初始換算流量為60 kg/s的條件下錘激波經(jīng)過(guò)進(jìn)氣道時(shí)壁面的壓力變化和內(nèi)部的三維流動(dòng)情況。并進(jìn)一步分析了不同換算流量、超壓比和進(jìn)氣道中心線(xiàn)幾何形狀對(duì)進(jìn)氣道所受錘激波壓力的影響。結(jié)果表明：

1) 在馬赫數(shù)為0.65，進(jìn)氣道初始換算流量為60 kg/s(出口馬赫數(shù)Mae=0.44)的條件下，錘激波以相對(duì)進(jìn)氣道大約306 m/s(Ma=0.9)的速度向上游傳播，相對(duì)激波上游來(lái)流的速度約481 m/s(Ma=1.42)。在錘激波傳播過(guò)程中伴隨三維復(fù)雜流動(dòng)結(jié)構(gòu)和回流的產(chǎn)生。通過(guò)記錄進(jìn)氣道上下壁面各點(diǎn)的壓力隨時(shí)間變化情況，發(fā)現(xiàn)進(jìn)氣道上下壁面所受峰值動(dòng)態(tài)壓力達(dá)到了來(lái)流靜壓的2.3～2.65倍。錘激波經(jīng)過(guò)進(jìn)氣道彎道處時(shí)，彎道外側(cè)的壁面壓力會(huì)大于內(nèi)側(cè)，尤其靠近入口的彎折處，上壁面所受的峰值壓力為整個(gè)進(jìn)氣道所受錘激波壓力最高的部位。在錘激波傳出進(jìn)氣道后，進(jìn)氣道內(nèi)部流動(dòng)不會(huì)立刻恢復(fù)到初始流動(dòng)，而是經(jīng)歷了由進(jìn)氣道中部同時(shí)向入口和出口流出以及隨后變?yōu)橥耆亓鞯攘鲃?dòng)狀態(tài)，經(jīng)過(guò)了大約0.185 s后才恢復(fù)為正常狀態(tài)。

2) 通過(guò)對(duì)不同換算流量條件下錘激波演化過(guò)程的對(duì)比，發(fā)現(xiàn)較低流量時(shí)進(jìn)氣道所受峰值壓力更大，錘激波相對(duì)進(jìn)氣道傳播速度更快。對(duì)比不同超壓比時(shí)錘激波演化過(guò)程，發(fā)現(xiàn)較高超壓比時(shí)進(jìn)氣道受到峰值壓力更大，錘激波傳播速度更快。

3) 對(duì)具有不同中心線(xiàn)幾何形狀的進(jìn)氣道進(jìn)行錘激波計(jì)算，表明了中心線(xiàn)形狀對(duì)錘激波壓力的大小和分布有重要影響。不同中心線(xiàn)形狀時(shí)，S形進(jìn)氣道前半段的彎道所受的錘激波壓力產(chǎn)生了明顯的改變，曲率較大的地方易產(chǎn)生大的錘激波峰值壓力，而后半段受到的影響較小。

[1] YOUNG L C, BEAULIEU W D. Review of hammershock pressures in aircraft inlets[J]. Journal of Aircraft, 1975, 12(4): 210-216.

[2] MARSHALL F L. Prediction of inlet duct overpressures resulting from engine surge[J]. Journal of Aircraft, 1973, 10(5): 274-278.

[3] KURKOV A P, SOEDER R H, MOSS J E. Investigation of the stall hammershock at the engine inlet[J]. Journal of Aircraft, 1975, 12(4): 198-204.

[4] BELLMAN D, HUGHES D. The flight investigation of pressure phenomena in the air intake of an F-111A airplane: AIAA-1969-0488[R]. Reston, VA: AIAA, 1969.

[5] NUGENT J, HOLZMAN J K. Flight-measured inlet pressure transients accompanying engine compressor surges on the F-111A airplane: NASA TN D-7696[R]. Washington, D.C.: NASA, 1974.

[6] EVANS P J, TRUAX P. YF-16 air induction system design loads associated with engine surge[J]. Journal of Aircraft, 1975, 12(4):205-209.

[7] HINDASH I, BUSH R, COSNER R. Computational modeling of inlet hammershock wave generation: AIAA-1990-2005[R]. Reston, VA: AIAA, 1990.

[8] PAYNTER G, MAYER D, TJONNELAND E. Flow stability issues in supersonic inlet flow analyses: AIAA-1993-0290[R]. Reston, VA: AIAA, 1993.

[9] GOBLE B D, KING S, TERRY J, et al. Inlet hammershock analysis using a 3-D unsteady Euler/Navier-Stokes code: AIAA-1996-2547[R]. Reston, VA: AIAA, 1996.

[10] YTTERSTRCM A, AXELSON E. Hammershock calculations in the air intake of JAS 39 GFUPEN using dual timestepping: AIAA-1999-3113[R]. Reston, VA: AIAA, 1999.

[11] BECKER J. Dynamic hammershock effects on the air intake design of supersonic aircraft[M]∥Structures under Shock and Impact Ⅲ. Southampton: WIT Press, 1994: 419-426.

[12] HSIEH T, WARDLAW A B, COLLINS P, et al. Numerical investigation of unsteady inlet flowfields[J]. AIAA Journal, 1987, 25(1): 75-81.

[13] GRIDLEY M, SYLVESTER T, TRUAX P. Impact of a probabilistic approach on inlet hammershock design loads: AIAA-1999-2114[R]. Reston, VA: AIAA, 1999.

[14] 朱宇, 沈天榮. 飛機(jī)進(jìn)氣道錘擊波載荷評(píng)估方法研究[J]. 航空發(fā)動(dòng)機(jī), 2015， 41(3): 6-11.

ZHU Y, SHEN T R. Evaluation approach of hammer shock loading for aircraft inlet[J]. Aeroengine, 2015, 41(3): 6-11 (in Chinese).

[15] ANDERSON B H, REDDY D R, KAPOOR K. Study on computing separating flows within a diffusion inlet S-duct[J]. Journal of Propulsion & Power, 1994, 10(5): 661-667.

[16] NICHOLS R. Calculation of the flow in a circular S-duct inlet: AIAA-1991-0174 [R]. Reston, VA: AIAA, 1991.

[17] STANLEY R, MOHLER J. Wind-US flow calculations for the M2129 S-duct using structured and unstructured grids: AIAA-2004-0525[R]. Reston, VA: AIAA, 2004.

[18] DEBIASI M, HERBERG M R, YAN Z, et al. Control of flow separation in S-ducts via flow injection and suction: AIAA-2008-0074[R]. Reston, VA: AIAA, 2008.

[19] LEE C, BOEDICKER C. Subsonic diffuser design and performance for advanced fighter aircraft: AIAA-1985-3073[R]. Reston, VA: AIAA, 1985.