支板增強混合超聲速燃燒的大渦模擬研究

黃志偉，何國強，秦 飛，魏祥庚

(西北工業大學 燃燒、熱結構與內流場重點實驗室，西安 710072)

?

支板增強混合超聲速燃燒的大渦模擬研究

黃志偉，何國強，秦 飛，魏祥庚

(西北工業大學 燃燒、熱結構與內流場重點實驗室，西安 710072)

基于開放源代碼軟件OpenFOAM，建立了三維可壓縮反應流動大渦模擬求解器，采用了PaSR亞格子燃燒模型和27步的氫氣-空氣反應動力學機理，開展了支板增強混合的超聲速燃燒大渦模擬研究，對比了滑移和粘性2種不同壁面邊界條件的影響。計算結果表明，不同截面上的平均軸向速度和溫度與實驗數據吻合良好，較好捕捉了超聲速擴散燃燒的火焰空間發展過程。詳細討論了剪切層增長、發展和破碎對燃燒過程的影響，揭示了支板后旋渦脫落與燃燒過程的耦合作用，區分了支板下游亞聲速區和超聲速區內不同的摻混模式。利用化學爆炸模式分析方法，獲取了爆炸化學過程及其特征時間尺度，得到了詳細的火焰結構及其穩定機制。

超聲速燃燒；化學爆炸模式分析；大渦模擬；PaSR燃燒模型；OpenFOAM

0 引言

雙模態超燃沖壓發動機和組合循環發動機是天地往返輸運系統和臨近空間飛行器的潛在高效動力系統。這類吸氣式發動機工作過程中涉及諸多基礎性科學問題。例如，高超聲速和高溫氣體動力學，超聲速條件下的點火與火焰穩定，復合防熱材料與主動熱防護，激波/附面層相互作用和吸熱型碳氫燃料技術等[1]。在這些前沿技術挑戰中，最基本的問題之一是燃燒室內超聲速流動條件下的湍流混合與高效燃燒。隨著計算機和數值計算研究的發展，具有較高時/空分辨率的湍流燃燒大渦模擬分析，成為超聲速燃燒和吸氣式發動機研究的有效手段[2]。

國內，清華大學的周建興等采用RANS方法和單步反應模型，通過隱式耦合求解可壓縮N-S方程，對氫氣的超聲速燃燒進行了三維數值模擬[3]。國防科技大學的范周琴等利用Flamelet燃燒模型和混合LES/RANS方法對支板構型超燃沖壓發動機進行模擬，氫氣動力學模型為9組分19方程，對比了三維與二維計算的區別[4]。中國科學院力學研究所的李曉鵬等分析了超聲速燃燒中的特征尺度及其影響因素，討論了Flamelet模型的適用性[5]。國外，Genin和Menon等使用LES-LEM和LES-EBU模型結合兩步氫氣動力學機理，對比研究了德國宇航院模型超燃沖壓發動機燃燒室，驗證了大渦模擬在超聲速流動與燃燒方面的潛在適用性[6]。Fureby等利用OpenFOAM平臺開展了超聲速燃燒大渦模擬方面的研究[7-9]，通過與大量試驗對比，驗證了該軟件在超燃沖壓發動機設計和超聲速燃燒機理研究方面的適用性。

國內外針對氫燃料超聲速燃燒的大渦模擬，大多采用Flamelet或EBU亞格子燃燒模型結合單步或簡化的動力學模型。本文基于近年來極具吸引力的開放源代碼計算平臺OpenFOAM，采用PaSR亞格子燃燒模型結合氫氣-空氣的詳細動力學機理，分析了湍流與詳細動力學的耦合作用，討論了超聲速燃燒流場火焰結構的演化，分析了壁面粘性對流場參數分布的影響。

1 數值方法與計算實現

1.1 OpenFOAM計算平臺與數值方法

OpenFOAM是由倫敦帝國理工學院開發的求解連續介質力學問題的開放源代碼軟件平臺[10]。其本質是基于面向對象程序設計的C++庫，具有優異的可移植性，包含豐富的物理和數學模型，可根據需要進行開發擴展，與其它軟件具有良好接口。近年來，已經在基礎研究和工程領域得到了廣泛的校驗和應用[7-9]。

本文的數值模擬基于OpenFOAM平臺，建立了適用于任意非結構化網格的可壓縮反應流動大渦模擬求解器。對所有對流通量的重建采用了基于非線性vanLeer限制器的二階精度TVD格式，對粘性和亞格子通量的空間離散化采用了中心差分格式。計算的最大CFL數不超過0.3，對應本文計算的物理時間步長為10-8s量級。當二階統計矩收斂時，認為求解收斂。

1.2 大渦模擬反應流的控制方程

1.3 亞格子流動與燃燒模型

(1)

(2)

該燃燒模型同時考慮了Kolmogorov尺度與求解尺度的綜合效應，避免了將湍流與火焰的解耦處理，較深刻地刻畫了湍流-化學相互作用過程，并得到了DNS和實驗數據的廣泛校驗，顯示了對于高速流動燃燒計算較優的預測能力[7-9,13]。

1.4 化學爆炸模式分析方法

典型反應流的微分方程可表達成如下離散形式[14]：

(3)

其中，D/Dt為隨流導數(物質導數);y為依賴變量矢量，如溫度、組分濃度等,對于空間離散化的流動方程，不同網格點上的化學組分濃度對應著不同的y入口值;ω表示化學源項;s包含所有的非化學源項，如火焰中的擴散、攪拌反應器中的均質混合等。當地化學信息完全包含在化學源項ω的Jacobian矩陣Jω中，化學模式即定義為Jω的特征模式，即一個特征值及其相應的左、右特征向量。如果化學模式特征值的實部λexp為正，則可進一步定義其為化學爆炸模式(CEM)[15]。因此，當非化學反應源項s在方程(3)中可忽略時，CEM的存在表征了混合物的本質是爆炸性的，即化學反應速率沿著其特征向量的方向呈指數型增長趨勢。

1.5 化學動力學模型

在湍流燃燒的大渦模擬中，選擇恰當的化學反應機理是關鍵的問題之一。為了更深入地研究反應機理與湍流的耦合作用，本文使用了根據文獻[16]中的氫氣詳細動力學機理得到的9組分、27步反應的氫氣-空氣化學動力學模型，在不過多增加計算量的同時，較好地反映氫氣燃燒的主要動力學過程，為獲得湍流火焰穩定機制提供了保證。表1給出了該動力學模型每個反應速率常數表達式中的參數。

表1 氫氣9組分、27步反應機理的速率常數Table1 Rate constants for 9 species 27-step chemical kinetics of hydrogen (cm3 ·mol·s·cal·K)

2 計算構型與邊界條件

圖1給出了本文所模擬的三維模型超聲速燃燒室尺寸示意圖[7]。計算網格包含590萬結構化的六面體單元，沿燃燒室上下壁面、支板壁面、支板前緣和尾跡區及剪切層區域進行了局部加密，以保證大渦模擬在這些區域內具有較好的網格分辨率。支板底端開有15個直徑為1 mm的小孔，每個小孔間隔2.4 mm。由于在寬度方向具有對稱性，同時為了節約計算成本，本文的計算區域只包含3個燃料噴孔。氫氣-空氣的等值比為0.034。

本文采用不同類型的壁面邊界條件，計算了2種工況。在空氣和燃料入口，由于流動是超聲速或聲速的，根據特征理論，所有的變量應用Dirichlet邊界條件；在燃燒室出口，采用Neumann邊界；在寬度方向使用周期性邊界條件。

圖1 模型超燃沖壓發動機燃燒室示意圖Fig.1 Schematic of the model scramjet combustor

在燃燒室的上下壁面及支板壁面上，第一個工況應用了無粘、絕熱邊界條件，第二個工況則應用了粘性壁面。

表2給出了空氣和燃料入口的狀態參數，來流空氣的湍流強度為0.5%，氫氣為5%。

表2 空氣來流與氫氣射流的入口條件Table2 Inflow conditions for air stream and hydrogen jet

3 結果與討論

3.1 LES結果與實驗數據的校驗

圖2給出了不同截面上軸向速度的分布曲線。

(a)x=0.120 m

(b)x=0.249 m

圖2(a)的位置非常靠近支板底端，由于燃燒釋熱，使得回流區內壓力升高，同時溫度升高，導致氣體粘性增大。因此，氫氣自支板噴出后受到十分明顯的減速作用。但由于氫氣出口速度很高，導致射流核心區相比于外圍氣體仍具有略高的速度，因而形成雙峰狀的速度分布。大渦模擬的結果在回流區內總體偏差較大，Menon和Fureby等指出，由于該位置恰好處于高度湍流的尾跡區，實驗中很難對速度進行精確地測量[6-7]，是造成偏差較大的原因之一。圖2(b)在距離支板很遠的下游位置，大尺度旋渦結構逐步轉化為小尺度漩渦，在燃燒釋熱膨脹和外層空氣流剪切力的雙重作用下，核心流逐步加速至接近主流的速度，同一截面上速度分布較為平緩。

圖3給出了不同位置上時間平均的溫度與實驗數據的對比曲線。圖3(a)中，反應較為劇烈，燃燒比較充分，釋放出大量的熱，形成高溫燃燒區，射流中心溫度達到局部最大值。圖3(b)中，由于氫氣和中間活性物質逐漸消耗完畢，反應強度減弱，流場中心溫度也有所降低。同時，由于分子的湍流輸運和輻射傳熱等作用，流場高溫區向上下壁面擴展，表現為曲線在高度方向上拓寬而峰值降低。無粘壁面的結果在第二個位置氣流核心區內偏差較大，可能由于假設了支板壁面無粘而導致燃料與空氣混合更充分，因而在燃燒室后段仍具有較強的熱釋放。

(a)x=0.167 m

(b)x=0.275 m

通過對比顯示，本文的大渦模擬計算結果與實驗數據測量結果吻合較好，有效地捕捉到了超聲速/亞聲速混合燃燒流場的演化過程。粘性壁面得到的速度和溫度分布準確性高于無粘壁面假設的結果。

3.2 無粘壁面與粘性壁面的結果對比

為了研究不同壁面邊界條件對流場結構發展、兩股氣流摻混及燃燒過程的影響，定義了標量場λ，它表征了流場不同的流動類型[17]：

(4)

其中

(5)

當λ=-1時，表示有旋流動；λ=0表示簡單剪切流動；λ=1表示平面拉伸流動。圖4(a)和(b)分別給出了粘性和無粘壁面條件下，燃燒室中心截面上λ的分布。考慮壁面粘性時，附面層內仍保持為簡單剪切流動，而附面層外迅速轉變成平面拉伸流動。同時，由于流體往下游發展，并不存在完全意義上的“純旋轉流動”，流體發生旋轉運動的同時，總伴隨著往下游的平動，因而并不存在λ=-1的情況。支板后的簡單剪切流動區域加寬，往下游發展時，表現出明顯向外擴展的趨勢。圖4(a)與圖4(b)的對比，清楚地表明了壁面粘性效應的影響范圍十分有限，當主要關心核心流的發展時，假設壁面無粘具有一定的合理性。

圖4 粘性與無粘壁面下λ分布的對比Fig.4 Distribution of λ for comparison of viscous and non-viscous wall

3.3 支板摻混增強與火焰穩定的機理分析

圖5給出了標量耗散率(Scalar Dissipation Rate)的等值面圖和中間活性物質OH質量分數的切面，以區分支板下游亞聲速和超聲速不同區域內的混合與燃燒模式。標量耗散率χ表征了燃料與氧化劑的混合速率，而OH基團可表征火焰面位置[18]。在近支板底端區域及初始反應混合層外邊界，由于兩股氣流之間存在較大參數梯度，混合主要通過擴散實現；而在更遠的下游，大尺度旋渦結構的卷吸與拉伸對摻混過程起決定性作用。在燃燒室末端，燃料與氧化劑已經達到了充分的混合，并實現了火焰的保持。從圖5可看到，該部分存在著大量的OH自由基。燃燒首先發生在支板尾跡區較薄的剪切層內，該區域聚集了大量的富燃燃氣和活性自由基，并對空氣和燃料進行預加熱，從而有利于反應剪切層內的燃燒與火焰穩定。

圖5 標量耗散率等值面與OH質量分數切面圖Fig.5 An iso-surface of scalar dissipation rate and contours of OH

來流空氣與氫氣射流作用形成反應剪切層，發生劇烈的湍流燃燒，并實現質量、動量和能量的交換。圖6給出了表征旋渦強度的速度梯度張量第二不變量 (Second Invariant of Velocity Gradient Tensor)等值面圖[19]。可看到，剪切層內含有豐富的旋渦結構。支板底部射流-尾跡區內大的旋渦結構往下游發展時，破碎成很多小旋渦。一方面，促進了反應物之間的摻混及隨后的燃燒；另一方面，又給熱釋放帶來較大的擾動，使火焰在三維空間內產生褶皺，并發生扭曲。支板前緣形成的兩道壓縮波以及支板末端形成的兩道膨脹波與燃燒室的上下壁面相互作用，形成復雜的反射波結構。由于超聲速反應剪切層較為穩定，且厚度較大，反射波無法直接穿透中心射流，從而在壁面與剪切層之間形成新的反射波系。燃燒釋熱還增加了回流區內的壓力，使得剪切層開始向外擴展，回流區擴大，有利于火焰的保持和穩定。

圖6 速度梯度張量第二不變量等值面與溫度切面圖Fig.6 An iso-surface of second invariant of velocity gradient tensor and contours of temperature

為進一步研究燃燒室內的湍流-化學相互作用過程，圖7給出了燃燒產物H2O作為混合物分數函數的散點圖分布，揭示燃燒室不同位置上摻混好壞與反應進行程度之間的關系。可看到，H2O的質量分數與混合分數隨著位置的后移表現出愈加明顯的線性關系，且分布更趨于集中化。在x=0.167 m位置處，相對較為分散的分布表明，該位置主要的反應活性強烈地受到湍流的影響，空氣與燃料通過大的湍流脈動強制進入到混合區的反應混合層內。摻混越充分，相應的燃燒反應越劇烈，產物H2O的質量分數也越大，但由于湍流脈動自身的隨機性，導致水的質量分數與混合物分數之間并不是嚴格的線性函數。而在更遠的下游位置，如x=0.199 m和x=0.275 m，隨著大尺度旋渦結構的破碎及湍動能的耗散，湍流的作用相對減弱，隨著位置的后移趨勢更加明顯。對于H2-O2燃燒系統，當量比混合物分數為Zst=0.028 5，在燃燒室后段，混合物分數已經低于該值，說明主燃燒反應區結束，在支板下游，實現了高效燃燒。

(a)x=0.167 m

(b)x=0.199 m

(c)x=0.275 m

圖8給出了放熱率作為混合物分數函數的散點圖分布。在x=0.167 m處，隨著混合物分數增加，放熱率首先緩慢上升，超過某一臨界值后迅速下降。更高的混合物分數意味著更充分的摻混，進而預示著更有效的燃燒。恰當的混合對于高效燃燒是必要的，但過大的混合速率會導致局部熱量和活性基團流失，化學反應反而受到抑制，甚至導致火焰不能自持穩定。在x=0.180 m處，放熱率隨混合物分數增加只保持緩慢增加，且分布極為分散，表明該位置放熱已趨于飽和，摻混對放熱率不起決定性作用。放熱率表現出強烈的脈動表明該位置的燃燒反應主要受湍流脈動影響。在x=0.199 m處，放熱率普遍較低，主放熱區結束。較低的混合物分數表明，H2已基本消耗完畢，中間活性物質對該位置的化學反應進程起主導作用。

(a)x=0.167 m

(b)x=0.180 m

(c)x=0.199 m

3.4 化學爆炸模式分析

利用化學爆炸模式分析方法，對每個網格節點進行了處理，圖9給出了化學源項的特征值，即爆炸模式時間尺度倒數的分布(對數坐標)。由于主流的流通時間不超過1 ms量級，故對時間尺度大于1 s的模式進行了截斷(該時間尺度相比于流通時間幾乎為“靜止的”，認為是非爆炸性的)。從圖9中可看到，整個燃燒室可劃分成兩個不同類型的區域：即藍色背景的非爆炸區域和彩色標識的核心流爆炸區域。其中，非爆炸區域又包含3個子區域：氫氣核心射流外圍的“貧燃”來流空氣流動區，氫氣核心射流緊鄰出口下游的低溫“富燃”湍流區，以及富含平衡燃燒產物的上、下混合層區域內。在前兩個子區域內主要是因為混合物超過了其貧/富可燃極限。

在近支板底端區域，由于主要通過擴散來實現摻混，因此在非爆炸和爆炸混合物之間存在較為顯著的界限，CEM特征時間尺度表現出較大的梯度。與圖6中相對低溫的高溫混合物相對應的兩個混合層外邊界將爆炸混合物與非爆炸混合物分開，即在該層內形成了所謂的預混燃燒模式。在燃燒室的下游位置，主要通過大尺度旋渦的卷曲拉伸來實現混合。相對而言，爆炸/非爆炸混合物的混合更為均勻，二者間形成較厚的中間過渡剪切層。

圖9 當地化學爆炸模式時間尺度lgλexp的空間分布Fig.9 Spatial distribution of the time scale of local chemical explosive mode

4 結論

(1)跨越支板回流區后大尺度旋渦的破碎對燃燒過程具有雙重作用：一方面，由于形成大量的小旋渦結構，有效促進了燃料-氧化劑的摻混，進而促進了化學反應進程；另一方面，也給熱釋放率帶來很大擾動。PaSR亞格子燃燒模型能捕獲超聲速燃燒流場的復雜物理化學過程。

(2)在近支板底端區域內的摻混主要通過擴散來實現；而在較遠的下游位置，湍流混合主要是由大尺度旋渦結構的卷吸與拉伸作用來完成。通過化學爆炸模式分析，獲得了燃燒室火焰結構及其穩定機理。

(3)在支板下游不同位置處，放熱速率存在不同的影響因素：隨著離開支板軸向距離的增加，首先是燃料與氧化劑的摻混，然后是湍流脈動，再者是中間活性物質，分別對放熱率起主導作用。

(4)粘性壁面相比于無粘壁面最主要的影響在壁面附面層很窄的區域內，對核心流場的發展影響很小。在主要關心核心流的參數變化時，假定壁面無粘具有一定的合理性。

(5)大渦模擬的平均軸向速度和溫度均與試驗數據吻合較好，表明OpenFOAM求解器有較高的數值計算精度。其開放源代碼的特性為后續開展燃燒流場計算與發動機設計、多學科優化和不確定性分析方法的耦合提供了充分的拓展空間。

[1] 蔡國飚, 徐大軍. 高超聲速飛行器技術[M]. 北京:科學出版社, 2012.

[2] Clark Ryan J, Shrestha SO Bade. A review of numerical simulation and modeling of combustion in scramjets[J]. Proceedings of the Institution of Mechanical Engineers, Part G: Journal of Aerospace Engineering, 2014:1-23.

[3] 周建興, 樸英, 豈興明, 等. 超聲速燃燒室內氫氣燃燒的三維數值研究[J]. 航空動力學報, 2007, 22(12):1984-1988.

[4] 范周琴, 孫明波, 劉衛東. 支板噴射超聲速燃燒流場三維大渦模擬[J]. 國防科技大學學報, 2011, 33(1):1-6.

[5] 李曉鵬, 張泰昌, 齊力, 等. 超聲速燃燒中的特征尺度及影響因素[J]. 航空動力學報, 2013, 28(7):1458-1466.

[6] Genin Franklin, Menon Suresh. Simulation of turbulent mixing behind a strut injector in supersonic flow[J]. AIAA Journal, 2009, 48(3):526-539

[7] Fureby C, Fedina E, Teger J. A computational study of supersonic combustion behind a wedge-shaped flameholder[J]. Shock Waves, 2014, 24(1):41-50.

[8] Fureby C, Chapuis M, Fedina E, et al. CFD analysis of the HyShot II scramjet combustor[J]. Proceedings of the Combustion Institute, 2011, 33(2):2399-2405.

[9] Berglund M, Fedina E, Fureby C, et al. Finite rate chemistry large-eddy simulation of self-ignition in a supersonic combustion ramjet[J]. AIAA Journal, 2010, 48(3):540-550.

[10] www.openfoam.com[DB/OL].

[11] Poinsot T J, Veynante D. Theoretical and Numerical Combustion, Third Edition[M]. US: R. T. Edwards, 2012.

[12] Bensow R E, Fureby C. On the justification and extension of mixed models in LES[J]. Journal of Turbulence, 2008, 54(8):1-17.

[13] Sabelnikov V, Fureby C. Extended LES-PaSR model for simulation of turbulent combustion[J]. Progress in Propulsion Physics, 2013, 4(1):539-568.

[14] Luo Z Y, Yoo C S, Richardson E S, et al. Chemical explosive mode analysis for a turbulent lifted ethylene jet flame in highly-heated coflow[J]. Combustion and Flame, 2012, 159(1):265-274.

[15] Lu T F, Yoo C S, Chen J H, et al. Three-dimensional direct numerical simulation of a turbulent lifted hydrogen jet flame in heated coflow: a chemical explosive mode analysis[J]. Journal of Fluid Mechanics, 2010, 652(1):45-64.

[16] Marinov N M, Westbrook C K, Pitz W J. Detailed and global chemical kinetics model for hydrogen[C]//8th International Symposium on Transport Processes, 1995.

[17] http://openfoamwiki.net/index.php/FlowType[DB/OL].

[18] Hawkes, Evatt R, Sankaran Ramanan, et al. Direct numerical simulation of turbulent combustion: fundamental insights towards predictive models[J]. Journal of Physics: Conference Series, 2005, 16(1):65-79.

[19] Khashehchi M, Elsinga G E, Bornstein H, et al. Studying the invariants of the velocity-gradient tensor of a round turbulent jet using Tomo-PIV[C]//8th International Symposium on Particle Image Velocimetry, 2009.

(編輯：崔賢彬)

Large eddy simulation of strut enhanced mixing for supersonic combustion

HUANG Zhi-wei, HE Guo-qiang, QIN Fei, WEI Xiang-geng

(Science and Technology on Combustion, Internal Flow and Thermal-Structure Laboratory, Northwestern Polytechnical University, Xi'an 710072, China)

Large Eddy Simulation (LES)of supersonic combustion in a model scramjet combustor based on an Open Source Field Operation and Manipulation (OpenOAM)computing platform was established, with two different wall boundary conditions, i.e. slip and viscous walls applied. The three-dimensional LES solver, which adopts a Partially Stirred Reactor (PaSR)sub-grid combustion model along with a skeleton 27 steps hydrogen chemical kinetics, was used to study strut-enhanced mixing and combustion. LES results show that mean axial velocity and temperature at different cross sections match well with experimental data, and spatial evolution of the supersonic diffusion flame is well captured. Effects of shear layers growth, development and breaking down on combustion processes were discussed in detail, and the coupling effects with vortex shedding at the strut base were revealed. Different mixing modes were recognized after the strut where subsonic and supersonic flows co-exist. Explosive chemical processes and their characteristic time scales were acquired by the Chemical Explosive Mode Analysis (CEMA)method, and the detailed structure and stabilization mechanism of the flame was identified.

supersonic combustion；chemical explosive mode；large eddy simulation；PaSR combustion model；OpenFOAM

2014-10-17；

：2015-03-26。

黃志偉(1989—)，男，博士生，研究方向為航空宇航推進理論與工程。E-mail:huangzhiwei504@mail.nwpu. edu.cn

V435

A

1006-2793(2015)05-0664-07

10.7673/j.issn.1006-2793.2015.05.012