文章快速检索     高级检索
  气体物理  2018, Vol. 3 Issue (1): 20-27   DOI: 10.19527/j.cnki.2096-1642.2018.01.003
0

引用本文  

刘朋欣, 李沁, 张涵信. 初始点火条件引起的三维旋转爆震波单/双波传播模式[J]. 气体物理, 2018, 3(1): 20-27.
Liu P X, Li Q, Zhang H X. Single-wave mode and double-wave mode caused by initiation condition in three-dimensional rotating detonation wave[J]. Physics of Gases, 2018, 3(1): 20-27.

基金项目

国家重点基础研究发展计划(2014CB744100);国家自然科学基金(91541105)

作者简介

刘朋欣(1990-)男, 硕士, 研究实习员, 主要研究方向为高阶差分格式、爆震燃烧数值模拟.E-mail:liupengxin2008@163.com

通信联系人

李沁(1969-)男, 博士, 研究员, 主要研究方向为计算空气动力学.E-mail:qin-li@vip.tom.com

文章历史

收稿日期:2017-09-04
修回日期:2017-10-23
初始点火条件引起的三维旋转爆震波单/双波传播模式
刘朋欣1,2, 李沁1,2, 张涵信1,2     
1. 中国空气动力研究与发展中心空气动力学国家重点实验室,四川绵阳 621000;
2. 中国空气动力研究与发展中心计算空气动力研究所,四川绵阳 621000
摘要:在不同工况下,旋转爆震波能够以单波、双波、多波模式进行传播.但在同一工况下,是否存在不同模式的稳定传播爆震波还有待进一步研究.基于Euler方程,耦合氢气/空气的有限化学反应速率模型,并采用高分辨率的5阶有限差分格式WENO-PPM5离散对流项,对三维旋转爆震波进行了数值模拟.计算结果表明,在同一特定工况下,旋转爆震波能够以两种不同的传播模式稳定传播,即单波模式和双波模式.详细地对比了两种传播模式下的流场特征、爆震波传播特性、推力性能等.在同一工况下,两种传播模式的爆震波周向传播速度相差不多,但双波模式的频率约为单波模式的2倍;双波模式下质量流量、比冲、推力的平均值均略高于单波模式;且双波模式的可燃混气层高度约为单波模式的1/2,这有助于缩小旋转爆震发动机的长度,使之更加紧凑.
关键词旋转爆震    点火条件    传播模式    传播特性    推力性能    
Single-Wave Mode and Double-Wave Mode Caused by Initiation Condition in Three-Dimensional Rotating Detonation Wave
LIU Peng-xin1,2, LI Qin1,2, ZHANG Han-xin1,2     
1. State Key Laboratory of Aerodynamics, China Aerodynamics Research and Development Center, Mianyang 621000, China;
2. Computational Aerodynamics Institute, China Aerodynamics Research and Development Center, Mianyang 621000, China
Abstract: Based on Euler equations associated with finite-rate chemical reaction mechanism of hydrogen/air mixtures, numerical calculations were performed by using a fifth-order WENO-PPM5 scheme for spatial discretization to simulate the rotating detonation in an annular chamber. Results show that under specific condition, the rotating detonation wave can propagate in two different modes, i.e. single-wave mode and double-wave mode. Comparison was made between the two modes, including the characteristic of flowfield, the feature of detonation propagation, the propulsion performance. Under the same condition, the circumferential velocity of detonation propagation in two modes is approximately equal, and the frequency of double-wave mode is twice as that of single-wave mode; the time-averaged mass flux, specific impulse and thrust of double-wave mode are all a little higher than that of single-wave mode. The height of fresh mixture in two-wave mode is about a half of that in single-wave mode, which is favorable to shorten the length of rotating detonation engine.
Key words: rotating detonation    initiation condition    propagation mode    propagation feature    propulsion performance    
引言

爆震燃烧是一种激波诱导的燃烧, 具有很多的优势, 如更高的热效率, 更快的反应速度.目前, 有多种基于爆震燃烧技术的爆震发动机形式, 其中旋转爆震发动机(rotating detonation engine, RDE)由于理论上只需要一次点火即可连续工作, 且推力相对稳定, 近年来受到了较多的关注.

Nicholls等[1]和Voitsekhovskii[2]较早地对旋转爆震开展了理论和实验研究.随后的几十年里, 各国学者进一步研究了影响RDE工作的各种因素[3-10], 包括点火方式、喷注方法、燃料和氧化剂、燃烧室形状等. RDE中的爆震波能够以不同的模式传播. Bykovskii等[3]在总结了大量实验数据并结合理论分析后, 给出了估算爆震波数目的经验公式. Wolański[4]则指出旋转爆震波的传播模式受到实验工况、燃烧室尺寸、可燃混气的喷注速率、可燃混气的成分等因素影响, 并从传播周期和燃料喷注速率的角度给出了决定爆震波数目的关系式. Zhou等[6]和Wu等[10]在不同工况条件下得到了单波、双波、多波的传播模式; 并通过数值模拟发现, 随着喷注总压增大, 质量流量随之增大, 单波模式的流场变得不稳定; 当总压增大到一定程度时, 单波模式会出现周期性震荡且强弱交替变化; 而在同一工况下, 双波模式则可以稳定传播.虽然已经发现旋转爆震波的传播模式受到诸多因素的影响, 但在同一工况条件下是否存在不同的稳定传播模式还须进一步研究.

本文通过对三维旋转爆震波进行数值模拟, 在同一工况下, 得到了爆震波的两种稳定传播模式(即单波模式和双波模式), 并对比分析了流场特征、推力性能等, 有助于进一步认识旋转爆震波的传播机制.

1 计算方法及模型 1.1 控制方程及数值方法

在爆震波传播过程中, 输运现象的影响一般小于对流作用, 因此大部分数值计算中忽略了输运特性, 包括黏性、热传导、物质扩散[4].假设混合气体为理想气体, 并忽略输运特性的影响, 且不考虑体积力及热辐射.使用含化学反应的三维Euler方程作为控制方程, 其在曲线坐标系下的形式为

$ \frac{{\partial \mathit{\boldsymbol{Q}}}}{{\partial t}} + \frac{{\partial \mathit{\boldsymbol{F}}}}{{\partial \xi }} + \frac{{\partial \mathit{\boldsymbol{G}}}}{{\partial \eta }} + \frac{{\partial \mathit{\boldsymbol{H}}}}{{\partial \zeta }} = \mathit{\boldsymbol{S}} $ (1)

其中,

$ \mathit{\boldsymbol{Q}} = {\mathit{\boldsymbol{J}}^{-1}}\left( {\begin{array}{*{20}{c}} \rho \\ {\rho u}\\ {\rho v}\\ {\rho w}\\ E\\ {\rho {f_i}} \end{array}} \right), \mathit{\boldsymbol{F}} = {\mathit{\boldsymbol{J}}^{-1}}\left( {\begin{array}{*{20}{c}} {\rho \bar U}\\ {\rho \bar Uu + P{\xi _x}}\\ {\rho \bar Uv + P{\xi _y}}\\ {\rho \bar Uw + P{\xi _z}}\\ {\left( {E + P} \right)\bar U}\\ {\rho \bar U{f_i}} \end{array}} \right) $
$ \mathit{\boldsymbol{G}} = {\mathit{\boldsymbol{J}}^{-1}}\left( {\begin{array}{*{20}{c}} {\rho \bar V}\\ {\rho \bar Vu + P{\eta _x}}\\ {\rho \bar Vv + P{\eta _y}}\\ {\rho \bar Vw + P{\eta _z}}\\ {\left( {E + P} \right)\bar V}\\ {\rho \bar V{f_i}} \end{array}} \right), \;\;\mathit{\boldsymbol{H}} = {\mathit{\boldsymbol{J}}^{-1}}\left( {\begin{array}{*{20}{c}} {\rho \bar W}\\ {\rho \bar Wu + P{\zeta _x}}\\ {\rho \bar Wv + P{\zeta _y}}\\ {\rho \bar Ww + P{\zeta _z}}\\ {\left( {E + P} \right)\bar W}\\ {\rho \bar W{f_i}} \end{array}} \right) $
$ \mathit{\boldsymbol{S}} = {\mathit{\boldsymbol{J}}^{-1}}\left( {\begin{array}{*{20}{c}} 0\\ 0\\ 0\\ 0\\ 0\\ {{S_i}} \end{array}} \right), {\rm{ }}{\mathit{\boldsymbol{J}}^{-1}} = \left| {\frac{{\partial \left( {x, y, z} \right)}}{{\partial \left( {\xi, \eta, \zeta } \right)}}} \right| $
$ \begin{array}{l} \bar U = u{\xi _x} + v{\xi _y} + w{\xi _z}\\ \bar V = u{\eta _x} + v{\eta _y} + w{\eta _z}\\ \bar W = u{\zeta _x} + v{\zeta _y} + w{\zeta _z}\\ \;\;\;\;\;\;P = \rho RT \end{array} $ (2)
$ E = \rho \sum\limits_i {{f_i}{h_i}}-P + \frac{1}{2}\rho ({u^2} + {v^2} + {w^2}) $ (3)
$ \begin{array}{l} {S_i} = {M_i}\sum\limits_{j = 1}^{nr} {({\beta _{ji}}- {\alpha _{ji}})} \left[{{K_j}\prod\limits_{i = 1}^{ns} {{{\left( {\frac{{{\rho _i}}}{{{M_i}}}} \right)}^{{\alpha _{ji}}}}{{\left( {\sum\limits_{k = 1}^{nr} {\frac{{{\rho _k}}}{{{M_k}}}{C_{jk}}} } \right)}^{{L_j}}}} } \right.-\\ \;\;\;\;\left. {{K_{-j}}\prod\limits_{i = 1}^{ns} {{{\left( {\frac{{{\rho _i}}}{{{M_i}}}} \right)}^{{\beta _{ji}}}}{{\left( {\sum\limits_{k = 1}^{nr} {\frac{{{\rho _k}}}{{{M_k}}}{C_{jk}}} } \right)}^{{L_j}}}} } \right] \end{array} $ (4)

其中, J-1为坐标变换矩阵的行列式; ρP分别为混合物的密度和压力; u, v, w为分速度; fi为第i种组分的质量分数(i=1…ns-1, ns为总组分数); E为单位体积的总内能; hi为比焓, 通过温度的多项式拟合得到[11]; Si为第i种组分的生成速率; Mi为分子量; αjiβji为第i种组分在第j个基元反应式中的等当量系数; nr为基元反应方程的总数; Lj代表第j个基元反应是否含有三体碰撞, Cjk代表组分的碰撞效率; KjK-j分别代表正逆反应的反应速率, 可以通过Arrhenius公式得到.

通过求解方程(1)~(4), 就得到耦合了ns组分nr基元反应的Euler方程.本研究采用氢气/空气混合物的7组分(H2, O2, H, O, OH, H2O, N2)8基元有限速率化学反应模型进行计算[12].

控制方程的时间离散采用具有总变差不增(total variation diminishing, TVD)性质的3阶Runge-Kutta法, 对流项采用Steger-Warming矢通量分裂, 并使用5阶WENO-PPM5格式[13]进行空间离散, 化学反应源项的求解采用耦合显式方法.

1.2 方法验证

首先通过一维爆震波行进问题来验证化学反应模型和数值方法的正确性.计算采用300 mm长的爆震管, 管内充满等当量比的氢气/空气混合气体(H2:O2:N2=2:1:3.76).初始时管内气体状态为P=2 atm, T=400 K, 在爆震管封闭端放置一小段高温高压区域作为起爆区: [0, 5 mm], P=30 atm, T=6 000 K.

计算采用了3种尺度的均匀网格, Δx=0.1, 0.3, 0.5 mm.时间步长为Δt=10-3 μs, 计算至t=120.0 μs的温度和压力如图 1所示.在所选网格尺度下, 压力和温度间断面是耦合在一起的, 准确地捕捉到了爆震波的基本特征.起爆区附近的震荡是由于采用高能起爆方式初始条件存在强间断引起的, 文献[16]中也有类似分布出现. 表 1将计算得到的爆震波整体性质与理论值[14]进行了对比, 包括爆震波传播速度UD, Chapman-Jouguet理论点处的压力PCJ和温度TCJ.可以看到, 计算值与理论值均符合良好, 也说明了计算程序能够准确模拟爆震波行进问题.

图 1 t=120.0 μs, 一维爆震波温度和压力分布 Fig.1 Temperature and pressure distributions of one-dimensional detonation wave at t=120.0 μs
下载CSV 表 1 一维爆震波整体性质对比 Tab.1 Comparison of properties of 1D detonation
1.3 计算模型及边界条件

计算模型采用旋转爆震波发动机典型的同轴圆环腔燃烧室.燃烧室内半径为40 mm, 外半径为50 mm, 长度为100 mm.燃料为等当量比的氢气/空气混合气体.喷注总温T0和总压P0分别为300 K和10 atm.可燃混合气体通过喷管喷入燃烧室, 形成一个或多个可燃混合气体层, 沿周向连续传播的爆震波则会燃烧消耗此可燃混合气体, 并将燃烧产物从出口一侧排出.

圆环内外壁采用无黏滑移壁面边界; 周向边界采用周期性边界处理; 计算域出口采用无反射边界条件, 出口背压为1.0 atm; 在实验中, 燃料和氧化剂通常是通过孔型或缝型喷注器分别进行填充的.数值模拟中忽略混合过程及喷注细节, 假定均匀可燃混合气体通过Laval喷管喷注进入燃烧室, 并设定Laval喷管的出口面积与喉道面积比为Aw/Athroat=10.当给定来流总压P0, 总温T0时, Laval喷管出口处的压力存在两个临界值:喷注达到壅塞状态时的临界压力Pcr1, 可燃混合气体以超声速喷入燃烧室时的临界压力Pcr2.根据入口边界处压力Pw的不同, 喷注可以分为以下4种情况:

(1) 堵塞(PwP0), 可燃混合气体不能进入燃烧室, 入口边界处理为无滑移壁面;

(2) 亚声速喷注(Pcr1PwP0), 可燃混合气体以亚声速喷入燃烧室, 各热力学参数利用等熵关系计算;

(3) 壅塞、亚声速喷注(Pcr2PwPcr1), 喷管内流动达到壅塞状态, 可燃混合气体以亚声速喷入燃烧室, 各热力学参数利用激波关系计算;

(4) 壅塞、超声速喷注(PwPcr2), 喷管内流动处于壅塞状态, 可燃混合气体以超声速喷入燃烧室, 各热力学参数利用等熵关系计算.

计算采用均匀网格, 总的网格数量为1 800×50×500(周向×径向×轴向), 网格尺度约为Δ=0.2 mm.初始时, 燃烧室内气体的状态为P=1 atm, T=300 K.在燃烧室的头部设置一定高度的等当量比的氢气/空气可燃混合气体层, 使用一段一维爆震波进行起爆, 并在其后方设置一段空气区作为隔离段.其余区域则全部填充为空气.整个燃烧室给定一个初始的轴向速度(200 m/s), 以便更快地形成旋转爆震波.本文设置了两个算例, 两者只是初始流场设置不同, 算例1的起爆区和隔离段的宽度约为算例2的一半, 初场设置的示意图见图 2.这样做的目的主要是改变初始爆震波前可燃混合气体层沿周向的长度; 考察不同波前长度对旋转爆震波的发展过程的影响.本文初始点火采用了直接起爆方式, 而没有考虑间接起爆方式中缓燃向爆震转变的过程(deflagration-to-detonation transition, DDT).为增强计算的稳定性且提高时间步长, 前期爆震波发展阶段使用低阶迎风格式进行计算; 当爆震波基本稳定后, 改用WENO-PPM5高精度格式.为建立旋转爆震流场, 不同计算可能会采用不同的计算格式和起爆方式[10, 16], 但数值起爆方式并不具有真实物理意义, 研究通常关注爆震建立稳定运行后的传播模式和流动结构.

图 2 算例初场设置 Fig.2 Initiation flowfield of two cases
2 结果与分析 2.1 稳定爆震波形成过程及流场结构

本文研究发现, 在特定工况下, 爆震波可以存在两种不同的稳定传播模式, 即单波模式和双波模式.由于初场设置的不同, 两算例中爆震波发展过程及稳定后的传播模式也不同.算例1中的稳定爆震波最终以单波模式传播, 而算例2则形成了稳定传播的双波爆震波.

算例1中, 稳定爆震波的形成过程较为简单:使用一维爆震波起爆后, 燃烧室内形成了一个单爆震波, 爆震波传播了几个周期以后即趋于稳定, 见图 3.而算例2中, 稳定爆震波的形成则较为复杂.使用一维爆震波起爆后, 形成了一个明显的爆震波锋面1, 且由于流场中的缓燃燃烧及压缩波向初始爆震波反方向的传播, 当能量累计到一定程度时, 会在局部形成热点(hot spot), 如图 4t1时刻的热点, 并发展成新的爆震波锋面; 在t2时刻可以清晰地看到新的爆震波锋面2.爆震波锋面1和2会发生碰撞, 并分别形成了新的爆震波锋面3和4;爆震波锋面4在传播过程中与局部高压区碰撞, 形成新的热点, 并进一步发展为新的爆震波锋面5, 见t3时刻.各爆震波锋面继续向前传播, 爆震波锋面3和4会发生碰撞; 碰撞后, 原爆震波锋面3形成新的热点1, 并逐渐发展成新的爆震波锋面6;而原爆震波锋面4碰撞后形成热点2, 接着继续与爆震波锋面5发生碰撞, 并被吸收, 见t4时刻和t5时刻.因此, 最后流场中的爆震波锋面只剩下5和6.初期, 爆震波锋面5和6之前的可燃混合气体层长度和高度都不相同; 经过几个传播周期的调整, 可燃混合气体层长度和高度趋于一致, 且爆震波以双波模式稳定传播, 见t6时刻.

图 3 算例1的稳定爆震波形成过程 Fig.3 Development processes of steady detonation in case 1
图 4 算例2的稳定爆震波形成过程 Fig.4 Development processes of steady detonation in case 2

刘世杰等[15]在实验中发现, 在同向传播模式下, 当燃烧室几何尺寸一定时, 爆震波数目主要受实验工况的影响, 随着推进剂总流量的增大而增多, 并通过改变氢气和空气的喷注质量流量, 得到了单波、混合单/双波、双波3种旋转爆震波传播模态.而本文研究则是关注同一工况下能否存在不同的稳定传播模式.在起爆初始条件设置方面, 有的数值模拟采用设置高温、高压区[16]或者采用两段一维爆震波[10]进行起爆, 也可以得到同向双波爆震波.这里不再对其他的点火条件进行计算测试和讨论, 而是关注同一工况下, 单、双波传播模式的异同.

图 5给出了两种传播模式在t=2 200.0 μs时的压力和温度分布.两模式均捕捉到了旋转爆震波流场的主要特征.以单波模式为例(图 5(b)), 可燃混合气体层D的存在, 使得爆震波A可以沿周向连续传播; 爆震波后产生的高温、高压产物沿轴向膨胀时, 与上一轮的燃烧产物相互作用, 就形成了斜激波B; 由于斜激波与爆震波后燃烧产物的热力学状态不同, 两者接触形成间断面C; 上一轮的高温燃烧产物会使接触面E处的可燃混气发生缓燃, 并消耗较少的可燃混合气体.

图 5 t=2 200.0 μs, 流场的压力和温度分布 Fig.5 Pressure and temperature distributions at t=2 200.0 μs

Bykovskii等[3]给出了估算爆震波数目的经验公式: Kdc/nh, 其中dc为圆环燃烧室外壁的直径, n为爆震波数目, h为波前可燃混合气体层的高度.对于所有使用气态氧化剂的同轴圆环形燃烧室, K=7±2.对于单波模式, n=1, h=36.9 mm, K=8.51;对于双波模式, n=2, h=18.3 mm, K=8.58;两者均符合经验公式.文献[3]还指出, 旋转爆震波前可燃混合气体层的高度h接近或超过临界值h*时, 爆震波才可以稳定传播, 并给出了使用胞格横向尺寸a来估算h*的关系式: h*≈(12±5)a.爆震波锋面是由许多三波点结构组成的, 当爆震波在壁面熏有烟膜的管道内传播时, 由于三波点附近的压力、涡量都较大, 三波点会在烟膜上留下运动痕迹, 其轨迹就是所谓的胞格结构[16].胞格结构的横向宽度即为胞格横向尺度a.本文通过记录各个点上在一个传播周期内的最大压力值[9], 得到了数值胞格.以单波模式为例, a≈3.7 mm, 则h*的最小值约为25.9 mm.当爆震波数继续增加时, 可燃混合气体层的高度会随之下降, 其偏离h*更大, 使本文计算工况下多波(n≥3)传播模式的存在更倾向于不可能.

2.2 爆震波整体传播特性

取采样点(θ=5.92 rad, r=45 mm, z=4 mm), 并记录采样点压力随时间的变化曲线, 结果如图 6所示.对于单波模式, 虽然前期(0~400.0 μs)流场中会存在局部的压力震荡, 但是波峰的周期性一直较好, 单波爆震波的传播很快趋于稳定.而双波模式中, 前期(0~400.0 μs)流场较为混乱, 峰值之间的时间间隔相差较大, 且存在局部压力震荡; 中期时(400.0~1 200.0 μs), 同向双波爆震波虽然已经形成, 但是由于波前可燃混合气体层的长度和高度不同, 而导致采样点出现峰值的时间间隔不一致; 经过几个周期的调整, 后期时两个爆震波的分布趋于一致, 采样点的峰值呈现出良好的周期性.

图 6 采样点压力随时间的变化 Fig.6 Temporal variations of pressure at sample point

取时间间隔1 500.0~2 200.0 μs, 计算爆震波的平均传播周期T和频率f, 并据此计算出爆震波的周向传播速度VD, 结果见表 2.两种模式下, 爆震波的周向传播速度十分接近; 且由于爆震波数目不同, 双波模式的频率约为单波模式的2倍.

下载CSV 表 2 两种传播模式下爆震波传播特征对比 Tab.2 Comparisons of propagation properties in two modes
2.3 推力性能

t=2 200.0 μs时, 观测入口中间线上(r=45 mm)的瞬时压力分布, 如图 7所示, 黑色线为喷管壅塞临界压力Pcr1(9.975 atm).当入口处的压力低于Pcr1时, 可燃混合气体喷入燃烧室, 并发生壅塞, 且质量流量为常数.单波模式发生壅塞流动的入口面积比约为0.917;双波模式则约为0.924;两者十分接近, 这与文献[10]中的结论一致.

图 7 入口中间线处的压力分布(t=2 200.0 μs, r=45 mm) Fig.7 Pressure distributions along middle radius at inlet

可燃混合气体质量流量mflux, 基于氢气的比冲Isp, 推力F随时间的变化见图 8.在观测的时间段内, 3个参数均围绕着一个固定值存在小幅度的震荡, 且双波模式的参数值基本上大于单波模式.

图 8 两种传播模式下参数随时间的变化 Fig.8 Temporal variations of parameters in two modes

计算此时间段内可燃混合气体的质量流量、比冲、推力的平均值, 结果见表 3.这3种参数的平均值相差不大, 双波模式均稍高于单波模式.双波模式的平均比冲比单波模式提高约1.86%, 可见双波模式的燃烧效率略高.两种模式下, 可燃混合气体的平均质量流量相差1.32%, 这也与两者相近的入口壅塞流动的面积比一致; 但是从图 5可以看出, 双波模式下, 可燃混合气体层的累积量大约只有单波模式的1/2, 且高度也只有后者的1/2.根据文献[3]的经验, 旋转爆震发动机的最小轴向长度应该大于波前可燃混合气体层高度的2倍.因此在同样的工况下, 双波模式所需的发动机轴向长度可以更小, 且燃料效率更高一些.

下载CSV 表 3 两种传播模式下参数平均值对比 Tab.3 Comparison of averaged parameters in two modes
3 结论与讨论

基于含化学反应源项的Euler方程, 采用氢气/空气的7组分8基元反应模型和高阶有限差分格式WENO-PPM5, 对三维旋转爆震波进行了数值模拟, 并得出如下结论:在同一工况下, (1)不同的点火条件可能引起旋转爆震波两种不同的稳定传播模式, 即单波模式和双波模式; (2)两模式中爆震波的周向传播速度相近, 且双波模式的频率约为单波模式的两倍; (3)双波模式下质量流量、比冲、推力的平均值均略大于单波模式; (4)双波模式中可燃混合气体层的高度约为单波模式的1/2, 这有助于减少旋转爆震发动机的轴向长度, 从而缩小发动机尺度, 使之更加紧凑.

虽然本文研究显示同一工况下旋转爆震波可能存在两种不同传播模式, 但导致具体模式的途径具有一定的不确定性.本文分析表明, 初始流场中预置爆震波前可燃混合气体层的长度不同, 是可能的因素之一; 是否还存在其他因素(如可燃混合气体层高度, 其他起爆方式)有待于进一步研究.

致谢: 此研究得到了国家自然科学基金(91541105), 国家重点基础研究发展计划(2014CB744100)的资助, 特此感谢.
参考文献
[1]
Nicholls J A, Dabora E K, Gealer R L. Studies in connection with stabilized gaseous detonation waves[J]. Symposium (International) on Combustion, 1958, 7(1): 766-772. DOI:10.1016/S0082-0784(58)80118-6
[2]
Voitsekhovskii B V. Maintained detonations[J]. Doklady Akademii Nauk UzSSR, 1959, 129(6): 1252-1256.
[3]
Bykovskii F A, Zhdan S A, Vedernikov E F. Continuous spin detonations[J]. Journal of Propulsion and Power, 2006, 22(6): 1204-1216. DOI:10.2514/1.17656
[4]
Wolański P. Detonative propulsion[J]. Proceedings of the Combustion Institute, 2013, 34(1): 125-158. DOI:10.1016/j.proci.2012.10.005
[5]
Lu F K, Braun E M. Rotating detonation wave propul-sion:experimental challenges, modeling, and engine con-cepts[J]. Journal of Propulsion and Power, 2014, 30(5): 1125-1142. DOI:10.2514/1.B34802
[6]
Zhou R, Wu D, Wang J P. Progress of continuously rotating detonation engines[J]. Chinese Journal of Aeronautics, 2016, 29(1): 15-29. DOI:10.1016/j.cja.2015.12.006
[7]
Bykovskii F A, Zhdan S A. Current status of research of continuous detonation in fuel-air mixtures (review)[J]. Combustion, Explosion, and Shock Waves, 2015, 51(1): 21-35. DOI:10.1134/S0010508215010025
[8]
Le Naour B, Falempin F, Coulon K. MBDA R & T effort regarding continuous detonation wave engine for propul-sion-status in 2016. Proceedings of the 21st AIAA International Space Planes and Hypersonics Technologies Conference, Xiamen, China: AIAA, 2017.
[9]
Cocks P A, Holley A T, Greene C B, et al. Develop-ment of a high fidelity RDE simulation capability. Kissimmee, Florida, USA: 53th AIAA Aerospace Scien-ces Meeting, 2015. http://arc.aiaa.org/doi/abs/10.2514/6.2015-1823
[10]
Wu D, Zhou R, Liu M, et al. Numerical investigation of the stability of rotating detonation engines[J]. Combus-tion Science and Technology, 2014, 186(10/11): 1699-1715.
[11]
Stull D R, Prophet H. JANAF thermochemical Tables. 2nd ed. Washington, DC, USA: U. S. Government Printing Office, 1971.
[12]
Evans J S, Schexnayder C J Jr. Influence of chemical kinetics and unmixedness on burning in supersonic hydro-gen flames[J]. AIAA Journal, 1980, 18(2): 188-193. DOI:10.2514/3.50747
[13]
Li Q, Liu P X, Zhang H X. Piecewise polynomial mapping method and corresponding WENO scheme with improved resolution[J]. Communications in Computational Physics, 2015, 18(5): 1417-1444. DOI:10.4208/cicp.150215.250515a
[14]
Yi T H, Lou J, Turangan C, et al. Propulsive perfor-mance of a continuously rotating detonation engine[J]. Journal of Propulsion and Power, 2011, 27(1): 171-181. DOI:10.2514/1.46686
[15]
刘世杰, 刘卫东, 林志勇, 等. 连续旋转爆震波传播过程研究(Ⅰ):同向传播模式[J]. 推进技术, 2014, 35(1): 138-144.
Liu S J, Liu W D, Lin Z Y, et al. Research on continuous rotating detonation wave propagation process (Ⅰ):one direction mode[J]. Journal of Propulsion Techno-logy, 2014, 35(1): 138-144. (in Chinese)
[16]
刘世杰. 连续旋转爆震波结构、传播模态及自持机理研究. 长沙: 国防科学技术大学, 2012.
Liu S J. Investigations on the structure, rotating mode and lasting mechanism of continuous rotating detonation wave. Changsha: National University of Defense Technology, 2012(in Chinese). http://cdmd.cnki.com.cn/Article/CDMD-90002-1013047748.htm