近臨界馬赫數(shù)條件下楔形誘導(dǎo)斜爆轟數(shù)值模擬

王成，吳京昌, 韓文虎, 楊同會(huì)

(北京理工大學(xué) 爆炸科學(xué)與技術(shù)重點(diǎn)實(shí)驗(yàn)室，北京 100081)

爆轟波是由激波和波后化學(xué)反應(yīng)放熱相耦合而成并以超聲速傳播的燃燒波. 由于爆轟波在高超聲速推進(jìn)系統(tǒng)中的巨大應(yīng)用潛力，其相關(guān)研究得到廣泛重視. 以爆轟波為基礎(chǔ)的推進(jìn)系統(tǒng)存在多種類別，如旋轉(zhuǎn)爆轟發(fā)動(dòng)機(jī)[1]及采用駐定斜爆轟波的沖壓加速器、斜爆轟發(fā)動(dòng)機(jī)[2]. 以斜爆轟推進(jìn)系統(tǒng)為基礎(chǔ)的相關(guān)研制取得了一定進(jìn)步，斜爆轟發(fā)動(dòng)機(jī)適應(yīng)高馬赫數(shù)飛行，且結(jié)構(gòu)簡(jiǎn)單，由于斜爆轟反應(yīng)區(qū)極短，因此燃燒室無須很長(zhǎng)，可大大減少發(fā)動(dòng)機(jī)的質(zhì)量以及壁面損失. 如駐定斜爆轟波穩(wěn)定，則系統(tǒng)無需采用火焰助穩(wěn). 然而斜爆轟波的不穩(wěn)定性，使得形成穩(wěn)定斜爆轟波較為困難. 因此，從內(nèi)在爆轟物理機(jī)制方面對(duì)斜爆轟波系結(jié)構(gòu)及不穩(wěn)定性進(jìn)行研究，對(duì)斜爆轟發(fā)動(dòng)機(jī)的設(shè)計(jì)具有重要參考意義.

為探索ODW的基本特性及其在推進(jìn)系統(tǒng)中的應(yīng)用，已開展大量研究. Lehr等[3]在觀測(cè)氫氧介質(zhì)中高速?gòu)椡栾w行時(shí)，發(fā)現(xiàn)彈頭部分產(chǎn)生斜激波誘導(dǎo)的燃燒現(xiàn)象，并同時(shí)觀測(cè)到不同燃料當(dāng)量比、彈丸結(jié)構(gòu)及飛行馬赫數(shù)等條件下，斜爆轟的結(jié)構(gòu)不同. Li等[4]通過數(shù)值模擬發(fā)現(xiàn)斜爆轟波由斜激波、誘導(dǎo)區(qū)、斜爆轟波波面和一組爆燃波組成. 斜激波下方氣體經(jīng)過壓縮與自燃，形成爆燃波. Viguier等[5]進(jìn)行了以氫氣為燃料的斜爆轟研究，實(shí)驗(yàn)所得斜爆轟結(jié)構(gòu)穩(wěn)定且實(shí)驗(yàn)重復(fù)性好，與Li的數(shù)值模擬結(jié)果基本一致. Palalexandris等[6]首次發(fā)現(xiàn)斜爆轟波波面上可能存在三波點(diǎn)及橫波，并發(fā)現(xiàn)這種結(jié)構(gòu)往往出現(xiàn)在較大化學(xué)反應(yīng)活化能的可燃?xì)怏w中. Choi等[7]通過高精度模擬，發(fā)現(xiàn)有三波點(diǎn)和橫波構(gòu)成的小尺度多維結(jié)構(gòu)確實(shí)存在，且與正爆轟波的胞格結(jié)構(gòu)明顯不同. Silva等[8]通過氫氧混合基元反應(yīng)模型，分析了不同楔面角、不同氫氧當(dāng)量比及不同馬赫數(shù)條件下斜爆轟波的發(fā)展演化，并對(duì)斜爆轟波誘導(dǎo)時(shí)間和誘導(dǎo)長(zhǎng)度進(jìn)行了研究. Liu和Wang等[9]從理論和數(shù)值模擬上發(fā)現(xiàn)Chapman-Jouguet 斜爆轟波(C-J ODW)在ODW結(jié)構(gòu)中起重要作用，并解釋了C-J ODW可以上行傳播的原因. Teng等[10-12]對(duì)斜激波向斜爆轟波轉(zhuǎn)變區(qū)域的爆燃波與激波耦合結(jié)構(gòu)進(jìn)行了研究，歸明月等[13]則進(jìn)一步解釋了斜爆轟精細(xì)結(jié)構(gòu)，并分別探討了尖劈角度及來流馬赫數(shù)對(duì)劈面壓縮性的影響. 綜上所述，關(guān)于斜爆轟的研究大多在無粘性反應(yīng)流體范疇，且集中在高馬赫數(shù)條件強(qiáng)過驅(qū)的駐定斜爆轟行為，對(duì)于近臨界馬赫數(shù)條件下粘性斜爆轟傳播特征的研究仍然很少. 盡管Li等[14]、Fusina等[15]研究認(rèn)為黏性對(duì)斜爆轟波的結(jié)構(gòu)影響僅局限于薄邊界層，對(duì)主流區(qū)結(jié)構(gòu)幾乎無影響，而近年來的研究表明，黏性擴(kuò)散對(duì)爆轟波后未反應(yīng)物質(zhì)的燃燒及斜爆轟上行速度等方面存在影響，不可忽略. 因此開展有、無黏性擴(kuò)散的斜爆轟數(shù)值模擬，揭示斜爆轟的陣面結(jié)構(gòu)及其臨界馬赫數(shù)條件下的傳播特征.

1 控制方程

控制方程采用具有黏性、熱傳導(dǎo)及分子擴(kuò)散，帶化學(xué)反應(yīng)的2維多組分可壓縮流體Navier-Stokes方程組，具體形式為

-AρYexp(-Ea/RpT)，

p=ρRpT/M.

式中：Y為可燃?xì)怏w質(zhì)量分?jǐn)?shù)；h=QY+CpT=QY+RpγT/(γ-1)為焓；e=QY+CVT=QY+RpT/(γ-1)為內(nèi)能；Q為反應(yīng)熱；CV和Cp分別表示定容和定壓比熱容；M為可燃?xì)怏w摩爾質(zhì)量；Rp為理想氣體常數(shù)；γ為可燃?xì)怏w絕熱指數(shù)；u、v分別為流體速度在x方向和y方向的分量；ρ、p、T、Y分別為密度、壓力、溫度、反應(yīng)物的質(zhì)量分?jǐn)?shù).

能量擴(kuò)散矢量σj和應(yīng)力張量ζij為

式中：μ=pv為動(dòng)黏性系數(shù)；Pr和Sc分別為Prandtl和Schmidt數(shù)；v為分子運(yùn)動(dòng)黏度；δij為Kronecker符號(hào).

上述物理量及均質(zhì)氣體混合物的初始狀態(tài)為量綱一的參數(shù)，因此初始密度、壓力和溫度是統(tǒng)一的. 根據(jù)未燃?xì)怏w的狀態(tài)，使方程量綱歸一化，

表1 量綱一化參數(shù)[16]

2 數(shù)值計(jì)算方法及收斂性

2.1 數(shù)值計(jì)算方法

在楔形斜面上設(shè)定一個(gè)量綱一化的大小800×600矩形域，由于氣流以超聲速傳播，會(huì)在楔面上方形成一道以頂點(diǎn)為起點(diǎn)的斜激波，斜激波誘導(dǎo)可燃混合氣體發(fā)生燃燒反應(yīng)，在下游形成復(fù)雜的斜爆轟結(jié)構(gòu). 為方便數(shù)值模擬，將坐標(biāo)系旋轉(zhuǎn)使之與楔面方向一致(如圖1). 通過改變?nèi)肓黢R赫數(shù)來模擬固定楔形角為30°的斜爆轟. 空間離散采用基于WENO保守有限差分格式的5階局部特征方法，時(shí)間離散采用3階TVD Runge-Kutta時(shí)間離散. 并用數(shù)值算例驗(yàn)證了相對(duì)網(wǎng)格分辨率和保正數(shù)值格式的正確性[17].

2.2 收斂性驗(yàn)證

采用20，40及80 pts/L1/2(L1/2為半反應(yīng)區(qū)厚度)3種網(wǎng)格精度對(duì)Ma=12,楔面角30°的2維無黏斜爆轟進(jìn)行收斂性驗(yàn)證. 由圖2可知，當(dāng)網(wǎng)格精度增加到80 pts/L1/2,發(fā)現(xiàn)爆轟結(jié)構(gòu)變化很小，表明40 pts/L1/2的模擬結(jié)果收斂于80 pts/L1/2時(shí)的結(jié)果. 因此，40 pts/L1/2網(wǎng)格精度能夠捕捉斜爆轟結(jié)構(gòu)特征及爆轟角. Choi與Kim等[18]采用3階MUSCL-TVD方案對(duì)對(duì)流通量進(jìn)行了離散化，網(wǎng)格精度為40 pts/L1/2能夠滿足斜爆轟模擬要求. 由于目前的數(shù)值模擬采用的格式為5階WENO計(jì)算格式，占用較少網(wǎng)格就能達(dá)到較高精度，且經(jīng)過驗(yàn)證其精度可以達(dá)到5階3維爆轟[19]. 高網(wǎng)格精度下接觸面會(huì)產(chǎn)生強(qiáng)渦，Choi及Kim等[18]進(jìn)行了詳細(xì)的討論. 本模擬對(duì)于高馬赫數(shù)條件下的斜爆轟進(jìn)行模擬，高過驅(qū)度會(huì)抑制接觸面上的不穩(wěn)定特征，因此小尺度渦不會(huì)出現(xiàn)[20].

3 數(shù)值結(jié)果分析和討論

3.1 馬赫數(shù)對(duì)斜爆轟結(jié)構(gòu)的影響

根據(jù)爆轟極理論推導(dǎo)得出，圖3為楔面角30°時(shí)不同馬赫數(shù)下的斜爆轟角與來流馬赫數(shù)間的關(guān)系(Ma=7,8,9.5,12),對(duì)比理論解給出臨界馬赫數(shù)為8.5. 觀察圖像可知，數(shù)值模擬結(jié)果預(yù)測(cè)的臨界馬赫

數(shù)與理論值吻合較好. 隨著馬赫數(shù)改變，斜爆轟呈現(xiàn)出駐定、非駐兩種不同的斜爆轟模態(tài). 當(dāng)來流馬赫數(shù)高于臨界值時(shí)，斜爆轟在下游起爆，上行一段距離后駐定；來流馬赫數(shù)低于臨界值時(shí)，斜爆轟波不斷上行，最終脫體. 當(dāng)馬赫數(shù)大于或接近臨界值8.5時(shí)，斜爆轟從下游起爆后上行至上游，最終駐定.

圖4給出馬赫數(shù)為7時(shí)，t=108.8，264.1，417.9時(shí)的流場(chǎng)溫度分布與壓力分布. 可以看出t=108.8時(shí)，斜爆轟處于起爆的初始階段，斜激波與斜爆轟波的過渡通過主三波點(diǎn)完成，同時(shí)爆轟下游存在C-J斜爆，推動(dòng)斜爆轟波面位置向上游移動(dòng). 對(duì)比圖4(a)中的溫度、壓力分布，爆轟壓力與火焰陣面之間解耦，不能形成駐定的斜爆轟；隨著時(shí)間的推進(jìn)，斜激波起爆區(qū)的長(zhǎng)度逐漸縮短，斜激波向斜爆轟過渡的三波點(diǎn)向上游推進(jìn)；在t=264.1時(shí)，可清晰地看到主三波點(diǎn)連接一激波，該激波經(jīng)楔面反射，并透過滑移間斷誘發(fā)下游流場(chǎng)失穩(wěn)；在t=417.9時(shí)，斜爆轟脫體，沒有形成駐定爆轟波. 主三波點(diǎn)的縱坐標(biāo)最終接近x軸，且跑出計(jì)算域外，斜爆轟波后大面積燃燒反應(yīng)進(jìn)行，且滑移間斷面部分出現(xiàn)一系列明顯的渦結(jié)構(gòu).

圖5給出了t為60.4，111.1，164.4時(shí)，來流馬赫數(shù)為8時(shí)的流場(chǎng)溫度分布與壓力分布數(shù)值模擬結(jié)果. 觀察模擬結(jié)果可知t=60.4時(shí)，起爆的斜激波區(qū)域長(zhǎng)度較圖4(a)短，斜激波向斜爆轟波的過渡較為平緩. 而主三波點(diǎn)附近的結(jié)構(gòu)較為復(fù)雜，斜爆轟波陣面下游結(jié)構(gòu)及陣面后橫波傳播均呈水平直線. 從流場(chǎng)整體溫度均值及壓強(qiáng)均值來看，t=60.4時(shí)的爆轟溫度及爆轟壓力還較低；t=111.1時(shí)，主三波點(diǎn)附近的結(jié)構(gòu)與圖4時(shí)相近時(shí)刻的爆轟場(chǎng)相比較而言，已基本穩(wěn)定，觀察圖5(b)中的溫度云圖可以發(fā)現(xiàn)，爆轟波陣面后方的反應(yīng)十分劇烈，并在爆轟右上部波陣面后方出現(xiàn)一些未反應(yīng)物質(zhì). 與馬赫數(shù)為7的模擬結(jié)果比較，斜爆轟角較小，且滑移線上渦結(jié)構(gòu)更為復(fù)雜；t=164.4時(shí)，斜爆轟仍沒有駐定，主三波點(diǎn)繼續(xù)沿x軸上行，滑移線下游的渦結(jié)構(gòu)發(fā)展充分，渦呈水平排列且尺度均勻，此時(shí)的未反應(yīng)物主要集中在波陣面后右上端. 爆轟最終脫體，沒有實(shí)現(xiàn)駐定，上行至計(jì)算域外.

3.2 與有粘狀態(tài)下的流場(chǎng)結(jié)構(gòu)比較

圖6給出了Ma=9.5的無黏斜爆轟在t為81.9，178.5，274.7及324.9時(shí)的溫度場(chǎng)模擬結(jié)果. 可以看出t=81.9時(shí)，斜爆轟起爆的上游斜激波區(qū)域長(zhǎng)度較圖4(a)和5(a)更短，斜爆轟角較大；在x=100～150之間區(qū)域，明顯觀察到三波點(diǎn)后連接一激波，在發(fā)生反射，并透過滑移間斷面與爆轟波陣面產(chǎn)生的橫波相互作用，使下游流場(chǎng)受到擾動(dòng)；爆轟波產(chǎn)生的橫波遇楔面后也發(fā)生較明顯反射，使得接觸間斷面下游出現(xiàn)一系列渦結(jié)構(gòu).t=178.5時(shí)，渦線下方的渦結(jié)構(gòu)進(jìn)一步發(fā)展，斜爆轟角進(jìn)一步增大，主三波點(diǎn)上行速度漸緩；t=274.7時(shí)，斜爆轟的波面后方出現(xiàn)較明顯的胞格結(jié)構(gòu)，三波點(diǎn)及滑移線的縱坐標(biāo)水平變得很低，斜爆轟上行過程中出現(xiàn)全局波動(dòng)現(xiàn)象；主三波點(diǎn)上行至x=35處，而后在t=324.9時(shí)退至x=50處；在t=324.9時(shí)，主三波點(diǎn)誘發(fā)的激波，及滑移線上的渦結(jié)構(gòu)變得更加明顯，爆轟波陣面后右上端緊貼波面處保留著部分未反應(yīng)物質(zhì)，此時(shí)斜爆轟已實(shí)現(xiàn)駐定，不再繼續(xù)上行傳播.

圖7給出了黏性斜爆轟在t=32.9，78.6，134.9，222.2，274.7及324.9時(shí)的溫度場(chǎng)云圖. 與圖4(a)及圖5(a)比較發(fā)現(xiàn)，在t=32.9時(shí)，斜爆轟角更小，誘導(dǎo)區(qū)的長(zhǎng)度更短，在t=78.6時(shí)，x=100～150之間渦線下方區(qū)域的激波及反射結(jié)構(gòu)與圖6(a)相似,但其爆轟波陣面后反應(yīng)區(qū)內(nèi)的小尺度結(jié)構(gòu)更加規(guī)則和均勻. 比較圖6(d)、(e)、(f)中的主三波點(diǎn)橫坐標(biāo)位置可知，黏性斜爆轟形成過程中也出現(xiàn)上行波的波動(dòng). 通過對(duì)比有、無粘性模擬結(jié)果，可以看出黏性對(duì)過驅(qū)情況的影響更為顯著，而對(duì)于接近C-J狀態(tài)的情況影響較小. 對(duì)于接近C-J狀態(tài)的斜爆轟，黏性C-J爆轟上行傳播速度要比非黏性情況快. 整體上，黏性不影響斜爆轟全局特征，但黏性抑制了滑移線上不穩(wěn)定誘發(fā)的小尺度渦結(jié)構(gòu)的產(chǎn)生.

4 結(jié) 論

基于一步總包化學(xué)反應(yīng)模型，采用NS及Euler控制方程組對(duì)2維楔形斜爆轟進(jìn)行了數(shù)值模擬，得到主要結(jié)論如下.

① 隨著馬赫數(shù)值增大，斜爆轟呈現(xiàn)出駐定、非駐斜爆轟模態(tài). 臨界馬赫數(shù)約為8.5，來流馬赫數(shù)高于臨界值時(shí)，斜爆轟可實(shí)現(xiàn)駐定；來流馬赫數(shù)低于臨界值時(shí)，斜爆轟波不斷上行，成為脫體爆轟.

② 通過對(duì)比有無黏性兩種情況下的ODW模擬，可以看出黏性對(duì)過驅(qū)情況的影響更為顯著. 整體上，黏性不影響整體特征，但黏度效應(yīng)抑制了斜爆轟結(jié)構(gòu)跡線上因不穩(wěn)定引起的小規(guī)模旋渦. 對(duì)于接近C-J狀態(tài)的ODW，爆轟出現(xiàn)上行波動(dòng)現(xiàn)象，爆轟在粘性情況下的上行傳播速度要快于非粘性情況.