高振勛,蔣崇文,李椿萱

(北京航空航天大學 國家計算流體力學實驗室,北京 10019)

0 引 言

超燃沖壓發(fā)動機的發(fā)展令超聲速湍流擴散燃燒的數(shù)值模擬技術受到了高度重視[1-3]?；诶字Z平均(RANS)以及大渦模擬(LES)方法對湍流燃燒流動進行數(shù)值模擬時,不可解尺度上發(fā)生的湍流與燃燒化學反應相互作用過程需要湍流燃燒模型考慮[4-7]。事實上,目前的湍流燃燒模型主要考慮的是湍流脈動效應對燃燒過程的影響,而燃燒放熱效應對湍流的影響則主要通過可解尺度上流場參數(shù)的變化來反映[8-11]。目前國外已發(fā)展了多種湍流燃燒模型,包括火焰面模型[12-13]、條件矩模型[14-15]、線性渦模型[16]、輸運PDF模型[17]等,而本文主要關注其中的火焰面模型。

火焰面模型是基于燃燒化學反應時間尺度足夠小的假設而建立的,該模型并不對組分輸運方程中的時間平均(RANS)或濾波后(LES)化學反應源項進行封閉,而是直接對反應標量(組分濃度和溫度)進行?；Ｔ撃Ｐ屠没旌戏謹?shù)變量,將不可解尺度上的火焰內(nèi)部結構模擬與可解尺度上的湍流混合模擬從形式上解耦?；鹧鎯?nèi)部結構可通過求解在混合分數(shù)空間內(nèi)基于漸近分析得到的火焰面模型方程得到[18],湍流混合過程則由混合分數(shù)和標量耗散率的聯(lián)合概率密度函數(shù)(PDF)描述。最后,通過PDF 積分則可完全確定湍流燃燒過程中組分濃度和溫度的平均特性。由于使用該模型時不需要求解組分的輸運方程,對于大規(guī)?；瘜W反應系統(tǒng)的模擬該模型在計算量上具有明顯優(yōu)勢。由于一般穩(wěn)態(tài)火焰面模型中只利用火焰面模型方程的穩(wěn)態(tài)解,因此無法考慮自點火和熄火過程,因此之后又進一步發(fā)展了火焰面/進度變量模型[19-21],該模型中所使用的火焰面數(shù)據(jù)庫同時包括了火焰面方程的穩(wěn)態(tài)和非穩(wěn)態(tài)解,并且其控制變量變?yōu)榛旌戏謹?shù)和進度變量(選定的某組分濃度),同時該模型需要求解被選定作為進度變量組分的輸運方程。總的來說,火焰面類模型在一定條件下才適用,其突出特點是計算量并不隨化學反應系統(tǒng)增加而明顯增大。

火焰面類模型最初均是針對低速湍流燃燒流動而建立的,模型中很多假設在超聲速可壓縮流動中并不成立,因此不能直接應用于超聲速湍流燃燒流動。國內(nèi)外部分學者已經(jīng)開展了火焰面類模型向超聲速湍流燃燒流動推廣的相關研究,本篇文章主要目的是對該方面的研究進展進行綜述,其中重點介紹本文作者所開展的相關工作。

1 火焰面類模型理論

1.1 火焰面模型

描述火焰內(nèi)部結構的火焰面模型方程是在一定物理假設下得到的,其物理圖像如圖1所示。該模型假設燃燒的化學反應足夠快,局部瞬時化學反應層的厚度(lr)小于湍流的Kolmogrove尺度(ηk),在該尺度下組分濃度和溫度的湍流脈動已由分子黏性所耗散,因此化學反應層內(nèi)部不受湍流脈動的影響,只由分子擴散控制。另外,大尺度湍流運動的應變率可引起化學反應層空間上的褶皺,但由于化學反應層厚度遠遠小于褶皺的曲率半徑(約等于湍流積分尺度,如圖1中的lt所示),以致組分濃度沿化學反應層切向的變化與沿法向的變化相比可以忽略。在以上假設基礎上所建立的火焰面模型相當于將湍流火焰看作是嵌入在湍流流場內(nèi)薄的具有一維結構的火焰面系綜。

圖1 火焰面模型的物理圖像Fig.1 Physical image of the flamelet model

通過引入擴散燃燒中的重要物理量混合分數(shù)Z,可將瞬時組分濃度輸運方程變換到混合分數(shù)空間,以混合分數(shù)空間內(nèi)的化學反應層厚度作為小參數(shù)進行奇異攝動分析[18],忽略非定常效應,可得到火焰面內(nèi)組分濃度的控制方程如式(1)所示,其中Ys為組分s的質(zhì)量分數(shù),ρ 為密度,χ 為標 量耗散率,為組分s的化學反應源項。另外,在低馬赫數(shù)條件下,忽略靜焓h 瞬時輸運方程中的黏性耗散項及壓力項,并在Lewis數(shù)為1的假設下,可得到瞬時靜焓h 與混合分數(shù)Z 的方程相似,即式(2)中的h 和Z 的線性關系。最后,在低馬赫數(shù)條件下假設整個燃燒空間為常壓空間,即式(3)中的p 為常數(shù)。于是,可最終導出火焰面模型的控制方程為:

其邊界條件為:

其中上標F 和O 分別表示燃料和氧化劑來流中的變量值。由于直接考慮火焰面模型方程中的非穩(wěn)態(tài)效應難度較大,因此一般在應用時忽略式(1)中的時間相關項,此時稱為穩(wěn)態(tài)火焰面模型。式(1)中的標量耗散率χ 是隨Z 變化的參數(shù),一般采用當量混合分數(shù)Zst對應的χ 值,即當量標量耗散率χst,作為湍流混合場對火焰內(nèi)部結構施加作用的控制參數(shù)。在給定不同的當量標量耗散率χst值條件下,求解火焰面模型方程(1)～(3)可得到Ys(Z,χst),T(Z,χst),即建立所謂的火焰面數(shù)據(jù)庫。之后,流場內(nèi)的平均組分濃度和溫度可由PDF(PDF形式采用預設方式得到)積分獲得:

基于隨機變量Z 和χst分別與湍流的大尺度和小尺度相關的觀點,可假設兩者統(tǒng)計獨立,則可表示為兩個邊緣PDF的乘積:

1.2 火焰面/進度變量模型

事實上,完整的火焰面模型方程的解包括穩(wěn)態(tài)和非穩(wěn)態(tài)兩個分支,即所謂的“S型曲線”。以上穩(wěn)態(tài)火焰面模型中以χst為控制參數(shù)建立的火焰面數(shù)據(jù)庫只包含了穩(wěn)態(tài)解,不能得到完整描述最大火焰溫度與當量標量耗散率χst關系的整個“S型曲線”,從而無法模擬自點火和熄火過程?；鹧婷?進度變量(Flamelet/Progress Variable,FPV)模型[19-21]的提出主要是解決火焰面模型對自點火模擬精度不準確的問題,方法是引入一個可描述化學反應進程的進度變量C,以進度變量C 作為火焰面方程組的控制參數(shù)。C 可定義為某種反應物或生成物的組分濃度,由此可得:

由于C 在χst一定的情況下與Z 一一對應,而其他組分濃度以及化學反應源項也與Z 一一對應,因此可以建立其他組分濃度以及化學反應源項與進度變量C 的一一對應關系,即建立包含Ys(Z,C)和關系的數(shù)據(jù)庫。由此建立的火焰面數(shù)據(jù)庫可以描述完整的“S型曲線”,從而準確描述自點火和熄火過程。

在獲得火焰面數(shù)據(jù)庫的基礎上,可得到平均的組分濃度和溫度在物理空間中的分布:

混合分數(shù)和進度變量的聯(lián)合概率密度分布為:

因此,還需要求解進度變量C 的輸運方程,其中的進度變量源項的計算公式為:

2 火焰面模型的可壓縮修正

要將火焰面模型用于超聲速燃燒流動中,首先面臨的一個問題是:火焰面類模型中的基本假設——燃燒化學反應的時間尺度足夠小而令局部瞬時化學反應層的厚度小于湍流的Kolmogrove尺度——在超聲速燃燒流動中是否能夠滿足? Bray等[25]分析認為超燃沖壓發(fā)動機燃燒室的燃燒滿足火焰面的假設。Balakrishnan等[26]在湍流積分尺度等于燃燒室尺寸以及湍流度為50%的假設下,估計了來流馬赫數(shù)1～25范圍內(nèi)超燃沖壓發(fā)動機燃燒室內(nèi)的湍流時間尺度的變化,并利用氫氣/空氣動力學模型估計了化學反應時間尺度,分析認為大部分情況下燃燒室內(nèi)超聲速湍流燃燒流動可滿足火焰面假設。Waidmann 等[27]對德國宇航局(DLR)氫燃料超燃沖壓發(fā)動機實驗數(shù)據(jù)的分析也顯示其燃燒處于火焰面模式。針對這一問題,國內(nèi)Gao等[28]和Fan等[29]分別開展了典型氫燃料超聲速燃燒室內(nèi)湍流和燃燒特征時間尺度的分析,利用Cuenot給出的判據(jù)[30]檢驗火焰面假設在超聲速燃燒流動中是否能夠滿足。其中,Gao等的分析基于RANS模擬,而Fan等則利用混合RANS/LES模擬,分析結果均顯示,火焰面模型的基本假設在所研究的超聲速燃燒場中均成立。然而,針對碳氫燃料的超聲速燃燒流場中火焰面假設成立情況的研究,目前尚未見到公開文獻報道。

在火焰面基本假設滿足的條件下,火焰面模型在超聲速燃燒流動中的應用還面臨其他問題。原始的火焰面模型在使用時,流場的平均溫度與組分濃度一樣,均是由火焰面數(shù)據(jù)庫的PDF積分獲得,這種方法在超聲速流動中無法考慮溫度與速度的強耦合,也即動能和內(nèi)能之間的轉(zhuǎn)化。Zheng和Bray[31]最早通過構造基于經(jīng)驗的平均溫度修正公式以考慮超聲速流動中動能向內(nèi)能的轉(zhuǎn)化。Sabel’nikov等[32]在此基礎上根據(jù)條件平均的思想,引入了速度關于混合分數(shù)的條件平均,從而將動能變化對火焰內(nèi)部結構的影響包含到火焰面數(shù)據(jù)庫中。以上兩種修正模型中溫度仍然采用火焰面數(shù)據(jù)庫積分獲得,無法有效考慮超聲速流動中速度和溫度的強耦合,故未得到廣泛應用。Overmann[33]首次將火焰面模型應用于實際構型超燃沖壓發(fā)動機燃燒室的二維數(shù)值模擬中,提出了僅利用火焰面數(shù)據(jù)庫PDF 積分獲得組分濃度,而溫度由能量方程隱式得到的策略,使得基于低馬赫數(shù)假設建立的火焰面模型可應用于復雜的超聲速流動中。該策略之后被其他一些學者采用,包括基于火焰面模型的RANS模擬[34-36]以及LES模擬[37-40]。

然而,Overmann的策略在建立組分濃度的火焰面數(shù)據(jù)庫過程中并未計及流動中的高馬赫數(shù)可壓縮效應,將基于該數(shù)據(jù)庫獲得的平均組分濃度直接用于反解溫度的方法在理論上是不完備的。在高馬赫數(shù)條件下分析1.1節(jié)中的火焰面模型理論:一方面,式(3)中將燃燒空間內(nèi)的壓強視為常數(shù),而在超燃沖壓發(fā)動機燃燒室流場中,由于存在復雜的波系,壓強的空間變化明顯;另一方面,式(2)中h 和Z 的線性關系在高馬赫數(shù)流動中顯然不能成立,此時動能變化引起的靜焓變化將對火焰內(nèi)部產(chǎn)生影響。

為了評估壓強的空間變化在火焰面內(nèi)部結構的影響,Gao等在文獻[41]中針對德宇航(DLR)的超聲速燃燒室[27]內(nèi)典型流動條件(如表1所示),在不同的p0值(1～1.6)×105Pa條件下求解了火焰面方程(1～3),結果由圖2 給出?？梢钥吹?在基準壓強1×105Pa基礎上增大60%后,混合分數(shù)空間內(nèi)的主要組分H2、O2和H2O 濃度變化非常小,且變化主要集中在Zst附近;而對于中間產(chǎn)物H、O 和OH 等,圖2中顯示壓強變化可誘導出較明顯的變化。接下來,將壓強值取1.6×105Pa時獲得的火焰面數(shù)據(jù)庫用于DLR 的燃燒流場模擬,與之前使用基準壓強火焰面數(shù)據(jù)庫所獲得溫度結果進行了對比,如圖3所示?？梢钥吹?之前p 取1×105Pa火焰面數(shù)據(jù)庫獲得的溫度峰值低于p 取1.6×105Pa時火焰面數(shù)據(jù)庫的結果,但誤差在3%以內(nèi)。可以預想,若生成火焰面數(shù)據(jù)時采用的壓強值高于實際燃燒流場壓強值時,所獲得的溫度結果將偏低。

表1 DLR 超聲速燃燒室實驗算例條件Table 1 Experimental conditions for DLR scramjet

圖2 不同壓強取值條件下求解火焰面方程獲得的組分濃度對比Fig.2 Species mass fraction distributions in the Z space under different values of p

圖3 不同壓強取值條件下求解火焰面方程獲得的溫度分布Fig.3 Temperature distributions obtained by different sets of flamelet library with different values of p

以上的研究發(fā)現(xiàn),火焰面數(shù)據(jù)庫中的主要組分濃度以及所獲得的流場溫度對數(shù)據(jù)庫建立時給定的環(huán)境壓強取值并不敏感,在實際應用時可以選擇燃燒室內(nèi)的典型壓強值作為火焰面數(shù)據(jù)庫建立時的環(huán)境壓強值。然而,若關心中間產(chǎn)物的濃度,則環(huán)境壓強的取值則會產(chǎn)生明顯影響。還有一點需要說明,如果實際燃燒室內(nèi)的壓強變化范圍很大,則壓強的變化也會對數(shù)據(jù)庫中的主要組分濃度產(chǎn)生顯著影響,此時則需要在火焰面數(shù)據(jù)庫中考慮壓強的影響。

圖4 DLR 燃燒室內(nèi)的溫度云圖Fig.4 Temperature contour in the DLR combustor

圖5 燃燒流場中不同站位的平均 和的分布對比Fig.5 Distributions of andat different streamwise locations

另一方面,分析高馬赫數(shù)可壓縮效應在火焰面方程中的另一個重要影響:式(2)中的瞬時靜焓h 與混合分數(shù)Z 的線性關系將不再成立。Gao 等在文獻[41]中利用DLR 燃燒室流場定量評估了這一誤差的影響。首先,提取了DLR 燃燒室流場(如圖4所示)不同流向站位上的、用空氣來流和氫氣射流靜焓值無量綱化后的分布,見圖5。若瞬時的Z 和h滿足線性關系,則應與相等,但圖5中顯示,不同站位上可以發(fā)現(xiàn)兩條曲線差別明顯,說明Z 和h的關系嚴重偏離線性關系。要定量評估這一誤差的大小,需要在火焰面模型方程中建立真實的瞬時h和Z 的關系,進而建立相應火焰面數(shù)據(jù)庫,將計算的溫度結果與之前的結果對比。然而,超聲速流動中瞬時h 和Z 關系并無確定的理論形式,因此以上做法無法實現(xiàn)。Gao等借鑒了Sabel’nikov等[32]提出的h 和Z 的二次函數(shù)關系。除了氧化劑和燃料來流確定的兩個邊界條件關系外,再提取每個站位上最大值與的對應關系,構造了不同站位上不同的h和Z 二次函數(shù)關系來生成不同的數(shù)據(jù)庫。如圖6給出了按以上方法獲得的四個站位上的h 和Z 二次關系曲線。

圖6 不同站位上h 和Z 的近似二次函數(shù)關系Fig.6 Different h-Z quadric relations for different streamwise locations

圖7 為不同的h 和Z 二次函數(shù)關系獲得的火焰面數(shù)據(jù)庫中典型組分(H2O 和OH)濃度在混合分數(shù)空間內(nèi)的分布,可以看到主要組分的濃度差別很小,而中間產(chǎn)物濃度則有明顯不同。應用以上火焰面數(shù)據(jù)庫分別獲得不同站位上的溫度分布(記為flamelet model h-Z),并與之前的溫度結果進行對比,如圖8所示?？梢钥吹?雖然圖5中顯示由于可壓縮效應導致的與的偏離非常大,但圖8中的數(shù)值實驗結果卻顯示此偏離所引起的火焰面模型獲得的溫度差別非常小,溫度峰值的偏差最大約為0.5%(x＝108 mm)。特別地,x＝108 mm 站位上,flamelet model h-Z 得到的溫度峰值低于火焰面模型結果,表明在動能向靜焓轉(zhuǎn)化的流場區(qū)域未考慮可壓縮效應的火焰面模型獲得的溫度場將偏高;x＝300mm 站位上,flamelet model h-Z 得到的溫度峰值均高于火焰面模型結果,表明在靜焓向動能轉(zhuǎn)化的流場區(qū)域未考慮可壓縮效應的火焰面模型獲得的溫度場將偏低。由此可見,式(2)中h 和Z 線性關系在超聲速流場中誤差很大,但該效應對于火焰面模型溫度場的模擬的影響很小。

圖7 不同的h 和Z 關系下求解火焰面模型方程獲得的組分濃度分布Fig.7 Mass fraction distributions in the Z space obtained under different h-Z relations

圖8 考慮動能與靜焓相互轉(zhuǎn)化后的火焰面模擬獲得的溫度分布Fig.8 Temperature distributions using different sets of flamelet library under different h-Z relations

為了進一步評估流場可壓縮效應對火焰面模型模擬的影響,高振勛等[42]建立了一種交互式的火焰面模型(Representive Interactive Flamelet model,RIF),即將流動方程與火焰面方程耦合求解,通過為燃燒區(qū)域的每點建立與當?shù)亓鲌鰠?shù)相匹配的火焰面數(shù)據(jù)庫,將流場當?shù)貕簭娮兓约癶-Z 關系的變化直接引入到火焰面數(shù)據(jù)庫中。利用RIF 模型對DLR 燃燒室的計算結果表明在該算例條件下,其與傳統(tǒng)火焰面計算結果差別很小,表明該流動中可壓縮效應對火焰面模型模擬結果影響很小。

3 火焰面/進度變量模型的可壓縮修正

火焰面/進度變量模型所采用的數(shù)據(jù)庫中包含了組分濃度和進度變量化學反應源項數(shù)據(jù),該模型應用于超聲速流動除了要考慮上文中流場當?shù)氐目蓧嚎s效應對數(shù)據(jù)庫中組分濃度的影響,還需要考慮進度變量化學反應源項受可壓縮效應的影響。Saghafian等[43]展示了火焰面數(shù)據(jù)庫中的進度變量化學反應源項值對于混合分數(shù)空間內(nèi)靜焓分布變化非常敏感,因此在超聲速流動中使用火焰面/進度變量模型時必須對數(shù)據(jù)庫中化學反應源項進行可壓縮修正。由于火焰面/進度變量數(shù)據(jù)庫的建立很復雜,因此在建立數(shù)據(jù)庫的過程中考慮流場當?shù)乜蓧嚎s效應非常困難。Saghafian[43]提出了另外一種思路:在建立數(shù)據(jù)庫的時候仍然采用給定的壓強值和靜焓分布,而對PDF積分后的進度變量源項值構造基于流場當?shù)刈兞康目蓧嚎s標度,其提出的可壓縮標度形式為:

高振勛等將其提出的可壓縮火焰面/進度變量模型應用于Cheng 的超聲速軸對稱射流燃燒實驗算例[44]的RANS模擬,其溫度場結果如圖9所示。

圖9 FPV模型所得Cheng算例的溫度場Fig.9 Temperature contour using FPV model for the Cheng case

圖10 式(13)可壓縮標度系數(shù)沿對稱軸的分布Fig.10 Distributions of the rescaling coefficient by Eq.(13)along axis

圖11 使用/不使用可壓縮標度的FPV模型得到對稱軸溫度分布對比Fig.11 Distributions of the temperature along axis using FPV model with and without the compressible rescaling

圖10 展示了利用式(13)計算得到的進度變量源項可壓縮修正系數(shù)在對稱軸上的分布,可以看到,該系數(shù)在大部分區(qū)域均大于1,表明流場中可壓縮效應顯著增大了進度變量的化學反應源項。圖11給出了使用式(13)可壓縮標度的火焰面/進度變量模型與不使用該修正所得到的對稱軸上溫度分布對比,與圖11中的結果相對應,式(13)的修正通過增大進度變量的源項從而所得到點火位置也提前,也與實驗數(shù)據(jù)更加吻合,初步展示了該可壓縮標度的有效性。該可壓縮修正在復雜超聲速燃燒流動中模擬的可靠性還待進一步開展研究。

4 結 論

本文主要針對火焰面模型以及火焰面/進度變量模型應用于超聲湍流燃燒流動中所開展的可壓縮修正相關研究進行了綜述,主要研究結論包括:

1)對于火焰面模型,其數(shù)據(jù)庫建立時壓強為常數(shù)以及靜焓與混合分數(shù)線性關系的假設在超聲速流動中均不成立。數(shù)值實驗發(fā)現(xiàn):數(shù)據(jù)庫中主要組分濃度對于壓強的變化并不敏感,因此所獲得的流場平均溫度也與建立數(shù)據(jù)庫所給定的壓強值不敏感;雖然靜焓與混合分數(shù)的關系在超聲速流動嚴重偏離線性關系,但這一點引起的火焰面模型獲得的溫度場差別卻很小。因此,嘗試為流場每點建立考慮可壓縮效應的數(shù)據(jù)庫的交互式火焰面模型,與傳統(tǒng)火焰面模型模擬的結果差別較小。

2)火焰面/進度變量模型數(shù)據(jù)庫中進度變量的化學反應源項對于壓強和靜焓-混合分數(shù)關系的變化很敏感,現(xiàn)有考慮這一影響的可壓縮修正主要是通過針對利用未經(jīng)修正的數(shù)據(jù)庫PDF積分后的平均源項進行可壓縮標度來實現(xiàn)。數(shù)值實驗研究發(fā)現(xiàn),可壓縮標度方法可有效考慮高馬赫數(shù)效應對進度變量源項的影響,從而改善火焰面/進度變量模型對超聲速燃燒流動的模擬精度。