卓長(zhǎng)飛,武曉松,封 鋒
(南京理工大學(xué)機(jī)械工程學(xué)院,南京 210094)
超燃沖壓發(fā)動(dòng)機(jī)是實(shí)現(xiàn)高超聲速飛行的首要關(guān)鍵技術(shù),是目前世界各國(guó)競(jìng)相發(fā)展的熱點(diǎn)領(lǐng)域之一。常規(guī)超燃沖壓發(fā)動(dòng)機(jī)中,超聲速氣流在燃燒室內(nèi)駐留時(shí)間很短,要在燃燒室內(nèi)進(jìn)行燃料和空氣的有效混合,并實(shí)現(xiàn)完全燃燒,需要較長(zhǎng)的燃燒室,這就導(dǎo)致發(fā)動(dòng)機(jī)的結(jié)構(gòu)重量急劇增大,增加了飛行器重量和阻力,從而極大制約了飛行器性能的提高[1]。激波引燃沖壓發(fā)動(dòng)機(jī)采用爆轟形式組織燃燒,燃燒距離較短,燃燒室長(zhǎng)度可大大縮短,可有效減輕發(fā)動(dòng)機(jī)本身及壁面冷卻系統(tǒng)的結(jié)構(gòu)質(zhì)量,同時(shí)爆轟過(guò)程的燃燒效率較高,使得激波引燃沖壓發(fā)動(dòng)機(jī)逐漸受到重視。目前,國(guó)外只有少數(shù)國(guó)家開(kāi)展激波引燃沖壓發(fā)動(dòng)機(jī)的研究工作。Sislian[2]和Redford[3]研究了燃料/空氣不完全混合和非設(shè)計(jì)工作條件對(duì)激波引燃沖壓發(fā)動(dòng)機(jī)的影響;Alexander和Sislian[4]數(shù)值研究了來(lái)流馬赫數(shù)為11、高度為34.5 km 條件下內(nèi)噴注型的激波引燃沖壓發(fā)動(dòng)機(jī)流場(chǎng)結(jié)構(gòu);Ess[5]對(duì)激波引燃沖壓發(fā)動(dòng)機(jī)等截面燃燒室中鈍頭體誘導(dǎo)燃燒和爆轟進(jìn)行了數(shù)值模擬,驗(yàn)證了發(fā)動(dòng)機(jī)的燃燒放熱和推進(jìn)性能,也為發(fā)動(dòng)機(jī)設(shè)計(jì)提供了依據(jù)。國(guó)內(nèi)對(duì)激波引燃沖壓發(fā)動(dòng)機(jī)研究較少,黃偉等[1]指出了發(fā)展這種新型推進(jìn)系統(tǒng)的關(guān)鍵技術(shù),對(duì)國(guó)內(nèi)的研究思路提出了建議;盛德林等[6]指出激波引燃沖壓發(fā)動(dòng)機(jī)將是一種很有前途的新型高超聲速推進(jìn)方法。
本文采用AUSMPW+迎風(fēng)格式,在以非結(jié)構(gòu)網(wǎng)格為控制體中,采用有限體積法求解二維帶化學(xué)反應(yīng)的多組分氣體Euler方程,對(duì)內(nèi)噴注型激波引燃沖壓發(fā)動(dòng)機(jī)一體化流場(chǎng)進(jìn)行數(shù)值模擬,研究了臺(tái)階長(zhǎng)度、斜劈尖角度對(duì)發(fā)動(dòng)機(jī)流場(chǎng)結(jié)構(gòu)和發(fā)動(dòng)機(jī)性能參數(shù)影響規(guī)律。
笛卡爾坐標(biāo)系下,微分守恒形式的二維化學(xué)非平衡流Euler方程:
其中
式中 Q為守恒變量;F、G為無(wú)粘通量;S為化學(xué)反應(yīng)源項(xiàng);ρ為混合密度;u、v分別為2個(gè)方向速度;p為混合壓力;E為混合氣體單位質(zhì)量的總能;fi(i=1,...,N-1)為i組分質(zhì)量分?jǐn)?shù);ωi為i組分質(zhì)量生成速率;N為總組分?jǐn)?shù)。
AUSMPW+迎風(fēng)格式執(zhí)行效率較高,在高壓梯度區(qū)域近似于 Van Leer格式,在低壓梯度區(qū)近似于AUSM+格式。該格式在高超聲速流計(jì)算中不會(huì)出現(xiàn)“carbuncle”現(xiàn)象,在間斷和邊界層內(nèi)都能得到精確解。AUSMPW+格式的數(shù)值通量項(xiàng)如下形式:
式中 w為壓力加權(quán)函數(shù);c為界面聲速;M為界面間斷處左右馬赫數(shù);其余變量含義參見(jiàn)文獻(xiàn)[7]。
為提高空間離散精度,在進(jìn)行對(duì)流通量迎風(fēng)格式離散之間,首先應(yīng)進(jìn)行物理量重構(gòu)。本文利用變量線性重構(gòu)技術(shù),以獲得空間二階精度。物理量重構(gòu)表達(dá)式如下:
式中 U表示原始變量;i,j分別為左、右單元中心點(diǎn)編號(hào);k為單元邊界編號(hào);r為兩點(diǎn)位置矢量;?U表示原始變量的梯度;φ表示抑制振蕩的限制器。
氫氣和氧氣的化學(xué)反應(yīng)機(jī)理采用7組分8步基元反應(yīng)模型,詳細(xì)反應(yīng)表達(dá)式與參數(shù)如表1所示。表1中,n為基元反應(yīng)級(jí)數(shù);阿烏累里斯公式表達(dá)式為Kf=ATbexp(-E/R0/T);第三碰撞體系數(shù):H2為 2.0,H2O為 6.0,其余組分均取 1.0。
表1 H2和O2基元反應(yīng)模型Table 1 Detail reaction model for H2+O2
化學(xué)非平衡流控制方程組分為流動(dòng)部分和化學(xué)反應(yīng)部分,兩者相互耦合,并會(huì)產(chǎn)生剛性問(wèn)題。本文采用時(shí)間算子分裂算法處理這種耦合過(guò)程,即首先計(jì)算流體流動(dòng)效應(yīng),得到物理變量的過(guò)渡值。然后繼續(xù)計(jì)算,將化學(xué)反應(yīng)的貢獻(xiàn)疊加到物理變量過(guò)渡值,最終得到下一時(shí)刻體現(xiàn)整體效應(yīng)的物理量值。其中,在計(jì)算化學(xué)反應(yīng)對(duì)流場(chǎng)的貢獻(xiàn)時(shí),需要把求解流動(dòng)偏微分方程時(shí)采用的時(shí)間步長(zhǎng)進(jìn)一步細(xì)分,作為求解化學(xué)反應(yīng)剛性常微分方程的步長(zhǎng)。具體做法是先凍結(jié)化學(xué)反應(yīng),求解得到流場(chǎng)參數(shù);然后,將化學(xué)反應(yīng)看做等容吸熱或放熱過(guò)程,保持內(nèi)能、速度參數(shù)不變,計(jì)算各組分的質(zhì)量變化率;最后,迭代求解溫度。
圖1為內(nèi)部噴注型激波引燃沖壓發(fā)動(dòng)機(jī)[1]。其工作原理:首先,高超聲速來(lái)流進(jìn)入燃燒室,與在燃燒室進(jìn)口噴注的燃料進(jìn)行混合;然后,高超聲速的燃料-空氣混合氣體在燃燒室末端的斜坡上產(chǎn)生較強(qiáng)的斜激波,波后溫度直接引燃預(yù)混氣體;最后,通過(guò)斜爆轟后充分燃燒的超聲速氣流通過(guò)尾噴管排出,并產(chǎn)生推力。
圖1 內(nèi)噴注型激波引燃沖壓發(fā)動(dòng)機(jī)典型結(jié)構(gòu)Fig.1 Typical structure of the shock-induced combustion ramjet engine based internal jet
本文計(jì)算采用文獻(xiàn)[4]給出的基本模型尺寸,由于本文計(jì)算的是二維問(wèn)題,沒(méi)考慮Z方向流動(dòng),假設(shè)本文模型的Z方向?yàn)?.1 m(即發(fā)動(dòng)機(jī)寬度),這方便對(duì)發(fā)動(dòng)機(jī)整體進(jìn)行推力積分計(jì)算。采用非結(jié)構(gòu)網(wǎng)格離散計(jì)算區(qū)域,如圖2所示。高超聲速來(lái)流條件:來(lái)流馬赫數(shù)11.0(設(shè)計(jì)馬赫數(shù)),來(lái)流速度為3 391 m/s,溫度為235 K,壓力為601 Pa。燃料噴注條件:燃料為 H2,噴注位置為燃燒室進(jìn)口擴(kuò)張位置的上、下壁面處,噴注溫度300 K,噴注壓力為15 000 Pa,噴注速度為1 500 m/s。燃燒室末端誘導(dǎo)爆轟燃燒的斜劈尖(記斜劈尖角度為α)和臺(tái)階(記臺(tái)階長(zhǎng)度為L(zhǎng))直接影響發(fā)動(dòng)機(jī)性能,本文針對(duì)這2個(gè)重要結(jié)構(gòu)參數(shù),計(jì)算了表2中所示的工況。
圖2 激波引燃沖壓發(fā)動(dòng)機(jī)計(jì)算模型與網(wǎng)格分布Fig.2 Model and grid of shock-induced combustion ramjet engine
表2 激波引燃沖壓發(fā)動(dòng)機(jī)計(jì)算工況Table 2 Calculation condition of shock-induced combustion ramjet engine
圖3~圖6給出了不同臺(tái)階長(zhǎng)度的激波引燃沖壓發(fā)動(dòng)機(jī)流場(chǎng)溫度云圖、密度紋影圖及爆轟產(chǎn)物H2O質(zhì)量分?jǐn)?shù)。
圖3 工況1發(fā)動(dòng)機(jī)流場(chǎng)與參數(shù)分布Fig.3 Flow field and parameter distribution of engine under case 1
圖4 工況2發(fā)動(dòng)機(jī)流場(chǎng)與參數(shù)分布Fig.4 Flow field and parameter distribution of engine under case 2
圖5 工況3發(fā)動(dòng)機(jī)流場(chǎng)與參數(shù)分布Fig.5 Flow field and parameter distribution of engine under case 3
圖6 工況4發(fā)動(dòng)機(jī)流場(chǎng)與參數(shù)分布Fig.6 Flow field and parameter distribution of engine under case 4
從圖3(a)可看出,激波引燃沖壓發(fā)動(dòng)機(jī)一體化流場(chǎng)的一些基本結(jié)構(gòu):高超聲速來(lái)流在外壓縮進(jìn)氣道斜面產(chǎn)生第一道斜激波,斜激波剛好相交于機(jī)身下端的整流罩,斜激波后溫度、壓力升高;被斜激波壓縮的來(lái)流經(jīng)過(guò)一定長(zhǎng)度的平直段,再通過(guò)一定的擴(kuò)張斜坡進(jìn)入燃燒室,上、下兩斜坡處存在噴注燃料H2的入口,H2通過(guò)垂直于斜坡的方式噴入來(lái)流中,燃燒室內(nèi)波系較為復(fù)雜;來(lái)流空氣與燃料預(yù)混后通過(guò)燃燒室末端的斜劈尖上產(chǎn)生斜爆轟波,混合氣體通過(guò)斜爆轟波時(shí)被迅速點(diǎn)燃,并充分燃燒,爆轟燃燒產(chǎn)物進(jìn)入尾噴管膨脹加速,產(chǎn)生推力;欠膨脹的爆轟燃燒產(chǎn)物流出尾噴管后繼續(xù)膨脹,與外流相互擠壓,并形成一道激波。
從圖3(a)、(b)可看出,混合氣體在斜劈尖前形成斜激波,斜激波強(qiáng)度較強(qiáng),波后溫度直接引燃混合氣體,燃燒波與斜激波耦合形成斜爆轟波。由于來(lái)流空氣過(guò)剩,可斷定,通過(guò)斜激波后H2基本燃燒完畢。斜激波剛好在燃燒室與尾噴管上壁面交點(diǎn)處發(fā)生反射,但由于反射激波前、后氣體均處于尾噴管中,超聲速氣流發(fā)生膨脹加速,這一過(guò)程嚴(yán)重干擾反射激波的形成,導(dǎo)致反射激波較為模糊。從圖3(c)可看出,由于斜劈尖產(chǎn)生的斜激波前的預(yù)混氣體并不均勻,由燃料噴注形式、位置決定了H2主要集中在燃燒室上、下壁面處,主流中心分布較少。因此,預(yù)混氣體通過(guò)斜激波引燃,生成的H2O主要集中在尾噴管上、下表面。
從圖4(a)、(b)可看出,斜激波第一反射點(diǎn)處于燃燒室內(nèi),因此產(chǎn)生正常的第一反射激波。當(dāng)?shù)谝环瓷浼げㄟ€沒(méi)達(dá)到第二反射壁面點(diǎn)時(shí),已經(jīng)處于尾噴管中,第一反射激波前、后氣體受尾噴管膨脹加速的影響,導(dǎo)致部分第一反射激波發(fā)生完全變形,并在尾噴管下壁面上發(fā)生第二次反射,但第二反射激波強(qiáng)度較弱,形狀較為彎曲。從圖4(c)可見(jiàn),與case 1中分布有較大不同,由于本工況的臺(tái)階較長(zhǎng),產(chǎn)生反射激波,流動(dòng)發(fā)生偏轉(zhuǎn),主要集中在燃燒室上、下壁面的爆轟產(chǎn)物H2O與主流中剩余的空氣混合更加充分,然后再流向尾噴管。因此,尾噴管中H2O分布更加廣泛,并不像case 1中主要集中于尾噴管上、下壁面處。
從圖5(a)、(b)和圖6(a)、(b)可看出,隨著臺(tái)階長(zhǎng)度的增加,第一反射激波不受尾噴管的影響,第二次反射點(diǎn)出現(xiàn)在臺(tái)階上,經(jīng)2次反射后,燃燒室末端即為噴管進(jìn)口處的氣體溫度大大升高,這樣將爆轟產(chǎn)物氣體的能量充分轉(zhuǎn)為熱能,有利于氣流在尾噴管內(nèi)的膨脹、加速,產(chǎn)生更大的推力。從圖5(c)、圖6(c)可見(jiàn),趨勢(shì)與case 2基本相同,均是由于臺(tái)階較長(zhǎng),產(chǎn)生多次反射激波,而導(dǎo)致主要集中在燃燒室上、下壁面的爆轟產(chǎn)物H2O與主流中剩余的空氣混合更加充分,尾噴管中H2O分布更加廣泛、均勻。
良渚玉器中也有龍的圖案,主要出現(xiàn)在玉鐲形器、玉圓牌、玉璜等器的邊緣上?,幧?號(hào)墓出土的玉鐲形器外緣刻有四個(gè)凸面,分別刻出四個(gè)龍首。又如瑤山2號(hào)墓出土的玉圓牌,見(jiàn)圖6,其邊緣也有三具龍首紋。
圖7給出了4種工況下燃燒室中H2質(zhì)量分?jǐn)?shù)分布。H2與空氣完全反應(yīng)時(shí)所需 H2的質(zhì)量分?jǐn)?shù)為0.028左右,因此把整個(gè)流場(chǎng)中H2質(zhì)量分?jǐn)?shù)顯示分布范圍調(diào)整為0~0.028,H2質(zhì)量分?jǐn)?shù)大于0.028的區(qū)域均用紅色顯示,這樣利于了解燃料與空氣的預(yù)混程度??煽闯?,H2從燃燒室進(jìn)口的上、下壁面垂直噴進(jìn)入燃燒室,與來(lái)流空氣進(jìn)行混合,但H2從壁面附近噴入,而來(lái)流速度較高,在燃燒室滯留時(shí)間較短,因而H2和來(lái)流混合并不均勻,H2主要集中在燃燒室上、下壁面處。隨著H2與來(lái)流空氣的流動(dòng)與擴(kuò)散混合,燃燒室流場(chǎng)中主流中開(kāi)始出現(xiàn)H2,但質(zhì)量分?jǐn)?shù)較低。從斜劈尖前的激波面看出,斜激波前H2仍集中于燃燒室上、下壁面,質(zhì)量分?jǐn)?shù)沒(méi)達(dá)到化學(xué)恰當(dāng)比所需要求,空氣是過(guò)量的,特別是在斜激波中間,H2所占的質(zhì)量分?jǐn)?shù)遠(yuǎn)低于完全反應(yīng)時(shí)所需H2的質(zhì)量分?jǐn)?shù)。因此,激波引燃沖壓發(fā)動(dòng)機(jī)技術(shù)關(guān)鍵問(wèn)題之一是使通過(guò)斜劈尖前的混合氣體應(yīng)盡量混合均勻。
圖7 燃燒室內(nèi)H2質(zhì)量分?jǐn)?shù)分布Fig.7 Distribution of H2mass fraction in combustion chamber
表3給出了4種工況下對(duì)應(yīng)的推力和基于燃料的比沖??煽闯?,隨著臺(tái)階長(zhǎng)度的增大,發(fā)動(dòng)機(jī)的推力也增大。但隨著臺(tái)階長(zhǎng)度的增大,推力增大的程度降低。由此推測(cè),應(yīng)存在一極限值,當(dāng)臺(tái)階長(zhǎng)度超過(guò)某一個(gè)值后,不再影響發(fā)動(dòng)機(jī)的推力。但在燃燒室長(zhǎng)度一定的情況下,臺(tái)階長(zhǎng)度過(guò)長(zhǎng),H2與空氣混合的距離將變短,不利于均勻混合。臺(tái)階長(zhǎng)度對(duì)發(fā)動(dòng)機(jī)性能參數(shù)影響的機(jī)理是臺(tái)階長(zhǎng)度影響斜激波的反射次數(shù),從而影響爆轟波后氣體溫度和推力。
表3 不同臺(tái)階長(zhǎng)度的發(fā)動(dòng)機(jī)性能參數(shù)Table 3 Performance of engine variation with bench length L
圖8、圖9分別給出了α=21°、41°的流場(chǎng)中溫度云圖、密度紋影圖及H2O質(zhì)量分?jǐn)?shù)分布。結(jié)合圖8中α=31°流場(chǎng)參數(shù)分布可看出,當(dāng)角度較小時(shí),預(yù)混氣體經(jīng)斜劈尖壓縮后產(chǎn)生的斜激波強(qiáng)度較弱,斜激波后溫度達(dá)不到直接引燃預(yù)混氣體的條件,經(jīng)過(guò)一定的誘導(dǎo)距離后發(fā)生著火,預(yù)混氣體發(fā)生燃燒。該狀態(tài)下預(yù)混氣體燃燒效率也較低,預(yù)混氣體的能量沒(méi)有得到充分釋放。當(dāng)α=31°時(shí),預(yù)混氣體經(jīng)斜劈尖壓縮后產(chǎn)生的斜激波強(qiáng)度較強(qiáng),波后溫度達(dá)到直接引燃預(yù)混氣體的條件,斜激波與燃燒波耦合形成斜爆轟波。當(dāng)α=41°時(shí),斜激波強(qiáng)度較強(qiáng)而形成斜爆轟波,斜爆轟波后溫度較高,約為3 300 K。同時(shí),斜爆轟波在燃燒室上壁面發(fā)生馬赫反射形成馬赫桿,預(yù)混氣體通過(guò)馬赫桿后的爆轟溫度遠(yuǎn)高于通過(guò)斜激波后的爆轟溫度,最高溫度能達(dá)到近4 000 K,預(yù)混氣體的能量得到充分釋放,再通過(guò)尾噴管的膨脹加速,能產(chǎn)生更大的推力。
表4給出了3種不同斜劈尖角度下發(fā)動(dòng)機(jī)性能參數(shù)。從表4可看出,隨斜劈尖角度的增大,發(fā)動(dòng)機(jī)推力減少。對(duì)于基于燃料的比沖而言,也是隨斜劈尖角度的增大而減少。實(shí)際上,增加斜劈尖角度會(huì)有兩方面作用:第一,斜劈尖表面壓力增大,即會(huì)增加發(fā)動(dòng)機(jī)的阻力;第二,預(yù)混氣體燃燒效率大大提高,能量充分釋放,通過(guò)尾噴管而產(chǎn)生更大的推力,而新增加的推力小于增加的阻力。因此,隨斜劈尖角度的增加,發(fā)動(dòng)機(jī)的凈推力和比沖是減少的。
圖8 工況5發(fā)動(dòng)機(jī)流場(chǎng)與參數(shù)分布Fig.8 Flow field and parameter distribution of engine under case 5
圖9 工況6發(fā)動(dòng)機(jī)流場(chǎng)與參數(shù)分布Fig.9 Flow field and parameter distribution of engine under case 6
表4 不同斜劈尖角度的發(fā)動(dòng)機(jī)性能參數(shù)Table 4 Performance of engine variation with oblique wedge angle α
(1)本文發(fā)展的數(shù)值方法能用于激波引燃沖壓發(fā)動(dòng)機(jī)的一體化流場(chǎng)和性能預(yù)示計(jì)算。
(2)臺(tái)階長(zhǎng)度影響激波反射次數(shù)和位置,從而影響爆轟產(chǎn)物H2O與剩余空氣的混合程度,導(dǎo)致爆轟產(chǎn)物H2O分布位置有較大差別;推力和燃料比沖均隨著臺(tái)階長(zhǎng)度的增加而呈增大的趨勢(shì),且增大程度逐漸降低。
(3)斜劈尖角度影響初始斜激波強(qiáng)度,從而影響爆轟波產(chǎn)生的位置。斜劈尖角度較小時(shí),預(yù)混氣體通過(guò)斜激波后,經(jīng)過(guò)一定誘導(dǎo)距離開(kāi)始燃燒放熱,斜劈尖角度較大時(shí),激波與燃燒波耦合形成較強(qiáng)爆轟波,同時(shí)在燃燒室上面形成較強(qiáng)的馬赫反射,波后溫度達(dá)4 000 K,預(yù)混氣體能量得到充分釋放;推力和燃料比沖均隨斜劈尖角度的增加而呈減小的趨勢(shì)。
[1] 黃偉,覃慧,羅世彬,等.激波誘燃沖壓發(fā)動(dòng)機(jī)關(guān)鍵技術(shù)研究[J].中國(guó)科學(xué):技術(shù)科學(xué),2010,40(1):64-70.
[2] Sislian J P,Dudebout R,Schumacher J,et al.Incomplete mixing and off-design effects on shock-induced combustion ramjet performance[J].Journal of Propulsion and Power,2000,16(1):41-48.
[3] Redford T.Effects of incomplete fuel-air mixing on the performance characteristics of mixed compression,shock-induced combustion ramjet(shcramjet)engines[D].Dissertation of Masteral Degree.Toronto:University of Toronto Institute,1998.
[4] Alexander D C,Sislian J P.Computational study of the propulsive characteristics of a shcramjet engine[J].Journal of Propulsion and Power,2008,24(1):34-44.
[5] Ess P R,Sislian J P,Allen C B.Blunt-body generated detonation in viscous hypersonic ducted flows[J].Journal of Propulsion and Power,2005,21(4):667-680.
[6] 盛德林,李文杰.一種新概念高超聲速推進(jìn)方案-激波引燃沖壓發(fā)動(dòng)機(jī)[J].飛航導(dǎo)彈,2010(6):44-46.
[7] Kim K H.Accurate computation of hypersonic flows using AUSMPW+ scheme and shock-induced grid technique[R].AIAA 98-2442.