兩相旋轉爆震發動機出口流場均勻性數值研究

馬 元，鄒 剛，孔光明，孫 靜，王文鵬，張宇坤

( 海軍航空大學 青島校區，山東 青島 266100)

相比于傳統的等壓燃燒方式，爆震燃燒具有熵增低、自增壓以及更高的熱循環效率等優點。旋轉爆震發動機(Rotating Detonation Engine，RDE)是基于爆震燃燒的一種動力裝置，它利用爆震波在環形燃燒室內的一端連續旋轉傳播，爆震產物從另一端開口處高速排出并產生連續推力。由于RDE沒有壓氣機、渦輪等轉動部件，故相比于傳統發動機，具有結構簡單、研制成本低的特點。此外，RDE還具有大比沖、高推重比、寬工作范圍等性能優勢，這對于開發新型高效航空航天發動機極具吸引力[1-2]。

近年來，國內外學者采用數值模擬方法對爆震燃燒室的流場特性和爆震旋轉機理等進行了探索和研究并取得了一定的進展[3-18]。姜宗林等[3]總結了氣相規則胞格爆震波傳播與起爆現象的相關研究進展，系統研究了氣相規則胞格爆震波起爆和傳播統一框架的幾個基本要素。滕宏輝等[4]采用數值模擬研究了氫氣-空氣混合物中二維斜爆震波的起始特性，分析了入口壓力和馬赫數對起始結構和長度的影響。王兵等[5]總結了有關爆炸、爆震及DDT的研究，討論了抑制管道內爆炸和爆震的被動/主動或混合方法。謝橋峰等[6]通過實驗研究了氫氣-空氣旋轉爆震燃燒室內貧油條件下的燃燒特性，得到了快速爆燃，不穩定爆震、準穩定爆震和穩定爆震等傳播模態，并給出了各傳播模態的影響因素和特點，對旋轉爆震燃燒室的設計和實際應用具有重要指導意義。Yi等[7]研究了噴管形狀、擴張角和長度對旋轉爆震發動機性能的影響。Schwer等[8-9]首次在數值模擬中將噴注腔與噴注孔考慮在內，對比研究了理想噴注模型與噴孔噴注模型RDE內流場結構及發動機性能的差異，并通過改變微型噴注孔的形狀及角度探索燃燒室爆震波對噴注腔壓力的影響規律。劉世杰等[10-11]模擬研究了旋轉爆震波的詳細結構和自持機理，并且發現當發動機尺寸小于臨界值時無法成功起爆。姜孝海等[12]采用9組分19步的詳細化學反應模型進行數值模擬，模擬結果表明爆震波作用后，受高溫膨脹及離心力等作用，產物擠向外壁，形成了一定的密度或壓力梯度。邵業濤等[13]采用兩步反應模型在矩形計算區域內模擬了爆震波的多個循環過程。馬虎等[14-15]基于Fluent軟件分析了旋轉爆震波的二維結構、入口總壓、外界反壓和燃燒室長度對RDE性能的影響。周蕊[16-17]通過對噴入流場中的粒子進行追蹤，對影響粒子軌跡及其物理參量變化的主要因素進行了分析，得到的RDE熱效率約為35.4%，并發現燃燒室半徑過大不利于爆震波的穩定傳播，波頭高度會周期性變化，采用多個波頭起爆能夠降低或消除大尺寸燃燒室內爆震波不穩定對發動機性能的不利影響。李寶星等[18]通過二維數值模擬研究了進氣總壓對旋轉爆震波的影響，獲得了燃燒室內部流場結構和爆震波傳播特性，對不同進氣總壓條件下發動機的爆震性能進行了計算。

早期研究大都針對氣態燃料，隨著研究的深入，液態燃料的旋轉爆震逐漸引起了人們的關注。在實際工程應用中受到重量和空間的限制，應用液態燃料更有優勢。與氣態燃料相比，液態燃料具有更高的密度比沖，相同體積下具有更高的沖量。鄭權等[19]對汽油為燃料，富氧空氣為氧化劑的連續旋轉爆震發動機進行了試驗，得到了不同當量比和不同汽油噴注壓力下爆震波的傳播特性，其中爆震波傳播頻率為2.1～2.4 kHz，傳播速度為1 022.2～1 171.8 m/s。王迪等[20-21]針對兩相連續旋轉爆震發動機的噴注器霧化混合、兩相爆震起爆與發動機熱試車試驗情況進行了總結，對俄羅斯、法國、波蘭在兩相連續旋轉爆震發動機有關的試驗研究情況進行了論述，對不同構型噴注器工作過程與性能、發動機點火起爆以及不同工質的發動機結構參數及試驗情況進行總結與分析。通過試驗得到爆震波的時域、頻域特征，對兩相連續旋轉爆震燃燒室中爆震波的起爆過程和穩定后的傳播過程進行了研究。馬元[22]以汽油為燃料，空氣為氧化劑，開展吸氣式氣液兩相RDE的試驗研究，分析了當量比、空氣噴注總溫、空氣噴注環縫、噴管等對氣液兩相旋轉爆震波傳播的影響。

上述學者的數值研究主要考察了氣態燃料下的旋轉爆震發動機，國內外對液態燃料旋轉爆震發動機數值模擬較少。兩相RDE出口流場均勻性對后續加裝的渦輪或噴管的設計具有重要的意義，而目前尚未有相關研究發表。本文對以辛烷為燃料、空氣為氧化劑的氣液兩相旋轉爆震發動機進行數值模擬，采用DPM模型計算辛烷液滴的軌跡，模擬非預混條件下的旋轉爆震波工作過程。對非預混噴注下兩相RDE的影響因素進行了研究，得到了空氣總溫與噴注均勻性對RDE出口流場均勻性的影響，為以后兩相RDE的深入研究提供參考。

1 數值方法和物理模型

本文利用商業軟件FLUENT，基于密度基求解器求解二維非穩態雷諾時均的N-S控制方程；對流項采用三階MUSCL格式離散，該格式對激波的捕捉具有較高的精度；物理通量采用AUSM矢通量分裂法進行分解；時間項采用四階龍格-庫塔法；采用標準k-ε湍流模型。以辛烷為燃料，空氣為氧化劑，化學反應方程式為C8H18+12.5(O2+3.76N2)→8CO2+9H2O+47N2。DPM模型考慮顆粒的耦合傳熱/傳質，考慮薩夫曼升力和壓力梯度力等附加力，考慮液滴的破碎、蒸發等物理過程，射流源采用組射流源。

1.1 模型和邊界條件

旋轉爆震燃燒室為柱狀環型燃燒室，相比于燃燒室的直徑，環形厚度相對較小，故將三維環形域簡化成二維計算域是合理的。因此為減少計算量節約成本，本文將燃燒室沿母線展開燃燒室內部流場簡化成二維矩形計算區域，尺寸為200 mm×50 mm，圖1中給出了監測點M1(x=100 mm，y=49 mm)的位置。本文網格尺度為0.5 mm×0.5 mm，為均勻正交網格。

圖1 旋轉爆震燃燒室二維計算域示意圖

本文計算域入口邊界的質量通量為100 kg/m2/s，空氣總溫為300～800 K，見表1所示。計算域的下邊界為壓力出口邊界，分兩種情況：當出口為亞音速時，邊界點壓力等于外界反壓，而其他守恒變量由內部流場外推得到，外界反壓為0.1 MPa。當出口為超音速時，所有守恒變量由內部區域外推得到。左右邊界定義為周期邊界，進行數據交換。

辛烷液滴的直徑為0.02 mm，初始溫度300 K，噴射速度20 m/s，辛烷的射流間距見表2所示，其中Case 15為理想均勻噴注。辛烷的總質量通量為6.667 kg/m2/s，總當量比為1。

表1 不同空氣噴注總溫(辛烷射流間距為2 mm)

表2 不同辛烷射流間距(空氣噴注總溫300 K)

1.2 數值驗證

圖2(a)為Case 5下爆震波穩定傳播時的溫度分布云圖，其中A是接近入口附近沿周向傳播的橫向爆震波，B是爆震波頭部產生的斜激波，C是滑移線，即新的爆震產物與上一循環的爆震產物形成的接觸間斷面，D是爆震產物與新鮮反應物的接觸面，E是新噴入的新鮮反應物，為旋轉爆震的傳播提供燃料。圖2(b)為 Bykovskii等[23]實驗獲得的旋轉爆震波結構，可以看出數值模擬的燃燒室流場結構與實驗結果定性一致。

圖2 Case 5的溫度云圖與Bykovskii實驗照片

圖3為Case 15監測點M1的壓力和溫度時程曲線。由圖3可知，在爆震波穿過監測點時，壓力均呈現急升緩降的特點，每個周期的壓力峰值基本保持不變，表明爆震波處于穩定傳播狀態，且壓力曲線與溫度曲線吻合較好，表明爆震陣面上激波與化學反應區高度耦合。計算得到平均峰值壓力為4.89 MPa，平均峰值溫度為 2 898 K，爆震波的平均傳播速度為1 649 m/s。與用Chemical Equilibrium with Application(簡稱CEA)計算得到的爆震波理論C-J壓力、C-J溫度和C-J速度的相對誤差分別為 -5.48%、-0.67%、-8.84%。

圖3 監測點M1的壓力和溫度曲線

1.4 網格無關性驗證

在數值方法準確性的基礎上，本文對二維RDE流場進行了網格無關性驗證。計算域為圖4所示二維計算域，網格尺度分別為1.0 mm、0.5 mm和0.25 mm，工況條件為 Case 5。由上節爆震波參數驗證表明0.5 mm的網格計算精度可以接受，本文不再對爆震波參數方面的網格無關性進行贅述。

圖4 不同網格尺度下燃燒室內溫度云圖

圖4為不同網格尺度下的燃燒室內溫度云圖，由圖4(a)可以看出，在網格尺度1.0 mm下，流場特征非常不明顯，如燃燒產物與新鮮反應物的接觸斷面、滑移線、斜激波等均不穩定，無法滿足流場分析的要求。由圖4(b)和圖4(c)可知，網格尺度0.5 mm和網格尺度0.25 mm下流場特征均非常明顯，爆震波、斜激波及接觸斷面等都很清晰。與 0.5 mm 下的流場相比，0.25 mm下的流場細節特征更加明顯，如爆震波后辛烷射流對流場的影響更加清晰，爆震產物與新鮮反應物的接觸面及滑移線等更加穩定，但采用DPM模型增加網格量會極大的減緩計算速度，由于本分不考慮流場內特別精細的結構，0.5 mm網格下的流場精度也能滿足要求。

綜上所述，1.0 mm網格計算精度很差，本文不采用，0.5 mm 網格下的精度和計算速度均能滿足要求，0.25 mm網格計算精度最高，但計算速度太慢，受實際計算資源所限，無法滿足大量工況的計算，因此本文數值計算均采用0.5 mm的網格。

可見，本文采用的數值方法可行，計算精度可以接受，后續的算例均采用該方法。

2 計算結果與分析

2.1 均勻性評價指標

為表征兩相燃燒室某截面的壓力分布均勻性，引入相對標準偏差CV，均勻性評價指數γP，克里斯琴森均勻系數CU，現分別介紹這些均勻性評價指標。

CV表征相對變異量的度量，為無量綱值，可以用來衡量均值顯著不同的總體離散性，也可用來比較流場均勻性的改善程度：

(1)

(2)

(3)

克里斯琴森均勻系數CU是克里斯琴森1942年提出的，是基于平均偏差的統計量，能夠直接反應被測參數與平均值的偏差程度，CU值越大，表明被測截面的壓力分布均勻性越好。其表達式為

(4)

為表征燃燒室出口溫度分布均勻性，現引入燃燒室出口溫度分布系數(Outlet Temperature Distribution Factor，簡稱OTDF)，其是衡量燃燒室出口溫度分布好壞的重要標志。它定義為燃燒室出口界面內最高燃氣總溫與燃氣平均總溫之差與燃燒室溫升的比值，定義式如式(5)所示，其中Tout_m為出口截面最高溫度，Tin_ave和Tout_ave分別為燃燒室進、出口平均溫度。

(5)

2.2 空氣噴注總溫

圖5給出了燃燒室出口處總壓的CU值、1-CV值隨著空氣噴注總溫的變化趨勢。由圖5可知，出口總壓的CU值、1-CV值均隨空氣噴注總溫的提高而提高，兩均勻性評價指標趨勢一致，均指明隨著空氣噴注總溫的提高，出口總壓的均勻性也隨之提高，CU值從0.363 2提高至0.477 8，提高了0.114 6，1-CV值從0.194 4提高至0.322 2，提高了0.127 8，提升數值接近。為探究空氣噴注總溫提高，出口壓力均勻性提高的原因，給出了不同噴注總溫下出口截面的總壓分布，如圖6所示。

圖5 不同空氣噴注總溫下燃燒室出口截面總壓CU值、1-CV值曲線

圖6 不同空氣噴注總溫下燃燒室出口截面處總壓分布曲線

由圖6可知，隨著空氣噴注總溫的提高，出口總壓峰值變化不大，但總壓上升至峰值的坡度卻不斷變緩，整體因此變得更加均勻。同時，單波模態下出口總壓存在兩個峰值，現給出了Case 1的總壓云圖，由圖7(a)可知：出口總壓的主峰值位于流場的滑移線附近，次峰值位于斜激波處。圖7(b)給出了該時刻下的氧氣組分分布，由圖可知，滑移線附近是未參與爆震反應的空氣的富集區，這是由于該工況為非均勻噴注，空氣不能充分參與爆震反應，剩余的空氣在斜激波的壓縮作用下流入滑移線附件區域。該區域溫度較低，而滑移線附近爆震產物流速是均勻變化的，因溫度降低而導致該區域馬赫數較高，如圖7(c)所示，滑移線附近區域為燃料室內馬赫數最高的區域，而高的馬赫數使得該區域燃燒產物總壓較大，即為出口總壓的主峰值區域。斜激波處的出口總壓較高是由于該處的反應產物受到了斜激波的壓縮作用，使其出口總壓增高。由圖6可知：無論是出口總壓的主峰值和次峰值，隨空氣噴注總溫的提高均變化不大，影響出口總壓均勻性的主要因素是主峰值之后，次峰值之前的燃燒產物出口總壓的均勻性。隨著空氣噴注總溫的提高，該部分燃燒產物出口總壓均勻性變好的主要原因可能是隨著空氣噴注總溫的提高，波前氣態辛烷含量更高且反應混合物活性更高，致使爆震波強度更高，相同條件下爆震產物的出口總壓更高，而出口總壓的峰值隨空氣噴注總溫的變化不大，因此出口總壓分布更加均勻。

圖 7 Case#1燃燒室內總壓、氧氣組分及馬赫數分布云圖

圖8為OTDF值以及出口平均總溫隨空氣噴注總溫的變化趨勢。由圖可知，隨著空氣噴注總溫從300 K提高至800 K，OTDF值從0.652 6降低至0.402 7，呈不斷下降趨勢，即燃燒室出口總溫分布越來越均勻。為分析出口總溫均勻性隨空氣噴注總溫提高而變好的原因，給出了不同空氣噴注總溫下燃燒室出口總溫分布，如圖9所示。由圖可知，不同空氣噴注總溫下的出口總溫分布曲線均存在一個峰值平臺及一個低洼，圖10給出了空氣噴注總溫300 K下的總溫云圖，兩圖對照可知，出口總溫分布曲線低洼處為滑移線附近區域，由圖7(b)可知，未參與爆震燃燒的空氣在斜激波的作用下在此區域匯集，導致該區域總溫較低，導致出口總溫分布曲線形成了一個低洼；出口總溫分布曲線的峰值平臺為斜激波與爆震燃燒產物相交的一塊三角形區域，爆震燃燒產物受斜激波的絕熱壓縮作用總溫迅速提高，致使出口總溫分布曲線出現了一個峰值平臺。圖8給出了出口總溫的平均值隨空氣噴注總溫的變化趨勢，可知隨著空氣噴注總溫的增加，出口總溫的平均值呈線性遞增的趨勢，由圖9可知，隨著空氣噴注總溫的提高，出口總溫的整體值均不斷提高，但出口總溫極大值與極小值的差值卻不斷降低，空氣噴注總溫300 K時極值差為1 711 K，空氣噴注總溫800 K時極值差降為1 233 K；同時，主峰值與次峰值的差值也在不斷降低，空氣噴注總溫為300 K、400 K、700 K、800 K下的出口總溫主次峰值差依次為649 K、540 K、411 K、245 K，隨著空氣噴注總溫的增加不斷降低。燃燒室出口總溫分布均勻性隨空氣噴注總溫的增加而增加的原因應為： 斜激波處三角形區域的總溫隨空氣噴注總溫的增加而增加，但增加幅度遠沒有平均總溫增加的幅度大；滑移線處低總溫區隨空氣噴注總溫的增加而增加，且增加幅度遠超平均總溫的增加幅度；斜激波處三角形區域的峰值總溫與滑移線邊緣處的次峰值總溫之差隨著空氣噴注總溫的提高而不斷減小。這些原因共同導致了提高空氣噴注總溫能有效提高出口總溫分布均勻性。

圖8 OTDF值及出口平均總溫隨空氣噴注總溫變化曲線

圖9 不同空氣噴注總溫下燃燒室出口總溫分布曲線

圖10 Case#1燃燒室內總溫云圖

2.3 辛烷射流間距

圖11給出不同辛烷射流間距下的燃燒室出口截面總壓CU值和1-CV值，由圖可知，兩均勻性評價指標趨勢一致，隨著辛烷射流間距增大，CU及1-CV均隨之增大，隨著辛烷射流間距從0 mm增加至8 mm，CU值從0.308 7增加至 0.446 7，增加了0.138，1-CV值從0.100 5增加至0.318 6，增加了0.218 1。兩指標均表征隨著辛烷射流間距增大，燃燒室出口總壓均勻性變好。圖12為不同辛烷射流間距下的燃燒室出口總壓分布，由圖可知，與上文分析類似，總壓分布曲線存在兩個峰值，分別對應著滑移線附近和斜激波附近的出口總壓，隨著辛烷射流間距的減小，總壓主峰值越來越大，曲線越來越陡峭，而次峰值變化較小，因此改變辛烷射流間距主要是影響主峰值，即滑移線附近的出口總壓而影響整個出口總壓的穩定性的。圖13給出了滑移線附近的馬赫數和總壓隨辛烷射流間距的變化趨勢，該區域的馬赫數和總壓也是整個出口截面的峰值馬赫數和峰值總壓，可以看到，峰值馬赫數和峰值總壓均隨辛烷射流間距增大而減小，兩者趨勢一致，在其余條件基本不變的情況下，馬赫數越高，總壓越大。因此，增大辛烷間距能改善出口總壓均勻性的原因應是：增大辛烷間距，滑移線附近馬赫數減小，導致總壓減小，而斜激波附近的總壓變化較小，進而使得整個出口截面的總壓均勻性增加。

圖11 不同辛烷射流間距下的燃燒室出口截面總壓CU值、1-CV值曲線

圖12 不同辛烷射流間距下的燃燒室出口總壓分布曲線

圖13 出口峰值馬赫數與峰值總壓隨辛烷射流間距變化曲線

圖14為OTDF值隨辛烷射流間距的變化趨勢。由圖可知，隨著辛烷射流間距從0提高至8 mm，OTDF值從0.583 5降低至0.703 0，呈不斷下降趨勢，即燃燒室出口總溫分布越來越均勻。由此可知，燃料噴注均勻性能影響燃燒室出口總溫均勻性，燃料噴注均勻性越好，燃燒室出口總溫均勻性越好。

圖14 OTDF值隨辛烷射流間距變化曲線

3 結論

采用均勻性評價指標CU，1-CV描述了空氣噴注總溫和辛烷射流間距對燃燒室出口總壓均勻性的影響，空氣噴注總溫和辛烷射流間距越大，燃燒室出口總壓均勻性越好。燃燒室出口截面總壓分布曲線存在著兩個峰值，分別對應著滑移線附近區域和斜激波附近區域。改變空氣噴注總溫，影響出口總壓均勻性的主要因素是主峰值之后，次峰值之前的燃燒產物出口總壓的均勻性；改變辛烷射流間距主要是影響主峰值，即滑移線附近的出口總壓進而影響整個出口總壓的穩定性。采用燃燒室出口溫度分布系數OTDF描述了空氣噴注總溫和射流間距對燃燒室出口總溫均勻性的影響，空氣噴注總溫越高，射流間距越小，燃燒室出口總溫均勻性越好，并得出了增大辛烷間距能改善出口總溫均勻性的原因應是：增大辛烷間距，滑移線附近馬赫數減小，導致總溫減小，而斜激波附近的總溫變化較小，進而使得整個出口截面的總溫均勻性增加。