汪洋 董剛
(南京理工大學(xué)瞬態(tài)物理重點(diǎn)實(shí)驗(yàn)室,南京 210094)
Richtmyer-Meshkov(RM)不穩(wěn)定是指不同密度且互不浸潤(rùn)的介質(zhì)間所形成的有一定曲率的界面,在脈沖加速(impulsive acceleration)作用下,其兩側(cè)介質(zhì)發(fā)生混合并導(dǎo)致其混合區(qū)(mixing zone,MZ)增長(zhǎng)的行為,其中脈沖加速通??捎杉げ_擊來(lái)實(shí)現(xiàn).RM不穩(wěn)定的機(jī)制來(lái)自斜壓效應(yīng)(baroclinic effect),即,當(dāng)激波與密度界面作用時(shí),界面形成的密度梯度與激波形成的壓力梯度的方向不一致會(huì)導(dǎo)致界面上的渦量沉積,在渦量的誘導(dǎo)下界面上的混合區(qū)逐漸增長(zhǎng),后期這種增長(zhǎng)還會(huì)逐漸發(fā)展為湍流混合.當(dāng)界面本身具有化學(xué)反應(yīng)活性時(shí)(如預(yù)混火焰界面),化學(xué)反應(yīng)會(huì)進(jìn)一步增加RM 不穩(wěn)定的復(fù)雜性.作為一種重要而又復(fù)雜的力學(xué)現(xiàn)象,帶化學(xué)反應(yīng)的RM 不穩(wěn)定頻繁地出現(xiàn)在自然現(xiàn)象(如宇宙超新星爆炸[1])和人類(lèi)生產(chǎn)生活(如航空發(fā)動(dòng)機(jī)超燃混合[2]、慣性約束聚變[3]、工業(yè)爆炸災(zāi)害[4]等)中.因此,開(kāi)展帶化學(xué)反應(yīng)的RM 不穩(wěn)定過(guò)程的研究,對(duì)認(rèn)識(shí)、控制和利用這一現(xiàn)象具有非常重要的現(xiàn)實(shí)意義.
對(duì)惰性界面條件下的RM 不穩(wěn)定過(guò)程,國(guó)內(nèi)外已經(jīng)開(kāi)展了大量的理論、實(shí)驗(yàn)和數(shù)值模擬研究.界面混合區(qū)增長(zhǎng)的理論預(yù)測(cè)模型包含了早期的線性模型[5]以及隨后針對(duì)激波單次作用下界面RM 不穩(wěn)定非線性發(fā)展過(guò)程的漸進(jìn)模型[6]和勢(shì)流模型[7-8]等,基于實(shí)驗(yàn)結(jié)果的適于激波多次作用下的界面混合區(qū)增長(zhǎng)模型也得到了進(jìn)一步研究[9-11].文獻(xiàn)[12-13]對(duì)這些預(yù)測(cè)模型進(jìn)行了很好的總結(jié).此外,有學(xué)者也采用實(shí)驗(yàn)測(cè)試的方法針對(duì)不同入射激波強(qiáng)度和形狀[14-17]、不同初始界面形態(tài)[18-21]條件下,RM 不穩(wěn)定過(guò)程的發(fā)展進(jìn)行了細(xì)致的研究.數(shù)值模擬方面,蔣華等[22-23]采用高精度數(shù)值模擬的方法研究了界面初始形態(tài)對(duì)混合區(qū)早期發(fā)展過(guò)程的影響,并采用理論模型對(duì)混合區(qū)的增長(zhǎng)進(jìn)行了預(yù)測(cè)[23];Schilling 等[24]采用類(lèi)似的數(shù)值格式模擬了激波及其反射激波作用下的單模形態(tài)界面RM 不穩(wěn)定過(guò)程; Thornber 等[25]采用三維數(shù)值模擬,考察了不同初始波長(zhǎng)下的界面混合區(qū)的增長(zhǎng)行為; Lombardini 等[26]利用三維大渦模擬方法研究了單次激波作用下,激波強(qiáng)度對(duì)界面混合區(qū)發(fā)展后期的影響,重點(diǎn)考察了后期界面混合發(fā)展特性;Tritschler 等[27]采用兩種不同的大渦模擬方法研究了激波及其反射波作用下界面混合區(qū)的發(fā)展,對(duì)界面發(fā)展后期小尺度湍流混合特性進(jìn)行了詳盡的分析.
盡管激波與惰性界面相互作用的研究得到了廣泛深入地開(kāi)展,但反應(yīng)性界面在激波誘導(dǎo)下的發(fā)展過(guò)程的研究卻遠(yuǎn)遠(yuǎn)不夠.由于特定形態(tài)的反應(yīng)性界面在實(shí)驗(yàn)控制上比較困難,因此實(shí)驗(yàn)研究幾乎未見(jiàn)報(bào)道,相關(guān)文獻(xiàn)僅限于少量的數(shù)值模擬研究.Khokhlov等[28]采用二維帶化學(xué)反應(yīng)的Navier-Stokes (NS) 方程,數(shù)值模擬了單次激波與單模態(tài)預(yù)混火焰界面的相互作用過(guò)程,指出RM 不穩(wěn)定可促進(jìn)燃燒放熱率,并增加火焰的面積;Kilchyk 等[29]采用帶單步化學(xué)反應(yīng)的二維NS 方程,研究了單次激波或膨脹波與單模預(yù)混火焰界面相互作用的過(guò)程,指出激波壓縮效應(yīng)和火焰長(zhǎng)度對(duì)反應(yīng)性RM 不穩(wěn)定有顯著影響.Massa和Jha[30]采用推導(dǎo)得到的二維線化擾動(dòng)方程,研究了單次激波作用下預(yù)混火焰界面的RM 不穩(wěn)定過(guò)程,指出化學(xué)反應(yīng)誘導(dǎo)的壓力梯度所形成的斜壓效應(yīng)可以減少火焰表面的小尺度結(jié)構(gòu)、減弱火焰界面的增長(zhǎng),因而很可能會(huì)弱化激波驅(qū)動(dòng)的小尺度混合過(guò)程.Dong 等采用帶單步化學(xué)反應(yīng)的NS 方程對(duì)激波及其反射激波與預(yù)混火焰界面的多次作用的二維[31-32]和三維[33]過(guò)程開(kāi)展了系列研究,研究了大尺度流動(dòng)、化學(xué)反應(yīng)以及小尺度擴(kuò)散等幾個(gè)因素對(duì)反應(yīng)性RM 不穩(wěn)定過(guò)程的影響,發(fā)現(xiàn)化學(xué)反應(yīng)在一定程度上可以削弱大尺度流動(dòng)帶來(lái)的界面混合效應(yīng),從而證實(shí)了文獻(xiàn)[30]的推測(cè).此外,Attal 和Ramaprabhu[34]還報(bào)道了非預(yù)混單模火焰界面的RM 不穩(wěn)定過(guò)程,考察了反射激波對(duì)界面的影響,其結(jié)果亦表明激波的多次作用明顯增加了燃燒效率,該工作同時(shí)也反映了反應(yīng)性RM不穩(wěn)定研究開(kāi)始向多樣性發(fā)展.
上述研究表明,化學(xué)反應(yīng)的存在進(jìn)一步增加了RM 不穩(wěn)定現(xiàn)象的復(fù)雜性.雖然目前已有眾多的理論模型可以預(yù)測(cè)激波誘導(dǎo)下惰性界面的混合區(qū)增長(zhǎng)過(guò)程,但對(duì)于反應(yīng)性界面的混合區(qū)在激波誘導(dǎo)下的增長(zhǎng)規(guī)律目前并不明確,例如,文獻(xiàn)[31-33]主要關(guān)注的是RM 不穩(wěn)定過(guò)程中流動(dòng)和化學(xué)反應(yīng)的時(shí)間尺度的變化規(guī)律及其對(duì)RM 不穩(wěn)定的影響,而化學(xué)反應(yīng)作用下混合區(qū)增長(zhǎng)的定量預(yù)測(cè)并未得到研究,也未見(jiàn)其它文獻(xiàn)報(bào)道.為此,本文采用高精度數(shù)值模擬的方法,研究平面入射激波及其反射激波與正弦形狀的預(yù)混火焰界面的相互作用過(guò)程,重點(diǎn)考察不同化學(xué)反應(yīng)活性條件下界面混合區(qū)的發(fā)展過(guò)程,對(duì)化學(xué)反應(yīng)條件下的界面混合區(qū)增長(zhǎng)速率進(jìn)行預(yù)測(cè),以期為反應(yīng)性RM 不穩(wěn)定的混合區(qū)發(fā)展模型的建立提供參考.
控制方程采用二維帶單步化學(xué)反應(yīng)的Navier-Stokes(NS)方程,表達(dá)如下
其中
以上各式中,t代表時(shí)間,x和y代表空間坐標(biāo)長(zhǎng)度,ρ代表體系的密度,u和v分別為x和y方向的速度分量,Y代表反應(yīng)物的質(zhì)量分?jǐn)?shù),p為體系壓力,e是單位質(zhì)量的總能,表達(dá)為
式(8)中,γ 為絕熱指數(shù),Q代表單位質(zhì)量的化學(xué)反應(yīng)放熱.qx和qy分別代表x和y方向的熱流量,表達(dá)為qx=?k?T/?x和qy=?k?T/?x,T是溫度.τij(i=x,y;j=x,y)為黏性應(yīng)力張量分量,表達(dá)為
其中,δi j為Kronecker 函數(shù).以上各式中出現(xiàn)的輸運(yùn)系數(shù),如導(dǎo)熱系數(shù)k,擴(kuò)散系數(shù)D和運(yùn)動(dòng)學(xué)黏度υ 可統(tǒng)一表達(dá)為[4]
本文計(jì)算使用的預(yù)混氣體參照了文獻(xiàn)[35] 中采用的C2H4+3O2+4N2氣體,同時(shí)假定該預(yù)混氣的Lewis數(shù)、Prandtl 數(shù)和Schemidt 數(shù)均為1,則式(10)中的k0,D0和υ0為相同的常數(shù),根據(jù)文獻(xiàn)[4],這些輸運(yùn)系數(shù)取值為k0=D0=υ0=3.2×10?7kg/(s·m·K0.7),n取值為0.7.這些取值能夠保證輸運(yùn)系數(shù)隨溫度的變化關(guān)系很好地符合實(shí)際C2H4+3O2+4N2氣體的這種變化.
本文采用了基于Arrhenius 形式的單步化學(xué)反應(yīng)描述燃燒現(xiàn)象,因此式(7) 中的反應(yīng)速率˙ω 可以表達(dá)為
式中β 代表指前因子,Ea為活化能,R代表通用氣體常數(shù).
對(duì)控制方程(1)采用了分裂算法進(jìn)行求解,其中對(duì)流通量的空間導(dǎo)數(shù)項(xiàng)?F/?x和?G/?y采用的是LF(Lax-Friedrichs)格式結(jié)合9 階WENO(weighted essentially non-oscillatory)模板的方法求解[36];輸運(yùn)通量的空間導(dǎo)數(shù)項(xiàng)?FD/?x和?GD/?y采用了10 階中心差分的格式進(jìn)行求解; 而化學(xué)反應(yīng)源項(xiàng)S和解U的時(shí)間導(dǎo)數(shù)項(xiàng)?U/?t則采用了3 階Runge-Kutta 方法進(jìn)行求解.上述控制方程及其計(jì)算方法在以前的研究中已經(jīng)得到了廣泛應(yīng)用[31-33].
本文模擬的激波與預(yù)混火焰界面相互作用的計(jì)算構(gòu)型如圖1(a) 所示,計(jì)算域?yàn)槎S矩形形狀.初始時(shí)刻,平面入射激波(incident shock wave,ISW) 位于計(jì)算域距下邊界為l1的位置并開(kāi)始向上運(yùn)動(dòng),在運(yùn)動(dòng)過(guò)程中與波長(zhǎng)為λ,初始振幅為a0的正弦形預(yù)混火焰界面(flame interface,FI)發(fā)生相互作用.令火焰界面的下方為密度較大的未燃預(yù)混氣(綠色部分),火焰上方為密度較小的等壓燃燒條件下的已燃?xì)怏w(藍(lán)色部分),則激波從大密度側(cè)向小密度側(cè)越過(guò)火焰界面后會(huì)形成繼續(xù)上行的激波和向下反射的稀疏波.上行激波行至計(jì)算域頂端的剛性壁面后會(huì)反射下行,下行激波則再次和火焰界面發(fā)生相互作用,作用后的激波會(huì)繼續(xù)反射向頂端壁面方向運(yùn)動(dòng),經(jīng)過(guò)再次反射后又會(huì)下行與界面發(fā)生作用.這個(gè)過(guò)程會(huì)一直進(jìn)行下去,直至反射激波強(qiáng)度逐漸減小至激波不再明顯為止.
圖1 正弦預(yù)混火焰界面的(a)初始形態(tài)和(b)典型的發(fā)展過(guò)程.ISW 為入射激波,FI 為火焰界面Fig.1 (a)Initially sinusoidal pattern and(b)developed pattern of flame interface induced by shock wave.ISW is incident shock wave,FI is flame interface
令初始的正弦形火焰界面數(shù)學(xué)形式如下
這里a(x)代表正弦火焰界面的流向(y方向)振幅沿x方向的變化.從界面的角度看,在入射激波的作用下,火焰界面經(jīng)RM 不穩(wěn)定過(guò)程發(fā)生變化,典型的形態(tài)如圖1(b)所示.可以看到,二維正弦形界面在激波的單次作用下演變?yōu)椤搬?spike)?泡(bubble)”結(jié)構(gòu),其中釘結(jié)構(gòu)比較復(fù)雜,它還包含了其上的“帽(cap)”結(jié)構(gòu)和兩側(cè)的“渦(vortex)” 結(jié)構(gòu).定義“帽” 的頂部到“泡” 的底部的流向長(zhǎng)度為界面混合區(qū)寬度WMZ,則該量是評(píng)價(jià)界面RM 不穩(wěn)定發(fā)展過(guò)程的重要物理量[12,24].
本文計(jì)算中,選取計(jì)算域流向長(zhǎng)度為l1+l2=0.6 m,其中l(wèi)1=0.1 m 以保證整個(gè)計(jì)算過(guò)程中入射激波與界面作用后形成的下行稀疏波不會(huì)走出下邊界而影響邊界條件; 計(jì)算域的橫向長(zhǎng)度為正弦界面波長(zhǎng)λ=0.012 8 m.計(jì)算域的左右邊界為周期性邊界條件,上邊界為絕熱剛性無(wú)滑移的壁面條件,下邊界為零梯度邊界條件.假設(shè)本文氣體符合量熱完全氣體的假設(shè),并令初始時(shí)刻波前未燃?xì)怏w(圖1(a)中綠色區(qū)域)的初始溫度和初始?jí)毫Ψ謩e為T(mén)0和p0,則波后未燃?xì)?圖1(a)中紅色區(qū)域)的狀態(tài)可由激波關(guān)系式求出.本文界面初始形態(tài)參數(shù)、預(yù)混氣的熱力學(xué)和動(dòng)力學(xué)參數(shù)的取值[33]均列入表1.
表1 計(jì)算參數(shù)[33]Table 1 Computational parameters[33]
此外,在本文計(jì)算中設(shè)初始平面入射激波Mach數(shù)為M0,則根據(jù)一維激波與界面相互作用的數(shù)值模擬,可以確定入射激波作用后的界面Atwood 數(shù)A=(ρ1?ρ2)/(ρ1+ρ2) (這里ρ1和ρ2分別代表界面兩側(cè)重和輕介質(zhì)的密度)和界面兩側(cè)流向速度差?V,以及第一次反射激波作用后的界面Atwood 數(shù)A1和界面兩側(cè)的流向速度差?V1,這些與一維流動(dòng)有關(guān)的參數(shù)值也列入表1.
為研究不同預(yù)混氣的化學(xué)活性對(duì)RM 不穩(wěn)定發(fā)展的影響,用不同活化能的數(shù)值來(lái)表征不同的化學(xué)反應(yīng)活性,活化能的變化的范圍列入表1 中,此外,還計(jì)算了化學(xué)凍結(jié)流條件下的預(yù)混火焰界面的發(fā)展過(guò)程.化學(xué)凍結(jié)流是一種保留了預(yù)混火焰熱力學(xué)狀態(tài)但將化學(xué)動(dòng)力學(xué)過(guò)程進(jìn)行凍結(jié)(設(shè)Ea無(wú)窮大,即反應(yīng)速率為零)的理想狀態(tài),考慮的目的是為了通過(guò)與反應(yīng)流計(jì)算結(jié)果的對(duì)比來(lái)分析活化能的影響.在本文計(jì)算中,計(jì)算了7 種不同活化能(Ea/RT0=30.0,31.2,31.6,32.0,32.5,38.2,50.0)的反應(yīng)流以及凍結(jié)流(Ea/RT0=∞)的界面混合區(qū)發(fā)展過(guò)程.但在后面的分析中,主要是以Ea/RT0=30.0,38.2,50.0 和凍結(jié)流等4 種典型條件的結(jié)果為主,其中活化能Ea/RT0=38.2的算例為基準(zhǔn)算例,其活化能取值來(lái)自文獻(xiàn)[37].
圖2 不同網(wǎng)格尺寸下計(jì)算的混合區(qū)寬度隨時(shí)間的變化Fig.2 Time histories of width of mixing zone with different grid resolutions
為考察本文計(jì)算結(jié)果的有效性,首先針對(duì)Ea=38.2RT0條件下的激波與預(yù)混火焰相互作用過(guò)程的數(shù)值模擬進(jìn)行了網(wǎng)格無(wú)關(guān)性驗(yàn)證.本文采用了均勻正交的正方形網(wǎng)格來(lái)劃分計(jì)算區(qū)域,選取了三種不同的網(wǎng)格尺寸(?x=?y=0.1 mm,0.05 mm 和0.025 mm)進(jìn)行激波與正弦形預(yù)混火焰界面相互作用過(guò)程的數(shù)值模擬.圖2 給出了3 種網(wǎng)格尺寸下計(jì)算得到的火焰界面發(fā)展過(guò)程中混合區(qū)寬度(WMZ,定義見(jiàn)圖1(b))隨時(shí)間的變化過(guò)程.圖中時(shí)間范圍I,II 分別為入射激波和第一次反射激波作用后WMZ的變化階段,可以發(fā)現(xiàn)階段II 的初期存在一個(gè)短暫的WMZ下降的過(guò)程,這是由于反射激波掠過(guò)火焰界面時(shí)對(duì)界面的壓縮所導(dǎo)致的.圖中結(jié)果表明,在階段I 和II 中,3 種網(wǎng)格尺寸計(jì)算的結(jié)果均相互吻合較好,但在階段II 之后的后繼反射激波與火焰界面的作用過(guò)程中,?x=?y=0.1 mm 尺寸下計(jì)算的結(jié)果要明顯低于?x=?y=0.05 mm 和0.025 mm 兩種尺寸下計(jì)算的結(jié)果,盡管后兩者的計(jì)算結(jié)果之間也略有差別.階段II 后計(jì)算結(jié)果的差異不僅體現(xiàn)了計(jì)算網(wǎng)格尺寸的影響,實(shí)際上也表明火焰界面的失穩(wěn)向著小尺度結(jié)構(gòu)甚至是湍流混合的方向發(fā)展,這超出了本文采用的控制方程和網(wǎng)格尺寸所研究的范圍,因此,本文研究?jī)H考慮階段I 和II 的相互作用過(guò)程.根據(jù)網(wǎng)格無(wú)關(guān)性測(cè)試的結(jié)果,本文采用?x=?y=0.05 mm 的網(wǎng)格尺寸對(duì)所有算例開(kāi)展數(shù)值模擬,每一個(gè)算例總網(wǎng)格數(shù)均為307.2 萬(wàn)個(gè),為提高計(jì)算效率采用了基于MPI的并行策略編制程序開(kāi)展計(jì)算.
圖3 不同化學(xué)活性下混合區(qū)寬度隨時(shí)間的變化. SI 為入射激波作用階段混合區(qū)增長(zhǎng)速率,SII 為反射激波作用階段混合區(qū)增長(zhǎng)速率Fig.3 Time histories of width of mixing zone with different chemical activities. SI is growth rate of mixing zone at incident shock wave stage,SII is growth rate of mixing zone at reflected shock wave stage
圖3 首先給出了3 個(gè)典型活化能(Ea=30.0RT0,38.2RT0和50.0RT0) 的反應(yīng)流以及凍結(jié)流下反應(yīng)區(qū)寬度(WMZ) 隨時(shí)間的變化.可以看到,階段I 中四個(gè)算例在經(jīng)歷了短暫的相同的發(fā)展過(guò)程(t< 0.1 ms)之后,其WMZ的增長(zhǎng)快慢開(kāi)始不同,活化能最小(Ea=30.0RT0)的算例增長(zhǎng)最快,而活化能最大(Ea=50.0RT0) 以及凍結(jié)流的算例增長(zhǎng)最慢.進(jìn)一步觀察可以發(fā)現(xiàn),階段I 的發(fā)展過(guò)程中盡管不同化學(xué)反應(yīng)活性下的界面混合區(qū)增長(zhǎng)快慢不同,但均呈現(xiàn)線性增長(zhǎng)的形式,通過(guò)對(duì)計(jì)算結(jié)果的線性擬合可以確定其WMZ的增長(zhǎng)速率SI.
在進(jìn)一步定量研究SI之前,首先考察了火焰界面在階段I 的演化規(guī)律和混合區(qū)寬度的增長(zhǎng)行為.圖4 顯示了階段I 中4 個(gè)典型時(shí)刻(t=0.2,0.3,0.4 和0.5 ms)下,用反應(yīng)物質(zhì)量分?jǐn)?shù)(Y)標(biāo)記的界面形態(tài)的變化,每個(gè)時(shí)刻下的4 個(gè)子圖從左到右分別代表了活化能為Ea=30.0RT0,38.2RT0和50.0RT0的反應(yīng)流以及凍結(jié)流等4 種不同反應(yīng)活性的情況.在t=0.2 ms時(shí)(圖4(a)),4 種條件下的界面結(jié)構(gòu)均呈現(xiàn)出“釘?泡” 結(jié)構(gòu),釘結(jié)構(gòu)的上方和兩側(cè)也分別呈現(xiàn)出“帽”結(jié)構(gòu)和渦結(jié)構(gòu)(見(jiàn)圖1(b)的說(shuō)明).但是,當(dāng)反應(yīng)活性強(qiáng)(Ea=30.0RT0)時(shí),渦結(jié)構(gòu)內(nèi)反應(yīng)物幾乎完全消耗,這說(shuō)明化學(xué)反應(yīng)的作用首先是發(fā)生在渦結(jié)構(gòu)中的,其原因主要是因?yàn)镽M 不穩(wěn)定過(guò)程中形成的渦結(jié)構(gòu)使得界面兩側(cè)的未燃?xì)馀c已燃?xì)膺M(jìn)行了充分的摻混,這種強(qiáng)化的傳熱和傳質(zhì)過(guò)程進(jìn)一步促進(jìn)了化學(xué)反應(yīng)的進(jìn)行.該時(shí)刻下的結(jié)果還表明,不同反應(yīng)活性下混合區(qū)寬度的變化差異不大.當(dāng)t=0.3 ms 時(shí)(圖4(b)),可以發(fā)現(xiàn)隨著反應(yīng)活性的增強(qiáng),渦結(jié)構(gòu)內(nèi)的混合過(guò)程逐漸被化學(xué)反應(yīng)過(guò)程所替代,其中,Ea=30.0RT0條件下的反應(yīng)流其渦結(jié)構(gòu)不再明顯,界面上方形態(tài)主要以“釘?帽” 結(jié)構(gòu)為主并明顯向上移動(dòng),從而導(dǎo)致其混合區(qū)寬度也較其他幾種情況更長(zhǎng).當(dāng)t=0.4 ms時(shí)(圖4(c)),除了界面上方的“釘?帽?渦” 結(jié)構(gòu)繼續(xù)發(fā)展外,還可以看到界面下方的泡結(jié)構(gòu)隨著反應(yīng)活性的增強(qiáng)也逐漸向下移動(dòng),活性越強(qiáng),“泡”結(jié)構(gòu)向下發(fā)展越快,因此導(dǎo)致整個(gè)混合區(qū)寬度隨著活性的增強(qiáng)而逐漸變大,進(jìn)而導(dǎo)致圖3 中階段I 的WMZ變化規(guī)律.到t=0.5 ms 時(shí)(圖4(d)),可以看到,反應(yīng)活性最強(qiáng)(Ea=30.0RT0) 的算例,其“帽” 結(jié)構(gòu)已經(jīng)由化學(xué)反應(yīng)所消耗,因而只留下釘結(jié)構(gòu),這也導(dǎo)致其反應(yīng)區(qū)寬度出現(xiàn)了突然的下降(見(jiàn)圖3);而活性次強(qiáng)(Ea=38.2RT0)的算例,其渦結(jié)構(gòu)不再明顯,但形成的“釘?帽”結(jié)構(gòu)開(kāi)始使界面混合區(qū)明顯向上發(fā)展,這與Ea=30.0RT0的算例在t=0.3 ms 時(shí)的結(jié)構(gòu)類(lèi)似.此外,隨著反應(yīng)活性的增強(qiáng),界面下方的泡結(jié)構(gòu)也進(jìn)一步向下發(fā)展.
圖4 階段I 中不同時(shí)刻的界面形態(tài)Fig.4 Patterns of flame interface at the different instants for stage I
根據(jù)以上的界面形態(tài)的演化過(guò)程,可以發(fā)現(xiàn),階段I 中混合區(qū)寬度隨時(shí)間的變化是由“釘?帽”結(jié)構(gòu)向上發(fā)展以及泡結(jié)構(gòu)向下發(fā)展所共同決定的.化學(xué)反應(yīng)活性對(duì)這一過(guò)程的影響十分顯著,化學(xué)反應(yīng)越強(qiáng)(活化能越小),“釘?帽”向上發(fā)展以及“泡”向下發(fā)展越顯著,因此使得反應(yīng)區(qū)寬度越大.
2.1 節(jié)首先觀察了不同反應(yīng)活性條件下界面形態(tài)及界面混合區(qū)寬度隨時(shí)間的演化規(guī)律,發(fā)現(xiàn)化學(xué)反應(yīng)對(duì)界面上方的“釘?帽”結(jié)構(gòu)以及下方的“泡”結(jié)構(gòu)均有顯著影響.為發(fā)展定量預(yù)測(cè)階段I 的混合區(qū)增長(zhǎng)速率SI的模型,圖5 進(jìn)一步考察了t=0.4 ms 時(shí)的界面結(jié)構(gòu)及流場(chǎng)特征.圖5(a) 將四種反應(yīng)活性條件下界面輪廓(用Y=0.9 的等值線表示)進(jìn)行了疊加對(duì)比,可以發(fā)現(xiàn),界面下方“泡”結(jié)構(gòu)下緣隨著反應(yīng)活性的增加不斷下移,這種下移符合預(yù)混火焰?zhèn)鞑サ奶卣?反應(yīng)活性越大,化學(xué)反應(yīng)速率越大,則火焰?zhèn)鞑ニ俣纫苍酱?因此可以得出結(jié)論:反應(yīng)性預(yù)混火焰界面混合區(qū)寬度的增長(zhǎng)其部分原因是由預(yù)混火焰?zhèn)鞑ニ暙I(xiàn)的.根據(jù)經(jīng)典的Frank-Kamenetsky 層流火焰?zhèn)鞑ダ碚揫38],層流預(yù)混火焰的傳播速度可以表達(dá)為
式中,Sb為“泡” 結(jié)構(gòu)的向下發(fā)展速度,也代表了層流預(yù)混火焰的傳播速度; α 為熱輸運(yùn)系數(shù),α=k0Tn,n=0.7;k0,β,Ea和R的含義及取值見(jiàn)表1;T代表入射波后的已燃區(qū)溫度,本文取T=7.0T0.
此外,觀察界面上方“釘?帽” 結(jié)構(gòu),可以看出,反應(yīng)活性越強(qiáng),則“帽” 結(jié)構(gòu)的上緣越向上移,但其結(jié)構(gòu)較為復(fù)雜,并不符合預(yù)混火焰的傳播規(guī)律.為此,圖5(b) 進(jìn)一步給出了界面附近的渦量幅值|ω|=|?0.5(?u/?y??v′/?x)|以及速度矢量V=(u,v′)的分布,這里v′=v?vp,vp為入射激波波后流場(chǎng)的平均速度.首先觀察渦量幅值分布可以看出,隨著反應(yīng)活性的增強(qiáng),流場(chǎng)的渦量明顯減弱,這說(shuō)明化學(xué)反應(yīng)的加入可以抑制渦量的形成,這與之前研究得到的結(jié)論相一致[30].接下來(lái)觀察速度矢量的分布可以看到,當(dāng)化學(xué)反應(yīng)較強(qiáng)(如Ea=30.0RT0) 時(shí),釘結(jié)構(gòu)兩側(cè)的渦量變?nèi)?漩渦結(jié)構(gòu)不明顯,因此對(duì)入射激波波后流動(dòng)速度抑制不夠明顯,因此“釘?帽” 結(jié)構(gòu)上方的流體流向速度較大.反之,當(dāng)化學(xué)反應(yīng)較弱(如Ea=50.0RT0) 或化學(xué)反應(yīng)凍結(jié)時(shí),由于釘兩側(cè)存在較強(qiáng)的互為反向的渦旋,渦旋外側(cè)誘導(dǎo)的向下的運(yùn)動(dòng)速度實(shí)際上抵消了波后流體向上的運(yùn)動(dòng).
圖5 t=0.4 ms 時(shí)不同化學(xué)活性條件下(a)反應(yīng)物質(zhì)量分?jǐn)?shù)Y=0.9 的等值線對(duì)比,(b)流場(chǎng)渦量幅度和速度矢量分布對(duì)比Fig.5 (a)Comparisons of(a)mass fraction contours with Y=0.9 and of(b)vorticity amplitude and velocity vector among cases with different chemical activities at t=0.4 ms
為驗(yàn)證上述分析,圖6 給出了計(jì)算域中心(x=λ/2) 處由于渦的旋轉(zhuǎn)所誘導(dǎo)的界面兩側(cè)的速度差?Va(本文稱(chēng)誘導(dǎo)速度)沿流向的分布,該速度差由界面上、下側(cè)中心處流向速度的差值所確定.圖6 的結(jié)果表明,雖然在較弱活性(Ea=50.0RT0)或凍結(jié)流的條件下,較強(qiáng)的漩渦對(duì)在其之間(y=0.031 6 m 左右)的流體所誘導(dǎo)的流向速度差較大,但在“釘?帽” 結(jié)構(gòu)的上方(y>0.032 5 m)其速度差迅速變小至反向流動(dòng),說(shuō)明漩渦外側(cè)的向下流動(dòng)對(duì)速度差的減小起著很大作用.相反,對(duì)活性較強(qiáng)的反應(yīng)流(Ea=30.0RT0),由于漩渦強(qiáng)度弱,漩渦外側(cè)誘導(dǎo)的反向流動(dòng)速度較小,因此渦對(duì)之間形成的正向速度差可以一直維持.通過(guò)以上分析可知,界面“釘?帽” 結(jié)構(gòu)受化學(xué)反應(yīng)的影響比較復(fù)雜,但是較強(qiáng)的反應(yīng)活性導(dǎo)致了其上方的具有較大的誘導(dǎo)速度差,進(jìn)而使得“釘?帽” 結(jié)構(gòu)向上發(fā)展得更快,由此可將其向上的發(fā)展速度表達(dá)為
圖6 t=0.4 ms 時(shí)不同化學(xué)活性條件下軸線上流向誘導(dǎo)速度分布Fig.6 Distributions of streamwise induced velocity at centerline for cases with different chemical activities at t=0.4 ms
式中?V代表一維構(gòu)型下入射激波波后的界面兩側(cè)流向速度差,其值參見(jiàn)表1,δ1為階段I 的經(jīng)驗(yàn)常數(shù).式(14)等號(hào)右端exp(?Ea/RT)反映了化學(xué)反應(yīng)活性對(duì)界面上緣“釘?帽” 結(jié)構(gòu)發(fā)展的影響,反應(yīng)活性越強(qiáng),該項(xiàng)越大,則Ss越大,因此符合圖5(b)和圖6 的分析.
本文假設(shè)所考慮的流場(chǎng)具有層流流動(dòng)特性,因此化學(xué)反應(yīng)和湍流的相互作用可以忽略,故綜合上述分析,可以采用線性疊加的方式給出階段I 中火焰界面的混合區(qū)寬度增長(zhǎng)速率SI的表達(dá)式
式中,SRM代表RM 不穩(wěn)定導(dǎo)致的混合區(qū)寬度增長(zhǎng)速率,Sc代表化學(xué)反應(yīng)導(dǎo)致的混合區(qū)寬度增長(zhǎng)速率.前者可以采用傳統(tǒng)的各類(lèi)基于RM 不穩(wěn)定的理論模型進(jìn)行預(yù)測(cè)[12-13],后者由于層流流動(dòng)的假設(shè)仍可采用如下線性疊加的方法進(jìn)行求解
式中Sb和Ss分別由式(13)和式(14)進(jìn)行求解,其中涉及到的變量?jī)H為活化能Ea,表明了化學(xué)反應(yīng)活性的影響,而其余各常量均可由初始條件和一維激波界面作用得到的參數(shù)來(lái)確定.
注意到式(14)中δ1為經(jīng)驗(yàn)參數(shù),為確定該值,本文對(duì)圖3 數(shù)值模擬得到的結(jié)果進(jìn)行線性擬合可以得到SI,令凍結(jié)流的SI=SRM(即凍結(jié)流的WMZ增長(zhǎng)速率不考慮化學(xué)反應(yīng),僅由RM 不穩(wěn)定所導(dǎo)致),則由式(15)可以求得數(shù)值模擬條件下的不同Ea的Sc值.圖7 給出了采用數(shù)值模擬得到的Sc隨Ea的變化關(guān)系,可以看出,隨著活化能的減小,Sc呈指數(shù)增大.進(jìn)一步采用式(13)、式(14)和式(16)進(jìn)行計(jì)算,發(fā)現(xiàn)當(dāng)δ1=9.3 時(shí),利用式(16)預(yù)測(cè)模型得到的Sc與數(shù)值模擬得到的Sc吻合良好.圖7 中還給出預(yù)測(cè)模型中Sb和Ss隨活化能的變化,可以發(fā)現(xiàn),隨著活化能的減小,火焰界面上方“釘?帽” 結(jié)構(gòu)的增長(zhǎng)速度比界面下方“泡”的增長(zhǎng)速度要更快.
圖7 數(shù)值模擬與理論預(yù)測(cè)的Sc 的比較Fig.7 Comparisons of Sc between numerical and predicted results
觀察圖3 中階段II(即第一次反射激波作用后的火焰界面發(fā)展階段)可以看出,反射激波作用后,不同化學(xué)反應(yīng)活性的反應(yīng)流和凍結(jié)流條件下,界面混合區(qū)寬度的增長(zhǎng)速率(用SII表示)趨于一致,其增長(zhǎng)速率也明顯大于階段I 的增長(zhǎng)速率.為了解WMZ的增長(zhǎng)規(guī)律,圖8 顯示了階段II 的兩個(gè)典型時(shí)刻(t=0.7 ms和0.8 ms)下,不同化學(xué)反應(yīng)活性下界面的發(fā)展形態(tài).圖8(a)的結(jié)果(t=0.7 ms)表明,反射激波與火焰界面作用后,混合區(qū)寬度仍然沿流向持續(xù)增長(zhǎng),當(dāng)反應(yīng)活性較強(qiáng)時(shí),界面呈現(xiàn)出沿流向的細(xì)長(zhǎng)結(jié)構(gòu),隨著反應(yīng)活性的下降,界面形狀逐漸變得復(fù)雜,已燃?xì)夂臀慈細(xì)庀嗷交煨纬傻亩S混合區(qū)變得越來(lái)越明顯.然而從混合區(qū)寬度的變化來(lái)看,不同活性下的混合區(qū)寬度變化并不明顯,即,隨著反應(yīng)活性的增強(qiáng),界面上方的“釘?帽”結(jié)構(gòu)和下方的“泡結(jié)構(gòu)”均同時(shí)下移,表明此時(shí)界面上下存在同樣的預(yù)混火焰?zhèn)鞑ミ^(guò)程,兩者傳播方向相同,因此對(duì)界面上、下方發(fā)展的影響相互抵消,從而導(dǎo)致混合區(qū)寬度對(duì)反應(yīng)活性的依賴(lài)不再變得明顯.圖8(b)的結(jié)果(t=0.8 ms)也顯示了類(lèi)似的規(guī)律.需要注意的是凍結(jié)流條件下,界面形態(tài)開(kāi)始失去對(duì)稱(chēng)性,這表明到第一次反射激波后期界面的小尺度流動(dòng)開(kāi)始增加,界面逐漸向湍流混合的階段發(fā)展,因此,本文的研究?jī)H限于第一次反射激波作用后的階段II.后繼的二次及多次激波與火焰界面的相互作用不在本文研究的范疇.
圖8 階段II 中不同時(shí)刻的界面形態(tài)Fig.8 Patterns of flame interface at the different instants for stage II
圖9 進(jìn)一步給出階段II 中t=0.7 ms 時(shí),不同反應(yīng)活性條件下中心線上流向誘導(dǎo)速度?Va沿流向的分布.由結(jié)果也可以看出界面下方反射激波作用之后的區(qū)域(y=0.39 ~0.41 m) 和界面上方的區(qū)域(y> 0.45 m)誘導(dǎo)速度均趨于零,這說(shuō)明入射激波和反射激波對(duì)火焰界面附近的渦旋運(yùn)動(dòng)的影響相互抵消,因而界面混合區(qū)的發(fā)展不再考慮渦旋運(yùn)動(dòng)對(duì)化學(xué)反應(yīng)的影響,這一點(diǎn)與階段I 中界面混合區(qū)的發(fā)展有所不同.
圖9 t=0.7 ms 時(shí)不同化學(xué)活性條件下中心線上流向誘導(dǎo)速度分布Fig.9 Distributions of streamwise induced velocity at centerline for cases with different chemical activities at t=0.7 ms
采用與式(15)相同的形式表達(dá)階段II 的混合區(qū)增長(zhǎng)速率SII
式中變量符號(hào)與式(15)相同,下標(biāo)II 代表階段II.根據(jù)前面分析,由于渦旋運(yùn)動(dòng)不再對(duì)火焰界面產(chǎn)生額外影響,界面下方“泡”結(jié)構(gòu)和上方的“釘?帽”結(jié)構(gòu)僅存在相似的火焰向下傳播的過(guò)程,因此導(dǎo)致兩者作用相互抵消,即Sc,II=0.因此,混合區(qū)寬度的變化僅由惰性RM 不穩(wěn)定過(guò)程所決定,采用反射波后的Mikaelian 模型[8]對(duì)其增長(zhǎng)速度進(jìn)行預(yù)測(cè),式(17)可表達(dá)為
式中δ2為階段II 的經(jīng)驗(yàn)常數(shù),?V1和A1的取值見(jiàn)表1.由7 種不同活化能條件下的反應(yīng)流以及凍結(jié)流的數(shù)值模擬結(jié)果可以擬合得到不同的SII值,由此得到的平均值為149.1 m/s,標(biāo)準(zhǔn)偏差為3.2 m/s,將平均增長(zhǎng)速度(149.1 m/s)代入式(17)可以得出δ2=0.79.需要說(shuō)明的是,原始的Mikaelian 模型[8]中經(jīng)驗(yàn)常數(shù)δ2=0.28,但進(jìn)一步的研究也表明[12],δ2的數(shù)值依賴(lài)于界面的形態(tài)與結(jié)構(gòu),比如對(duì)惰性界面和密度梯度很大的界面,二維和三維計(jì)算采用的δ2分別可取0.55和0.84.本文δ2取0.79 反映了預(yù)混火焰界面與原始Mikaelian 模型中不存在溫差的界面是有所區(qū)別的.
采用帶單步反應(yīng)的Navier-Stokes 方程,對(duì)入射激波及其反射激波與正弦形狀的預(yù)混火焰界面相互作用過(guò)程進(jìn)行了二維數(shù)值模擬,考察了介質(zhì)反應(yīng)活性對(duì)界面形態(tài)演化和混合區(qū)增長(zhǎng)過(guò)程的影響.得到的主要結(jié)論如下:
(1)正弦形預(yù)混火焰界面在入射激波作用后沿流向不斷拉伸,除常規(guī)RM 不穩(wěn)定誘導(dǎo)的界面混合區(qū)增長(zhǎng)外,與化學(xué)反應(yīng)相關(guān)的預(yù)混火焰?zhèn)鞑?duì)界面下方的“泡”結(jié)構(gòu)的發(fā)展起到了促進(jìn)作用,而化學(xué)反應(yīng)與渦結(jié)構(gòu)的相互作用減弱了渦的強(qiáng)度,促進(jìn)了界面上方的流體流向運(yùn)動(dòng)進(jìn)而使“釘?帽”結(jié)構(gòu)沿流向不斷發(fā)展.化學(xué)活性越強(qiáng),界面的流向增長(zhǎng)越快.
(2) 沿流向發(fā)展的火焰界面在反射激波的作用下,沿流向繼續(xù)發(fā)展,其中界面下方的“泡” 結(jié)構(gòu)和上方“釘?帽” 結(jié)構(gòu)均受到火焰向下傳播的影響,故化學(xué)反應(yīng)的作用抵消,從而導(dǎo)致反射激波作用后界面的發(fā)展基本不依賴(lài)于介質(zhì)的化學(xué)反應(yīng)活性.
(3)根據(jù)上述界面演化機(jī)理,一方面提出了入射激波作用階段,化學(xué)反應(yīng)導(dǎo)致的界面增長(zhǎng)速率Sc的表達(dá)式,該表達(dá)式僅依賴(lài)于介質(zhì)的反應(yīng)活化能和初始條件,具有較好的通用性;另一方面指出仍可以采用Mikaelian 的模型預(yù)測(cè)反射激波作用下火焰界面的混合區(qū)增長(zhǎng)速率,但模型常數(shù)要進(jìn)行修正.