高馬赫數(shù)下超聲速燃燒的自點火查表方法1)

張錦成 王振國 孫明波 汪洪波 王亞男 劉朝陽

(國防科技大學空天科學學院,長沙 410072)

引言

燃燒室是高超聲速沖壓發(fā)動機的核心部件,其內部超聲速氣流中燃料高效穩(wěn)定的燃燒一直是備受關注的課題[1].高馬赫數(shù)下,燃料的燃燒駐留時間極短難以充分預混,而燃料與氧化劑的摻混對點火及火焰穩(wěn)定的影響非常重要[2-3].化學反應和流動的強耦合性是超聲速反應流場的另一個關鍵問題.許多實驗和數(shù)值工作都以建立更準確精細的湍流化學反應模型為目標.本文借鑒火焰面的降維建庫查表思想,初步發(fā)展了具備描述高馬赫數(shù)自點火效應的查表方法.進一步與大渦模擬求解器耦合形成高馬赫數(shù)反應流求解器.有利于進一步理解自點火在火焰?zhèn)鞑ズ头€(wěn)定的作用并為高馬赫數(shù)工況沖壓發(fā)動機設計提供高效的仿真平臺.

高馬赫數(shù)下超聲速燃燒的自點火效應非常普遍,并且對燃燒室內火焰穩(wěn)定有重要作用[4-6].Cheng等[7]開展了一系列同軸射流實驗,觀察到了穩(wěn)定在下游的抬舉火焰.Boivin 等[8]應用氫氧三步簡化反應機理進行仿真,對比發(fā)現(xiàn)必須通過針對自點火效應的反應速率修正,可以得到與詳細機理較吻合的結果;捕捉流場自點火火焰基和描述湍流擴散火焰,并提出基于化學爆炸模型分析(CEMA)的自點火識別方法,結果表明上游“指狀”自點火核是流場最初反應階段.Moule 等[9]在文獻[7]實驗的基礎上進行數(shù)值仿真,對火焰穩(wěn)定機制進行分析.指出自點火效應和可壓縮效應共同控制化學反應過程,自點火過程通過生成水來開啟鏈式反應,下游激波后的高溫反應區(qū)主導釋熱,使抬舉火焰穩(wěn)定.Bouheraoua 等[10]進行相同工況的數(shù)值研究,發(fā)現(xiàn)流場中激波和點火的相互作用:入口激波提供了自點火發(fā)生的速度和混合條件,反應釋熱會在同一位置誘導出間歇的激波,其導致的壓力波動驅動自點火火核位置和強度產生不穩(wěn)定.Zhao 等[11]通過采用理想攪拌器模型考慮湍流影響結合化學反應建庫方法捕捉到了剪切渦中自點火火核,同時能夠對實驗的溫度場提供吻合很好的預測.

Ben-Yakar 等[12]開展了飛行Ma=10 的橫向射流實驗,高速紋影和OH-PLIF 結果顯示OH 基主要分布在射流迎風剪切層和近壁區(qū)域.Won 等[13]對Ben-Yakar 實驗構型和工況進行數(shù)值仿真,所用的DES (detached-eddy simulation)耦合有限速率模型,成功重現(xiàn)了高焓來流中橫向射流形成的三維非定常反應流場.計算結果說明源于上游回流區(qū)反向渦對(CVPs)相互作用持續(xù)地產生分離渦.Liu 等[14]也進行了相同條件的數(shù)值工作,采用高精度格式捕捉到了大尺度渦結構的產生和輸運;進行了化學反應特性分析,發(fā)現(xiàn)生成OH 基的鏈反應是吸熱的,并且發(fā)生在貧燃一側而釋熱明顯的鏈反應消耗OH 基同時生成水發(fā)生在富燃一側,總結了上游自點火維持下游擴散火焰的機制.Oliver 等[15]和Goldfeld[16]最新的工作揭示了自點火在橫向射流火焰時空演化中的作用.低釋熱的自點火發(fā)生在貧燃一側,主要分布在上游回流區(qū),釋熱的擴散導致自點火開始級聯(lián).發(fā)生在富燃側高釋熱的自點火,導致的膨脹擴展產生褶皺的火焰鋒面.

數(shù)值研究中,體現(xiàn)化學反應和湍流的交互,再現(xiàn)反應流場多時空尺度的演化非常困難.現(xiàn)有模型往往建立在一些假設上,只適用于一些條件下的仿真,許多針對性的模型改進由此開展.實際燃燒室的仿真模擬計算量大,所用燃料反應機理復雜,需要更加高效且適用于詳細化學反應機理的模型方法.化學反應的查表方法可以降低組分輸運帶來的計算量,同時可以使用詳細機理.

現(xiàn)有的處理湍流化學反應的方法主要有兩類:非降維方法和降維方法.典型的非降維方法為設定型PDF (assumed PDF)方法[17-18].該方法通過已知的高階矩信息設定概率密度函數(shù),可以給出反應速率在包含脈動量的狀態(tài)下的值.通過設定溫度高斯分布,組分聯(lián)合β分布[19-20]可以得到與實驗結果較為一致的仿真結果.但是當使用詳細化學反應機理時,設定型PDF 方法需要求解組分輸運方程,導致巨大的計算負荷.內在低維流形(ILDM)是最早也是最簡單的化學反應降維建庫方法,由Maas 和Pope[21]提出.由于ILDM 只采用了反應平衡后的狀態(tài)進行建庫,在反應預熱的低溫區(qū)域沒有被覆蓋,只能通過線性外推預測,會導致較大的偏差.擴展的火焰低維流形(FPI)方法[22]作為ILDM 擴展和改進被提出,能夠對低溫區(qū)反應得到更好的預測結果.FPI 使用了不同混合分數(shù)的一維層流預混火焰來描述低溫區(qū).自點火發(fā)生時反應速率非常低,特征時間大于局部流動特征時間,以上方法不能被準確描述非平衡過程.針對自點火現(xiàn)象,點火動力學查表(TKI)是一個被成功應用于低速壓燃發(fā)動機中的方法[23].TKI 采用一定初始溫度和壓力下的等壓預混反應器,在不同進度標量,生成燃料消耗速率和釋熱率數(shù)據(jù)庫[24].但是其在高超聲速反應流還沒有進一步的應用.經典的穩(wěn)態(tài)擴散燃燒火焰面模型引入薄反應層假設(Da? 1),建立火焰面方程.垂直于火焰面可以用混合分數(shù)標識,火焰面上可以用標量耗散率表示不同的狀態(tài),從而進行降維.穩(wěn)態(tài)擴散火焰面模型[25-26]不包含非穩(wěn)態(tài)燃燒分支,不能描述反應發(fā)生的瞬態(tài)過程.因此,火焰面進度變量模型被提出,增加一個定義的進度變量來標識反應進行的程度[27-28].但是自點火發(fā)生時,火焰面薄反應層假設不再成立,限制了火焰面方法在自點火效應明顯的高馬赫數(shù)反應流中的應用.

本文在火焰面進度變量模型的思路基礎上,提出在時間維度上降維以考慮自點火效應,建立一種超聲速燃燒的自點火查表方法.能夠與大渦模擬求解器耦合,適用于高馬赫數(shù)反應流的仿真應用.介紹了自點火數(shù)據(jù)庫的生成和映射過程,闡述流動求解器的耦合方式和參數(shù)更新方式,然后針對無流動的零維化學反應進行驗證.最后將耦合自點火查表方法的大渦模擬應用在兩種典型的高焓超聲速模型沖壓發(fā)動機中,證明該方法的有效性和適用性.

1 模型與方法

1.1 可壓縮大渦模擬方法

本文的模型實現(xiàn)和流場仿真工作以大渦模擬(LES)求解器為平臺.該求解器基于有限差分法求解超聲速反應流場.為保證數(shù)值方法的相容性和穩(wěn)定性,方程中不同數(shù)學物理性質的項采用不同的數(shù)值格式進行處理.對于空間離散項,采用五階WENO格式重構對流項,而二階中心格式離散黏性項.對于時間離散項,使用隱式LU-SGS 算法進行時間積分.有關方程和數(shù)值方法的詳細說明,請參見文獻[4,29].

1.2 自點火查表方法

借鑒于火焰面方法對火焰空間結構的降維,可以定義時間序列上一組熱力學狀態(tài)向量Vn.通過解如下的化學反應動力學常微分方程求得隨時間變化的熱力學參數(shù)Vn(如圖1 所示)

其中 ωi由基元反應模型和化學反應動力學公式[30]計算.圖1 展示了初始溫度1050 K,當量比為1,壓力100 kPA 的H2-O2反應曲線.通過在曲線上采集不同時刻的點,記錄諸如組分質量分數(shù),溫度,反應速率,釋熱率等信息組成一組原始數(shù)據(jù).計算不同初始狀態(tài)下的化學反應曲線,并進行數(shù)據(jù)采集則可以形成原始數(shù)據(jù)庫.可以通過初始狀態(tài)和不同初始狀態(tài)下的反應時間索引到每個時間序列上的所有數(shù)據(jù)點.考慮影響化學反應的主要因素,初始條件選定為混合分數(shù)Z,初始溫度T0,壓力p.初始條件設置的范圍和取值如表1,覆蓋了沖壓發(fā)動機燃燒室內涉及的大部分參數(shù)范圍,同時為了控制數(shù)據(jù)量對混合分數(shù)設置了不均勻的取值.

圖1 化學反應動力學計算出的溫度和H2 曲線Fig.1 Temperature and H2 curves calculated by chemical reaction kinetics

表1 初始參數(shù)的范圍和取值Table 1 Range and value of control parameters

值得注意的是,化學反應計算需要使用的參數(shù)隨時間變化是高度非線性的.圖2 展示了原始數(shù)據(jù)庫中產物H2O 反應速率和反應釋熱率隨時間變化的曲線,釋熱率和產物反應源項是化學反應與流動耦合的重要參數(shù);反應延遲時間,對流場中火焰分布位置和結構有重要影響.建立自點火模型的重點是:準確預測延遲時間,反應速率和釋熱率的極值;盡量還原反應速率和釋熱率隨反應進行的分布細節(jié).

圖2 H2O 的生成率和反應釋熱率Fig.2 Generation rate of H2O and reaction heat release rate

在流場仿真計算中,反應經歷時間和全局的物理時間是不一致的,通過其他方法也不容易求得,這給自點火查表耦合流動求解器帶來困難.通過組分構造一個隨時間單調變化的新變量,一方面可以代替時間進行數(shù)據(jù)庫索引,另一方面代替所有組分在流場中輸運,這是數(shù)據(jù)庫方法實現(xiàn)反應流動計算的關鍵.基于以上考慮在生成數(shù)據(jù)庫時要考慮將時間映射到一個新的變量,能夠通過組分反應反應進度.類似火焰面/進度變量模型,定義自點火進度變量

由于YC,max只取決于初始條件,這樣的歸一化是可行的.由于水的積累是單調的,所以歸一化過后的變量仍然保持對時間的單調性.進一步需要對進度變量的采樣點在點火延遲之前加密,以保證對自點火行為的精確捕捉.加密方式如圖3 所示,局部加密還可以減少數(shù)據(jù)庫的規(guī)模.

圖3 進度變量采樣點分布Fig.3 Distribution of value points for progress variables

圖4 展示了OH 基質量分數(shù)和溫度隨進度變量的變化,與圖1 對比可以看出進行映射的效果:反應過程在進度變量維度上變化的更均勻,這樣降低了化學反應在時間維度上的剛性.圖5 展示H2反應速率和進度變量反應速率曲線,進一步證明了選取進度變量的優(yōu)勢,反應速率曲線更光滑,能保證插值計算的精度.

圖4 OH 基質量分數(shù)和溫度隨進度變量的變化Fig.4 Variation of OH mass fraction and temperature with progress variables

圖5 H2 反應速率和進度變量生成率變化曲線Fig.5 H2 generation rate and progress variable generation rate curve

1.3 耦合查表與驗證

本節(jié)首先介紹數(shù)據(jù)庫的查表插值方法以及與大渦模擬求解器耦合方法.之后以零維的預混反應工況,驗證數(shù)據(jù)庫生成方式和查表插值方式.

生成的數(shù)據(jù)庫是一個四維 Φ (T0,p,Z,YC) 正交的大型數(shù)組,采用分段線性插值的方式進行查表輸出.如圖6 所示輸入當前的 (T0,p,Z,YC),通過查找并進行線性插值,得到進度變量速率,釋熱率hrr,組分Yi,然后進行溫度組分更新得到下一時間步的參數(shù).

圖6 數(shù)據(jù)庫查表插值示意圖Fig.6 Schematic diagram of lookup and interpolation from the auto-ignition database

在與大渦模擬耦合求解流場時,首先要在流場中輸運混合分數(shù)和進度變量用來代替組分的輸運,并模擬自點火火核的傳播和火焰擴散.輸運方程與火焰面/進度變量方法中構造的方程一致

其次定義數(shù)據(jù)庫查詢的輸入變量,根據(jù)求解方程組1 的定壓和守恒條件,從大渦模擬求解器輸入數(shù)據(jù)庫的 (T0,p,Z,YC) 可以如下定義:溫度根據(jù)當?shù)氐幕旌戏謹?shù)定義純混合溫度T0=Tox(1?Z)+TfZ,其中Tox,Tf分別為氧化劑測溫度和燃料一側溫度,Z為混合分數(shù).p,Z為當?shù)氐膲毫突旌戏謹?shù),YC由輸運方程求得.

然后通過查表插值后獲得的進度變量速率,釋熱率,組分,對流場進行更新.釋熱率,組分用于更新溫度,采用牛頓下山法進行迭代求解.進度變量速率用于更新進度變量,并添加到進度變量方程的源項上.

化學反應過程中遵循質量和原子守恒定律,如下定義的混合分數(shù)Z不隨時間變化

圖7 混合分數(shù)Z 守恒的驗證Fig.7 Verification of mixture fraction Z conservation

對初始溫度T=1200 K,壓力p=100 kPa,混合分數(shù)Z=0.0283 (當量比為1)工況進行計算,對比自點火數(shù)據(jù)庫插值方法和化學反應速率積分方法得到的結構來驗證組分求解的準確性.圖8 為插值和積分方式得到的H2和H2O 質量分數(shù)在反應中的變化,圖9 為不同初始溫度下兩種方式溫度隨時間和進度變量的變化.圖中紅色帶方塊標記的是查表進行插值得到的曲線,藍色圓圈標記的是反應速率積分的結果.積分結果可以當作驗證數(shù)據(jù)庫插值結果的基準.其中,作為進度變量的H2O 的質量分數(shù)吻合情況最好,說明進度變量的選取和進度變量速率的插值計算都是合適的.H2在快速反應階段吻合較好,在反應穩(wěn)定階段略偏低.總體上看自點火數(shù)據(jù)庫方法基本上可以實現(xiàn)化學反應的計算,可以成為有限速率方法的高效率替代方法.

圖8 組分H2 和H2O 質量分數(shù)在反應中的變化Fig.8 Variation of H2 and H2O mass fraction in the reaction

圖9 不同初始溫度下溫度隨時間變化Fig.9 Temperature variation with time at different initial temperatures

本小節(jié)主要驗證自點火數(shù)據(jù)庫方法的三個方面:混合分數(shù)守恒性,主要組分和溫度更新正確性,為該方法在超聲速反應流仿真中的驗證提供基礎.

2 算例驗證

當沖壓發(fā)動機燃燒室工作在高馬赫數(shù)工況下,入口來流空氣的總焓非常高.很多研究表明此時氫氣燃料能夠立即與熱空氣發(fā)生反應,無需借助任何穩(wěn)燃裝置能夠形成穩(wěn)定的火焰結構.在簡單的燃燒室構型產生的復雜射流燃燒流場,包含自點火與火焰之間的復雜的演化與競爭.本節(jié)基于 Gamba 的橫向射流燃燒實驗和Burrows-Kurkov 臺階射流實驗,針對高焓環(huán)境下的氫氣自點火主導的燃燒過程開展基于LES 的自點火查表方法的應用驗證.

2.1 高焓超聲速來流中橫向射流的自點火

Gamba 和Mungal[31]通過 OH 自發(fā)輻射和紋影技術研究了不同動量比(J=0.3～ 5.0)的聲速射流噴注進入超聲速來流誘導的混合燃燒特性.發(fā)現(xiàn)動量比J直接控制著下游的反應模式和分布區(qū)域.實驗條件下自由來流馬赫數(shù) 2.4,其滯止溫度可達到3000 K,相應的靜壓和靜溫分別為 40 kPa 和 1400 K.一個直徑 2 mm 的圓形噴注器安裝在距離平板前緣64 mm 位置.氫氣射流以聲速垂直噴注到超聲速流場中,其靜壓為 1074.5 kPa、靜溫為 250 K.計算域構型和參數(shù)如圖10,展向寬度為80 mm.

圖10 Gamba 實驗構型計算域示意圖 (單位:mm)Fig.10 Schematic diagram of the computational domain of the Gamba’s experiments (unit:mm)

燃料射流沿垂直壁面方向噴注到高焓超聲速氣流中,會形成非常復雜的反應流場.射流與來流之間的強相互干擾誘導出一系列激波和湍流渦等結構,在此過程中氫氣與熱空氣快速摻混并發(fā)生劇烈反應.

圖11 展示了實驗OH-PLIF 和計算得到的中心截面上OH 分布云圖.實驗和仿真均表明化學反應主要集中在兩個區(qū)域:一個靠近壁面,此處靜溫高點火延遲低;并且此處流動速度較低,燃料駐留時間長混合充分;因此可以維持燃燒穩(wěn)定.另一個反應區(qū)位于射流迎風剪切層,形成OH 基濃度很高的反應鋒面.盡管當?shù)厮俣群芨?但是在湍流渦的作用下自點火火核未被吹滅,并能夠輸運到流場下游.說明自點火查表方法能夠準確捕捉到流場中的自點火現(xiàn)象,并能夠反映自點火火核的傳播.

圖11 中心截面上OH 分布云圖Fig.11 OH distribution on the central section

2.2 臺階壁面燃燒室中射流誘導的抬舉火焰

Burrows-Kurkov 的實驗[32]采用了臺階壁面的燃燒室,燃料氫氣在臺階側壁以來流同向的速度注入由一定比例的H2O,N2,O2混合的污染空氣來流.臺階底壁略有擴展,具體的幾何構型見圖12.馬赫數(shù)2.4 的來流以1237.9 K 靜溫,96 kPa 靜壓在臺階上游矩形入口進入,根據(jù)一些基于此實驗的仿真研究入口處含有1 cm 厚度的邊界層;聲速的261.7 K的低溫射流以略高于主流的壓力噴注.來流和射流的參數(shù)見表2.

表2 Burrows-Kurkov 實驗射流和來流參數(shù)Table 2 Jet and inflow parameters in Burrows-Kurkov experiment

圖12 Burrows-Kurkov 實驗構型計算域示意圖 (單位:mm)Fig.12 Schematic diagram of the computational domain of the Burrows-Kurkov experimental configuration (unit:mm)

燃燒室內燃料和氧化劑剪切形成的渦的存在提供了混合的場所,自點火一般發(fā)生在剪切渦貧燃一側,局部溫度升高會降低點火延遲時間;自點火會發(fā)生級聯(lián)形成火焰基,向下游傳播形成穩(wěn)定的湍流火焰.由于存在點火延遲時間,火焰被抬舉維持在距射流出口一定位置.

圖13 為自點火數(shù)據(jù)庫查表法預測的中心截面上溫度和OH 基的分布,高溫區(qū)在上游成破碎的指狀,而在下游高溫區(qū)連片狀.這里由于上游的自點火火核的擴散和輸運形成了穩(wěn)定的火焰.氫氧基作為很重要的組分用來標記化學反應的鋒面,氫氧基出現(xiàn)的位置較反應升溫明顯的火焰區(qū)域提前,說明上游存在低溫自點火反應.自點火查表方法再現(xiàn)了自點火發(fā)展形成火焰的過程,證明其適用于火焰自點火共存的復雜燃燒流場,擴大了建庫類方法的應用范圍.但是細節(jié)上,本方法預測的火焰抬舉位置略低,可能是由于還沒有針對湍流對火焰的作用進行修正.

圖13 中心截面上溫度和OH 基的分布Fig.13 Distributions of temperature and mass fration of OH at the central slice

圖14 通過出口位置(x=356 mm)總溫剖面定量對比了實驗數(shù)據(jù)和仿真結果.可以看到在近壁的反應區(qū),總溫的上升與實驗數(shù)據(jù)吻合很好,遠離反應區(qū)的遠壁端總溫也基本與實驗一致.但是得到的反應區(qū)中心位置的總溫最大值略比實驗所得高.這是由于數(shù)據(jù)庫方法還是層流狀態(tài),使得預測的反映溫度偏高.定量的對比證實了自點火數(shù)據(jù)庫查表方法的可靠性.

圖14 出口位置(x=356 mm)總溫剖面Fig.14 Total temperature profile at the outlet (x=356 mm)

3 結論

超聲速高焓來流條件下,燃料室可能發(fā)生自點火行為對火焰穩(wěn)定有重要意義,并且給超聲速燃燒問題建模和仿真研究帶來了新的挑戰(zhàn).本文提出一種超聲速燃燒的自點火查表方法.首先介紹基于化學反應動力學建立自點火數(shù)據(jù)庫的思路并就零維預混火焰驗證了數(shù)據(jù)庫的建立和插值方法.然后對兩個三維射流誘導的燃燒流場進行仿真,對基于大渦模擬的自點火數(shù)據(jù)庫查表方法進行驗證,可以得到如下結論.

(1)選擇H2O 作為進度變量,降低了求解化學反應時的剛性.并且通過進度變量成功地耦合了流動求解和反應求解,可以正確地處理流動與化學反應的相互作用.

(2)算例驗證可以證明本文提出的方法初步具備對超聲速條件下包含自點火模式的反應流場進行仿真的能力.