圆球诱发斜爆轰波的数值研究1)

方宜申 胡宗民 滕宏辉姜宗林

(中国科学院力学研究所高温气体动力学国家重点实验室，北京100190)

(中国科学院大学工程科学学院，北京100049)

圆球诱发斜爆轰波的数值研究1)

方宜申 胡宗民 滕宏辉2)姜宗林

(中国科学院力学研究所高温气体动力学国家重点实验室，北京100190)

(中国科学院大学工程科学学院，北京100049)

斜爆轰发动机是飞行器在高马赫数飞行条件下的一种新型发动机，具有结构简单、成本低和比冲高等优点.但是斜爆轰发动机的来流马赫数范围广，来流条件复杂，为实现斜爆轰波的迅速、可靠引发，采用钝头体来诱发.利用Euler方程和氢氧基元反应模型，对超声速氢气/空气混合气体中圆球诱导的斜爆轰流场进行了数值研究.不同于楔面诱发的斜爆轰波，球体首先会在驻点附近诱发正激波/爆轰波，然后在稀疏波作用下发展为斜激波/爆轰波.模拟结果显示，经过钝头体压缩的预混气体达到自燃温度后，会出现两种流场：当马赫数较低时，由于稀疏波的影响，燃烧熄灭，钝头体下游不会出现燃烧情况；而当马赫数较高时，燃烧阵面能传到下游.分析表明，当钝头体的尺度较小时，驻点附近的能量不足以诱发爆轰波，只会形成明显的燃烧带与激波非耦合结构；当钝头体的尺度较大时，流场中不会出现燃烧带与激波的非耦合现象，且这一特征与马赫数无关.通过调整球体直径，获得了激波和燃烧带部分耦合的燃烧流场结构，这一流场结构在楔面诱发的斜爆轰波中并不存在，说明稀疏波与爆轰波面的相互作用是决定圆球诱发斜爆轰波的关键.

斜爆轰，氢氧，起爆，熄爆

引言

航空航天技术的飞速发展，要求飞行器的飞行速度越来越快，传统涡轮发动机难以满足要求.爆轰发动机是现代发动机的一种选择，主要有斜爆轰发动机、脉冲爆轰发动机[1]和旋转爆轰发动机[2]三种.斜爆轰发动机具有结构简单、质量小、成本低、比冲高和适应高速飞行的优点，在飞行马赫数高达10的情况下仍能为飞行器提供有效推力[3].

爆轰波是一种强激波诱发快速燃烧的耦合流场结构.爆轰波相对于预混气体是超声速传播的.在实验室坐标系下，如果预混气以超声速流动，理论上可能实现爆轰波的驻定.但实际上，驻定正爆轰波很难实现，而利用斜激波诱发的斜爆轰波已经成功实现了驻定[4]，并且相对正爆轰波具有较小的总压损失.而关于斜爆轰燃烧技术在高超声速飞行器动力装置中的应用，近年来已经开展了多方面相关的基础理论与工程应用的数值和实验研究[5-9].

斜爆轰发动机不仅具有超燃发动机的优点，且燃烧过程近似等容，相比于超燃发动机的等压燃烧具有更高的燃烧效率.Li等[10]通过数值分析发现，楔面诱发的斜爆轰流场结构包含斜激波、反应区、爆燃波以及之后的斜爆轰波.Viguier等[11]在1997年用实验证实了这一斜爆轰波结构的存在.Choi等[12]发现，斜爆轰波也有如正爆轰波一般的胞格结构，不过与正爆轰波不同的是：斜爆轰波三波点的轨迹不存在相互交叉.Teng等[13-15]对斜爆轰波胞格结构的变化规律给出了定量化的研究结论.斜爆轰发动机的来流马赫数范围广，来流条件复杂.Zhang等[16-17]研究发现来流当量比不同会影响斜爆轰波的特征长度.为了实现斜爆轰波的迅速、可靠引发，可以采用钝头体进行诱发.但是钝头体诱发斜爆轰波的过程比楔面诱发更加复杂，涉及到正激波、斜激波和稀疏波与燃烧反应的复杂相互作用，其中的波系结构和起爆机理还不清楚.

1972年，Lehr[18]利用弹状钝头体诱发爆轰波，试验结果表明钝头体诱发的爆轰波与楔面诱发的爆轰波在流场结构上有极大的不同.Kaneshige等[19]在1996年的氢气/氧气爆轰实验中发现，利用球形钝头体诱发爆轰波的过程出现了伴随着燃烧带的激波、爆燃转爆轰、直接起爆等三种流场结构.Ju等[20]通过理论和数值研究对起爆给出了定量化规律.Maeda等[21-25]在乙炔/氧气爆轰实验中，得到了燃烧带与激波非耦合、草帽状斜爆轰波、直接起爆爆轰波等三种典型的球形钝头体诱发斜爆轰波的流场结构.在非楔面斜爆轰中，由于稀疏波的存在，对斜爆轰波结构带来一定的影响[26].

本文在温度300K,压力50kPa的初始条件下，改变来流马赫数和球体直径，研究球形钝头体诱发爆轰波的过程以及不同流场结构之间的临界条件.通过数值模拟给出流场结构随马赫数与球体直径变化的分布图.

1 数学与物理模型

本研究对圆球诱发斜爆轰波过程进行数值模拟.预混可燃气体在超声速流动下，在钝头体顶部形成脱体激波，驻点附近为高温高压区.不同的来流马赫数及钝头体大小决定了驻点附近的流场参数，驻点附近的流场参数又决定了下游的流场结构.当驻点附近得到的能量高过诱发爆轰波的临界能量时，流场出现爆轰波；反之，则不能.

现有的数值研究表明，在超声速流动下，黏性对流场影响很小，数值研究主要采用无黏流假设.因此控制方程简化为二维轴对称Euler方程

其中，ωi为化学反应中第i种组分的质量生成率，它由化学反应模型决定.在上述方程中ρi是第i种组分的密度，总密度和v分别代表x方向和r方向上的速度.比内能e的计算公式为

其中，Ri代表第i种组分的气体常数，T为温度.

数值模拟中采用自适应无结构的四边形网格[27]和MUSCL-Hancock算法[28]，Rienmann问题的求解采用HLLC算子.时间项为由CFL数控制的显式格式，CFL数为0.4.化学反应采用时间分裂算法.氢/空气化学反应模型[29]含有11种组元(H2,O2,O,H,OH, HO2,H2O2,H2O,N2,N,NO)和23个基元化学反应，详细的机理见Chemkin程序说明，化学反应带来的刚性问题由DVODE软件包解决[30].氢气/空气混合物的化学当量比φ=1.0，即H2:O2:N2=2:1:4.球体边界为固壁边界，下边界为轴对称边界，左边界为来流，其他边界为自由边界.

2 数值结果与讨论

参考Lehr[18]在1972年的实验，给定流场的初始温度为300K，压力为50kPa.实验中的钝头体直径为15mm，为研究钝头体尺寸对于爆轰的影响，选取球形直径D=5～15mm.由于球形直径的不同，初始网格尺寸也不同，最大网格不大于1.2mm，加密后最小网格不大于0.1mm，此时继续加密网格流场结构不发生变化.图1所示为Ma=4.0,D=5mm，加密层数为3和4，即最小网格为0.125mm和0.06mm时的网格加密情况.可以看到加密处基本不变.图2为这两种加密情况下，球头前的压力温度分布曲线，图中的数字表示加密层数，两种加密网格的压力和温度曲线基本重合.

图1 Ma=4.0，D=5mm，加密层数3(左)，4(右)网格分布Fig.1Ma=4.0，D=5mm,encryption layer 3(left)and 4(right)grid distribution

图2 Ma=4.0，D=5mm，y=0时的压力和温度变化Fig.2Ma=4.0，D=5mm,the pressure and temperature at liney=0 variation diagram

2.1 算例与讨论

算例1 球体直径D=5mm

当球体直径为5mm时，马赫数的增大不会引发爆轰波.当马赫数为4.0的时候，驻点附近温度达到自燃温度，气体燃烧并向下游传播.由于钝头体为球形，气体先受到球体的压缩作用，后受到膨胀作用.由于稀疏波的影响，燃烧熄灭，钝头体下游不会出现燃烧带，如图3.

图3 Ma=4.0,D=5mm，压力(上)、温度(下)分布图Fig.3Ma=4.0,D=5mm,pressure(upper)and temperature(lower) distribution

当马赫数为5.0的时候，驻点附近的燃烧延续到下游流场，形成燃烧带与激波非耦合的现象，如图4.

图4 Ma=5.0,D=5mm，压力(上)、温度(下)分布图Fig.4Ma=5.0,D=5mm,pressure(upper)and temperature(lower) distribution

继续增大来流速度直至马赫数为15，流场结构都为燃烧带与激波非耦合的情况，马赫数较大时只会导致燃烧带变宽，流场结构不会出现变化.

为界定燃烧熄灭与非耦合情况的临界马赫数，分析马赫数为4.12和4.14时的流场结构.结果表明当马赫数大于4.1时，出现燃烧带与激波非耦合的情况，且最大温度也会有较大的变化，如图5所示.

可以看到，当马赫数从4.1变化到4.12的时候，不仅流场结构出现变化，驻点温度大幅升高，因此当精度为0.02时，临界马赫数等于4.1.

图5 最高温度随马赫数的变化Fig.5 Maximum temperature with di ff erent Mach numbers

算例2 球体直径D=15mm

当球体直径为 15mm时，流场结构与直径为5mm时的不完全一样.马赫数为4.0的时候，流场结构与直径为5mm时相似，但是在马赫数为5.0的时候，直径15mm钝头体可以直接引发爆轰波，如图6所示.

图6 D=15mm，压力(上)、温度(下)分布图Fig.6D=15mm,pressure(upper)and temperature(lower) distribution

为研究在马赫数4.0～5.0之间的流场情况，马赫数变化步长设为0.1.结果发现，在马赫数为4.1的时候，流场并不出现燃烧带与激波非耦合的情况，而是在驻点位置引发爆轰波，并诱发下游流场形成爆轰波.由于C-J爆速大于此时的来流马赫数，爆轰波无法驻定，向上游传播.

对比直径为5mm和15mm的结果可以看到，在爆轰燃烧中，在同样的流场条件下，钝头体的尺度是诱发爆轰波的一个重要参数.

2.2 流场结构随来流马赫数和球体直径变化的分布规律

来流马赫数和球体钝头体直径都是影响流场结构的关键参数.为了深入研究其对流场的影响，设置直径变化步长为2.5mm，马赫数变化步长为0.1.流场结构如图7所示.

从图7可以看到，马赫数的变化对于流场是否燃烧的影响较大.当马赫数从4.0变化到4.2之后，流场燃烧熄灭的情况都消失了.而球体直径的大小对于诱发爆轰波的影响较大.当球体直径为5mm时，钝头体无法诱发爆轰波；当球体直径增大到7.5mm时，钝头体才能诱发爆轰波.

图7 马赫数和圆球直径对起爆的影响Fig.7 The influenc of Mach number and ball diameter on initiation

对马赫数为4.2，直径为5～7.5mm进行进一步的数值研究发现，如果直径变化步长为0.5mm,当直径小于6.5mm时，流场结构为燃烧带与激波非耦合；当直径大于6.5mm时，预混气体在经过一段时间后起爆；当直径为6.5mm时，驻点位置的燃烧改变弓形激波在驻点附近的形状，但是不足以克服稀疏波的影响，对激波的影响不能传播到整个流场，这一流动是定常的.可以看到激波面由两道斜激波构成，在y=20mm处存在一个明显的波面拐点，如图8所示.拐点上游是斜爆轰波，下游是斜激波，燃烧面与激波面发生了明显的解耦.这一现象在楔面诱发的斜爆轰流场[13-16]中是不存在的.与楔面诱发斜爆轰波不同的是，气体在流经圆球时，不仅受到激波压缩作用，亦受到稀疏波的作用，圆球诱发的稀疏波是导致解耦的原因.因此可以得出结论，稀疏波与爆轰波面的相互作用是决定圆球诱发斜爆轰波的关键.

图8 Ma=4.2，D=6.5mm，压力(上)、温度(下)分布图Fig.8Ma=4.2,D=6.5mm,pressure(upper)and temperature(lower) distribution

3 结论

本文采用数值模拟对圆球诱发斜爆轰波进行了研究，重点分析了不同球体直径对起爆的影响.当球体直径为5mm时，流场结构不会出现爆轰现象.随着马赫数的增加，流场从燃烧熄灭变到燃烧带与激波非耦合的结构，临界马赫数为4.1；当球体直径为15mm时，流场中能够观测到爆轰现象.随着马赫数的增加，流场从燃烧熄灭变到诱发爆轰波，但是不会出现燃烧带与激波非耦合的情况；当球体直径在5～15mm之间，马赫数高于4.1的时候，燃烧熄灭的情况不再出现，球体直径对于爆轰波引发的影响明显，且当直径大于6.5mm的时候开始出现爆轰现象.

以前的学者对楔面诱发斜爆轰波进行了较多的研究，但是对球体诱发斜爆轰波还缺乏深入的研究.这种斜爆轰波的引发过程更加复杂，涉及到正激波/爆轰波、稀疏波与斜激波/爆轰波的相互作用.由于不同的直径会诱发不同强度的稀疏波，它们与燃烧的相互作用成为能否成功起爆的关键.本文的研究证实大直径圆球更容易诱发爆轰波，和以前的起爆实验与理论[17]结果也是定性一致的.下一步需要在此基础上，对临界状态下膨胀波与激波诱发燃烧的相互作用进行深入、系统的研究，从而建立适应工程应用的起爆模型.

1 董刚，范宝春，谢波.氢气--空气混合物中瞬态爆轰过程的二维数组模拟.高压物理学报，2004，18(1)：40-46(Dong Gang,Fan Baochun,Xie bo.Two-dimensional simulation of transient detonation process for H2-O2-N2mixture.Chinese Journal of High Pressure Physics,2004，18(1)：40-46(in Chinese))

2 LinW,ZhouJ,LiuS,etal.Experimentalstudyonpropagationmode of H2/air continuously rotating detonation wave.International Journal of Hydrogen Energy,2015,40(4):1980-1993

3 Voland RT,Huebner LD,McClinton CR.X-43A hypersonic vehicle technology development.Acta Astronautica，2006,59(1-5)：181-191

4 林志勇，周进，张继业等.预混超声速气流斜激波诱发脱体爆轰研究.航空动力学报，2009，24(1)：50-54(Lin Zhiyong,Zhou Jin, Zhang Jiye,et al.Investigation of detached detonation induced by oblique shock in premixed supersonic fl w.Journal of Aerospace Power,2009，24(1)：50-54(in Chinese))

5 Rudy W,Dziubanii K,Zbikowski M,et al.Experimental determination of critical conditions for hydrogen-air detonation propagation in partially confine geometry.International Journal of Hydrogen Energy,2016,1-8

6 Wang C,Dong XZ,Shu CW.Parallel adaptive mesh refinemen method based on WENO finit di ff erence scheme for the simulation ofmulti-dimensional detonation.Journalof Computational Physics, 2015,298:161-175

7 Zhang Z,Li Z,Dong G.Numerical studies of multi-cycle acetyleneair detonation induced by shock focusing.Procedia Engineering, 2015,99:327-331

8 王成，宁建国，雷娟.障碍物对氢氧预混气体爆轰波传播的影响//庆祝中国力学学会成立50周年暨中国力学学会学术大会，北京，2007(Wang Cheng,Ning Jianguo,Lei Juan.The influenc of obstaclesonthepropagationofH2-airdetonation//Celebrationofthe 50 Anniversary of the Founding of the Chinese Society of Mechanics and the Academic Conference of Chinese Society of Mechanics, Beijing,2007(in Chinese))

9 Wang C,Dong XZ,Shu CW.Parallel adaptive mesh refinemen method based on WENO finit di ff erence scheme for the simulation ofmulti-dimensional detonation.Journalof Computational Physics, 2015,298:161-175

10 Li C,Kailasanath K,Oran ES.Detonation structures behind oblique shocks.Phys Fluids,1994,6:1600-1611

11 Viguier C,Silva LFFD,Desbordes D,et al.Onset of oblique detonation waves:comparison between experimental and numerical results for hydrogen–air mixtures.Proc Combust Inst,1997,26:3023-3031

12 Choi JY,Kim DW,Jeung IS,et al.Cell-like structure of unstable oblique detonation wave from high-resolution numerical simulation.Proc Combust Inst,2007,31:2473-2480

13 Teng HH,Jiang ZL,Ng HD.Numerical study on unstable surfaces ofobliquedetonations.JournalofFluidMechanics,2014,744:111-128

14 Teng HH,Ng HD,Li K,et al.Evolution of cellular structures on oblique detonation surfaces.Combustion and Flame,2015,162(2): 470-477

15 Teng HH,Jiang ZL.On the transition pattern of the oblique detonation structure.Journal of Fluid Mechanics,2012,713:659-669

16 Zhang Y,Gong J,Wang T.Numerical study on initiation of oblique detonations in hydrogen–air mixtures with various equivalence ratios.Aerospace Science and Technology,2016,49:130-134

17 Wang T,Zhang YN,Teng HH,et al.Numerical study of oblique detonation wave initiation in a stoichiometric hydrogen-air mixture,Physics of Fluids,2015,27(9):096101

18 Lehr HF.Experiments on shock-Induced combustion.Astronautica Acta,1972,17:589-597

19 Kaneshige MJ,Shepherd JE.Oblique detonation stabilized on a hypervelocity projectile//26th Symp.(Int.)on Combustion,Pittsburgh, 1996.3015

20 Ju Y,Masuya G,Sasoh A.Numerical and theoretical studies on detonation initiation by a supersonic projectile.Symposium on Combustion,1998,27(2):2225-2231

21 Maeda,S,Inada,R,Kasahara,J,et al.Visualization of the nonsteady state oblique detonation wave phenomena around hypersonic spherical projectile.Proc Combust Inst,2011,33:2343-2349

22 Maeda S,Kasahara J,Matsuo A.Oblique detonation wave stability around a spherical projectile by a high time resolution optical observation.Combustion and Flame,2012,159(2):887-896

23 Maeda S,Sumiya S,Kasahara J,et al.Initiation and sustaining mechanisms of stabilized oblique detonation waves around projectiles.Proceedings of the Combustion Institute,2013,34(2):1973-1980

24 Maeda S,Sumiya S,Kasahara J,et al.Scale e ff ect of spherical projectiles for stabilization of oblique detonation waves.Shock Waves, 2015,25(2):141-150

25 Choi JY,Maeda S,Kasahara J,et al.Calculation of drag coefficients for hypersonic spherical projectiles initiating oblique detonation wave or shock-induced combustion//50th AIAA Aerospace Sciences Meeting Including the New Horizons Forum and Aerospace Exposition.2012

26 Liu Y,Han X,Yao S,et al.A numerical investigation of the prompt oblique detonation wave sustained by a finite-lengt wedge.Shock Waves,2016:1-11

27 SunM,TakayamaK.Conservativesmoothingonanadaptivequadrilateral grid.Journal of Computational Physics,1999,150:143-180

28 Toro EF.Riemann Solvers and Numerical Methods for Fluid Dynamics.(Second ed).Berlin:Springer,1999

29 Kee RJ,Rupley FM,Meeks E,et al.Chemkin-III:a fortran chemical kinetics package for the analysis of gas-phase chemical and plasma kinetics.UC-405,SAND96-8216,Sandia National Laboratories, 1996

30 BrownPN,ByrneGD,HindmarshAC.VODE:avariable-coefficient ODE solver.SIAM J Sci Stat Comput,1989,10(5):1038-1051

NUMERICAL STUDY OF THE OBLIQUE DETONATION INITIATION INDUCED BY SPHERES1)

Fang Yishen Hu Zongmin Teng Honghui2)Jiang Zonglin

(State Key Laboratory of High-Temperature Gasdynamics,Institute of Mechanics，Chinese Academy of Sciences，Beijing100190，China)

(School of Engineering Science，University of Chinese Academy of Sciences，Beijing100049，China)

The oblique detonation wave engine is a new kind of engine which has a simple structure,low cost,and high specifi impulse.In order to ensure the initiation,blunt body is used to induce the oblique detonation wave.The oblique detonation wave fl w fiel induced by spheres in supersonic hydrogen/air mixture is numerically simulated,based on the Euler equations and a detailed hydrogen-oxygen chemical reaction model.Unlike the oblique detonation wave induced by a wedge,the reacting fl w around a sphere is much more complex.First,a normal shock wave/detonation wave is formed, then oblique shock wave/detonation wave is developed in the presence of a rarefaction wave.The numerical simulation results show that after the gases being compressed by the blunt body and reaching the auto-ignition temperature,two kinds of fl wfiled will appear.When Mach numbers are low,the combustion will be quenched and can not appear downstream of the blunt body due to the influenc of the rarefaction wave.When Mach numbers are high,combustion can spreadto the downstream region.When the scales of blunt body are small,energy around the stationary point is not enough to induce detonation initiation and an obvious decoupling of combustion and shock wave is formed.As the sphere becomes large enough,decoupling of combustion and shock wave will not appear in the fl w and this feature is indpendent of the Mach number.By adjusting the spheric diameter,the fl w structures with partial coupling of shock wave and combustion zone was obtained which does not exist in a wedgy-induced oblique detonation.The present investigations suggest that the interaction between rarefaction wave and detonation wavefront is the key issue for detonation initiation induced by a spheric body.

oblique detonation,oxyhydrogen,initiation,quench

O381

A

10.6052/0459-1879-16-143

2016–05–25收稿，2016–12–09录用,2016–12–09网络版发表.

1)国家自然科学基金资助项目(11372333,91641130,11532014).

2)滕宏辉，副研究员，主要研究方向：激波与爆轰物理，爆轰燃烧推进，多相燃烧与爆炸.E-mail：hhteng@imech.ac.cn

方宜申，胡宗民，滕宏辉，姜宗林.圆球诱发斜爆轰波的数值研究.力学学报，2017,49(2):268-273

Fang Yishen,Hu Zongmin,Teng Honghui,Jiang Zonglin.Numerical study of the oblique detonation initiation induced by spheres.Chinese Journal of Theoretical and Applied Mechanics,2017,49(2):268-273