斜爆轰的胞格结构及横波传播*
2010-02-26王爱峰姜宗林
王爱峰,赵 伟,姜宗林
(中国科学院力学研究所高温气体动力学重点实验室,北京100190)
斜爆轰是超声速的可燃气体通过一定角度的楔面经过斜激波压缩、燃烧诱导直至快速放热的一种特殊的爆轰现象。由于斜爆轰波可以驻定,燃烧时间极短且拥有很高的燃烧效率,因此在超燃冲压推进系统(supersonic combustion ramjet,scramjet)和冲压加速器(ram accelerator)[1]等方面都有良好的应用前景。
D.T.Pratt 等[2]主要应用传统的爆轰理论、简化的物理模型,研究马赫数、放热量、楔面角对斜爆轰驻定的影响,结果表明:对于固定当量比的可燃气体(即固定放热量),不同的马赫数,形成驻定斜爆轰需要楔面角满足一定的范围,这个范围为CJ 线与极限角连线所夹的区域(如图1 所示)。H.F.Lehr[3]在观测高速弹丸在氢氧介质中飞行时,发现弹头产生斜激波诱导的燃烧现象,观测得到不同的燃料当量比、飞行马赫数以及弹丸结构等条件下出现的不同的斜爆轰波结构。T.Fujiw ara 等[4]模拟了无粘条件下的斜爆轰结构,得到斜爆轰结构中出现斜激波、爆燃波与斜爆轰波共存的复杂结构(如图2 所示),E.K.Dabora 等[5]则通过斜激波管实验观测到不同气体条件下的斜爆轰波的结构。由于楔面诱导的斜爆轰结构中涉及爆燃波与爆轰波之间复杂的相互作用,随着对斜爆轰现象认识的深化,近年来主要的研究集中在斜爆轰内部流场的精细研究。L.F.F.da Silva 等[6]利用氢氧混合基元反应模型,分析不同楔面角、不同马赫数、不同氢氧当量比条件下斜爆轰波的发展演化,并对诱导时间与诱导长度进行研究。C.Viguier 等[7]通过实验,发现斜爆轰中的胞格结构的横波传播特征(见图3),J.Y.Choi 等[8]的数值研究发现了高活化能状态下的特殊的胞格结构,同时发现横波的单向传播现象。
图1 斜爆轰驻定匹配关系(Q~=1.96)Fig.1 Matching relationship of standing ODW
图2 斜激波-斜爆轰波结构示意图Fig.2 Schematic diagram of OSW-ODW structure
为了探索斜爆轰形成机制并分析其中的胞格结构,本文中利用数值方法,分析在30°楔角条件下,临界爆轰马赫数附近6.8、7.0 与7.5 的3 种马赫数状态下的斜爆轰结构,观察斜爆轰波阵面产生的横波的传播特征。
1 控制方程与数值方法
1.1 控制方程
忽略粘性、热传导,可以采用二维欧拉方程耦合化学反应过程描述斜爆轰波。化学反应模型采用了总包化学反应模型,通过反应进度常数λ反应化学反应的完全度。
控制方程为
图3 斜爆轰波结构以及横波传播的纹影图[7]Fig.3 Schlieren photograph of ODW structure and propagation of transverse wave
式中:ρ为混合气体密度,u 和v 分别为x 和y 方向上的速度,E 为质量总能量,λ为反应进程度。p 为气体压力,满足以下关系
式中:γ为绝热指数,q 为化学反应生成热。考虑单步不可逆反应模型,化学反应生成源项表达式为
式中:Ea 为化学反应的活化能,k 为化学反应速率系数。
1.2 数值方法
控制方程的对流项采用三阶精度ENO 格式[9]离散,考虑迎风效应,采用Steger-Warming 流通量分裂[10]。时间采用三阶精度TVD Runge-Kutta 方法进行迭代,化学反应项采用解耦的方式处理。
所选计算区域为图4 中矩形区域,并按照楔面角度进行旋转。计算域大小为10 cm×4 cm,θ为楔面角,β 为斜爆轰角。为了研究斜爆轰波的复杂的细部结构,采用2 000×800 的网格进行计算,同时在入口边界增加3 个虚点,保证上游的入口超声速流体不会在楔面上形成数值反射,楔面采用滑移边界条件,出口边界为无扰动外插边界。
为了研究马赫数对斜爆轰胞格结构的影响,采用如下的物理参数:初始压力p0=101.325 kPa,初始温度T0=300 K,楔面角θ=30°,马赫数Ma=6.8,7.0,7.5,放热量Q=965 kJ/kg,Ea=4.794 MJ/kg,k=7.5 ns-1,γ=1.29。所设定的马赫数、放热量与楔面角关系均满足斜爆轰驻定条件[2]。
图4 斜爆轰结构示意图Fig.4 Schematic diagram of ODW structure
2 数值验证
2.1 斜爆轰角验证
为了验证不同状态下斜爆轰角的变化,假定平面斜爆轰:(1)来流为预混未燃气体,混合均匀,流动稳定;(2)化学反应层等效于放热,放热量为Q;(3)流动无粘绝热。根据上面的假设条件,利用流体守恒方程,可以得到下面关于放热量、楔角、马赫数和斜爆轰角之间的关系
式中:(Ma)1 为斜爆轰波前马赫数,(Ma)1n 为波前法向马赫数,θ为楔面角,β 为斜爆轰角,为量纲一放热量,cp为定压比热,T1为波前温度,γ为绝热指数。
楔面角θ=30°、马赫数不同、放热量不同时,斜爆轰角的理论解和数值结果及相对误差见表1。由表可见,数值解比理论解略大,最大相差低于1.25%,说明数值解的结果是可信的。
2.2 斜爆轰波流场验证
为了进一步验证斜爆轰结构中复杂的结构,将数值结果与实验结果进行对比,图5 是实验纹影图[5]与数值密度云图的对比结果。可以看到典型的斜激波、斜爆轰波以及斜激波、斜爆轰波与爆燃波交汇成的三波点。由于斜爆轰波与爆燃波波后的气体密度不同形成一条向下游发展的接触间断面。因此所采用的数值计算能够对斜爆轰波结构进行详细的解析。
表1 斜爆轰角Table 1 ODW angles
图5 实验纹影图[5]和数值密度云图(p0=101.325 kPa,T 0=293 K)Fig.5 Schlieren image of the flow field and density field
3 数值结果与讨论
3.1 不同马赫数条件下的温度流场
图6 为放热量Q=965 kJ/kg、马赫数为6.8、7.0 和7.5 的温度云图。为了保证流场定常,每种状态下的密度残差都小于10-4。从图中可以看出,马赫数为7.5 时,上游的斜激波区域长度较短,斜激波与斜爆轰波过渡平滑,流场中的胞格结构规则且胞格较小,波阵面呈现一条直线。而马赫数为7.0、6.8时,斜激波区域长度相对较长,波后流场存在较大的扰动,下游的胞格结构不再规则且胞格大小不一,斜爆轰波阵面总体趋向于直线,但局部出现弯曲现象。
马赫数的差别不大,但流场的结构存在较大的波动。进一步观察三波点的结构发现,马赫数为7.0、6.8 时斜激波区较长,且三波点附近的流场结构较为复杂;而马赫数为7.5 时,三波点附近的流场较为简单,因此有必要对主流三波点附近的流场结构进行细致的研究。
图7 为马赫数为7.0 时的三波点附近的流场图,从图中可以看出,三波点后会产生一道激波,这道激波传播到楔面发生反射,并能透过滑移间断面,与爆轰波阵面产生的横波相互作用后下游流场受到扰动,使得胞格结构呈现不规则的状态。同时爆轰波阵面产生的横波传播到楔面后发生了反射,使得接触间断面下游出现一系列的涡结构。
为了进一步分析,截取y=3 mm 的水平方向的压力和化学反应度曲线(见图8)。图中圆点部分为三波点后产生的一道激波,方点部分为激波在楔面上发生反射形成的反射激波。A 区压力经历突变,最终化学反应度为1,热量完全释放。B 区压力下降,是爆燃区域。C 区为激波与反射激波包围的混合区域。C 区之后由于爆轰前导反应段面上产生的横波在楔面上发生反射,形成振荡的D 区。
3.2 波阵面结构特征
前面的研究发现,在马赫数为7.0、6.8 时流场中出现不规则的胞格结构,而马赫数为7.5 时胞格结构却非常规则。三波点后产生的激波能够与波阵面上产生的横波相互作用,并产生扰动影响着下游流场的结构,同时横波传播到楔面也能发生反射,使滑移间断面呈现一系列的涡结构。这些现象的出现都与横波存在很大的关系,因此对波阵面上横波传播的研究是非常必要的。
分别对马赫数为7.5 和6.8 等2 种情况下的爆轰波阵面的结构进行细部分析,如图9 ~10 所示。当马赫数为7.5 时,斜爆轰波阵面上产生的横波均向同一个方向传播(如图11(a)所示),这与正爆轰的横波的双向传播特性不同,但这种单向传播的横波与文献[8]一致。当马赫数为6.8 时,由于上下游波阵面的结构相异,因此分为上游A 区和下游B 区(见图6(a)),研究表明,波后流场的上游横波仍然是单向传播(如图11(a)所示),但是到了下游,随着马赫干的出现,横波的传播不再单向,而是类似于正爆轰中双向传播的横波结构(如图11(b)所示),这种横波单向与双向传播的同时出现在以前没有发现过。
图6 相同放热量、不同马赫数条件下温度云图(Q=965 kJ/kg)Fig.6 Temperature contour at different Mach numbers
图7 三波点附近的流场图(Ma=7.0,Q=965 kJ/kg)Fig.7 Flow field near primary triple point area
图8 沿楔面方向压力与化学反应度分布(y=3 mm,Ma=7.0,Q=965 kJ/kg)Fig.8 Pressure and degree of chemical reaction profile along w edge direction
马赫数越大,加入的动能越大,通过压缩获得的温度越高,所得斜激波部分的长度越短,因此在一定尺度范围内,三波点后产生的激波对波后流场影响也很小,因此马赫数为7.5 时波后流场结构均匀,同时横波只出现单向传播的机制。但是马赫数为7.0、6.8 时,横波是由单向转为双向传播的,两种传播状态并存,如果上述的这种横波的传播方式是斜爆轰中一般的传播规律,那么马赫数为7.5 时必然也会出现。如果这个假定成立,那么影响这个方面的因素必然是长度尺度,因此增加长度尺度来考察马赫数为7.5 时的波阵面形成的横波的传播情况。
长度由10 cm 增加到15 cm、其余参数不变时的流场如图12 ~13 所示。由图13 可以看出,在11 cm后就开始出现马赫干,随后横波的双向传播越来越明显。因此横波由单向向双向传播的转变不是偶然的,它与长度尺度相关,是其中复杂波系逐渐削弱的一种体现。
图9 局部温度云图(Ma=7.5,Q=965 kJ/kg)Fig.9 Detailed temperature contour
图10 局部温度云图(Ma=6.8,Q=965 kJ/kg)Fig.10 Detailed temperature contour
图11 波阵面简图Fig.11 Schematic of w ave front structure
图12 温度云图(Ma=7.5,Q=965 kJ/kg,L=15 cm)Fig.12 Temperature contour
图13 局部温度云图(Ma=7.5,Q=965 kJ/kg,L=15 cm)Fig.13 Detailed temperature contour
4 结 论
研究了30°楔角条件下,临界爆轰马赫数附近的3 种来流状态时斜爆轰的胞格结构以及横波的传播机制,得到如下主要结论:
(1)马赫数为7.0、6.8 时,诱导长度较长,一定尺度范围内,斜爆轰的胞格结构呈现不规则的状态。
(2)主流三波点后形成一道激波,并在楔面上发生反射,穿越滑移间断面与斜爆轰波阵面产生的横波相互作用,使得横波传播发生偏向,对下游流场产生扰动。
(3)斜爆轰波阵面产生的横波在上游呈现单向传播状态,而在下游过渡到双向传播。增加长度尺度,马赫数为7.5 时,在11 cm 处也开始出现双向传播的横波,说明上述横波由单向向双向传播的转变不是偶然的,传播的状态依赖于长度尺度。
[1] Sasoh A.Laser-propelled ram accelerator[J].Journal de Physique Ⅳ,2000,10:11-14.
[2] Pratt D T,H umphrey J W.Morphology of standing oblique detonation w aves[J].Propulsion and Pow er, 1991,7(5):837-845.
[3] Lehr H F.Experiments in shock induced combustion[J].Astronautica Acta,1972,17(4-5):589-597.
[4] Fujiw ara T, Matsuo A.A tw o-dimensional detonation supported by a blunt body or a w edge[R].AIAA Paper 88-0098,1988.
[5] Dabora E K,Desbordes D,Wagner H G.Oblique detonation at hypersonic velocities[C]∥Borisov A A,Kuhl A L,Leyer J C,et al.Dynamics of Detonations and Explosions:Detonations,Progress in Astronautics and Aeronautics.Washington:AIAA,1991:187-201.
[6] da Silva L F F, Deshaies B.Stabilization of an oblique detonation wave by a wedge:A parametric numerical study[J].Combustion and Flame,2000,121(1-2):152-166.
[7] Viguier C,Gourara A,Desbordes D.Onset of oblique detonation waves:Comparison between experimental and numerical results for hydrogen-air mixtures[C]∥Proceedings of Twenty-Seventh Symposium(International)on Combustion.Pittsburgh:The Combustion Institute,1998:3023-3031.
[8] Choi J Y, Kim D W,Jeung I S, et al.Cell-like structure of unstable oblique detonation w ave f rom high-resolution numerical simulation[J].Proceeding of the Combustion Institute,2007,31(2):2473-2480.
[9] Shu C W,Osher S.Efficient implementation of essentially non-oscillatory shock capturing schemes Ⅱ[J].Computational Physics,1989,83(1):32-78.
[10] Steger J L, Warming R F.Flux vector splitting of the inviscid gasdynamic equations with application to finitedifference methods[J].Computational Physics,1981,40(2):263-293.