Analytical Models for the Thrust of a Rotating Detonation Engine

Joseph E. Shepherd
Professor of Aeronautics and Mechanical Engineering
joseph.e.shepherd@caltech.edu
Graduate Aerospace Laboratories
California Institute of Technology
Pasadena, CA 91125 USA
Jiro Kasahara
Professor of Aerospace Engineeringkasahara@nuae.nagoya-u.ac.jp
Nagoya University
Nagoya, Aichi, 464-8603, JAPANGALCIT Report FM2017.001
Graduate Aerospace Laboratories
California Institute of Technology
Pasadena CA 911125 USA

Submitted for Publication to the Journal of Propulsion and Power

September 11, 2017

Abstract
Two models are proposed for rotating detonation engine performance. The first model is motivated by models
of pulse detonation engine performance which are based on the pressure-time history within the detonation
tube. The present work extends those ideas to treat rotating detonation engines with a control volume analysis
that considers the forces within the combustion chamber. The key scaling parameters for this model are theheight of the reactant layer just ahead of the detonation wave and the computed Chapman-Jouguet pressureand velocity. The thrust can be estimated using these parameters and a simple functional form of the pressure
history on the injector surface. The second model is based on the approximation of mean axial flow anduses a conventional control volume analysis that focuses on exit conditions to evaluate thrust. The axial flow
speed and thrust are evaluated based on approximating the flow following the detonation as isentropic andconsidering a two-dimensional expansion that converts azimuthal motion into purely axial flow. Numericaland analytic evaluation of these models demonstrates that predicted thrust and specific impulse exhibit thesame scaling relationships with mixture properties as pulse detonation engines.

iContents
1

Introduction

1

2

Control Volume Analysis of Thrust

3

3

Pressure History Model
6
3.1 Evaluating the thrust . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11

4

Axial Flow Model
4.1 Modeling the expansion process . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.2 Analytical Expressions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.3 Axial Flow Model and Nozzles . . . . . . . . . . . . . . . . . . . . . . . . . . . .

13
14
18
20

5

Comparisons with simulations

22

6

Comparison to Test Results

22

7Comparison of RDE and PDE analysis

23

8Conclusion

25

A Equations of Motion

29

B Solutions to Triple-Point Configuration

32

C Existence of Attached Oblique Solutions

34

ii

List of Figures
1
2

3
4
5

6
7
8
9
10
11

12
13

a) A rotating detonation engine test fixture (Ishihara et al., 2015). b) Schematic of“unwrapped” flow field projected onto the θ - Z plane. . . . . . . . . . . . . . . .
Two approaches to formulating a control volume for a rotating detonation engine
to compute thrust base on a) chamber exit conditions. b) combustion chamber inlet
conditions. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Simplified model for detonation of a uniform layer of thickness H and detail of
triple point region (red box). . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Three examples of pressure-deflection diagrams for 3-5 interface mechanical matching solution. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Pressure at injector face (lines) Schwer and Kailasanath (2013) compared to the
exponential decay model (points) described in text. a) Stoichiometric C2 H4 -air andb) Stoichiometric C2 H4 -O2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
a) Specific thrust T / MÛ and b) fuel-based specific impulse Isp f of RDEs as predictedusing the pressure history model. . . . . . . . . . . . . . . . . . . . . . . . . . . .
Geometry of rotating detonation combustion chamber and associated axisymmetriccoordinate system. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Axial flow RDE model for C2 H4 -air: a) Pressure-enthalpy states; b) Specific thrust.
Û b) fuel-based specific
Sonic axial flow RDE model: a) Specific thrust T / M;
impulse Isp f . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Sonic axial flow RDE model: a) Effect of initial temperature T1 for P1 = 100 kPa;b) Effect of initial pressure P1 for T1 = 300 K. . . . . . . . . . . . . . . . . . . . .
a) Measured and model estimates of specific thrust T / MÛ for the C2 H4 -O2 RDE
tests of Kato et al. (2016). (a) Pressure history model and (b) sonic axial flow
model; lines are predictions, points are test data. . . . . . . . . . . . . . . . . . . .
Comparison of predicted RDE and PDE performance with measured performanceusing data from Fig. 10 of Rankin et al. (2017). a) H2 -air. b) C2 H4 -air. . . . . . . .
Pressure-flow deflection analyses. (a) Effect of air temperature (b) Effect of shear
flow for 3500 K air temperature. . . . . . . . . . . . . . . . . . . . . . . . . . . .

iii

2

3
6
8

10
12
15
17
18
18

24
25
35

List of Tables
1
2

Parameters for the cases reported in Schwer and Kailasanath (2013) and comparisonof numerical simulations with models of fuel-based specific impulse. . . . . . . . . 22
Triple-point and axial-flow solution states for representative fuel-air and fueloxygen mixtures. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33

iv

Shepherd and Kasahara - September 23, 2017

1

1

Introduction

The rotating detonation engine or RDE (Fig. 1a) is based on continuous flow through an annularchamber in which combustion takes place through a transversely rotating detonation. The fuel andoxidizer are injected at the entrance of the combustor and the high-pressure combustion products
behind the detonation wave expand and turn to create a primarily axial flow at the exit of the
combustion chamber. The reviews by Rankin et al. (2017) and Zhou et al. (2016) provide moredetails about RDEs and Kailasanath (2017) summarizes recent studies. Although variants ofthe combustor shown in Fig. 1a are being considered as components in air-breathing propulsion
systems, in the present study we examine the performance of the annular combustor in a rocket
motor configuration.1 For actual rocket motor applications, an obvious extension is to add an
exhaust nozzle; this is not the focus of the present study although the effect of an ideal nozzle isbriefly considered.
The purpose of the present article is to provide two simple analytical methodologies for computing the thrust of the configuration of Fig. 1a. The first methodology, termed the “pressure history
model”, is motivated by the success of ideas that were originally developed for pulse detonationengines (PDE) (Wintenberger et al., 2003; Endo et al., 2004) and the observations (Rankin et al.,2017; Fotia et al., 2017) that models developed to predict PDE performance appear to be quantitatively useful for RDEs. We shall demonstrate that although the principles of operation appear tobe quite different for RDE and PDE operation, it is possible to use very similar modeling ideas topredict RDE performance under certain conditions. The second methodology, termed the “axial
flow model”, is based on the approximation of mean axial flow and uses a conventional control
volume analysis that focuses on exit conditions to evaluate thrust. This model is motivated by the
results of numerical simulations that show the approach to axial flow with downstream distanceand the success of Kato et al. (2016) in using conventional rocket motor performance analyses
in correlating RDE test results. The axial flow speed and thrust are evaluated by approximating
the flow following the detonation as isentropic and considering a two-dimensional expansion thatconverts azimuthal motion into purely axial flow.
In order to make progress with analytical modeling, we will need to make a number of simplifying assumptions. First, we will assume that the fuel and oxidizer are supplied at a constant average
rate, are well mixed and the detonation wave moves at a constant speed UD . Under these conditions,
we can further assume that the flow within the engine will be periodic in time. Viewed from thelaboratory frame, as shown in Fig. 1a, the detonation wave appears to rotate with angular speed ω= UD /R in the azimuthal direction, which gives the device its name. The actual detonation speed,
UD in Fig. 1a, may differ from the Chapman-Jouguet (CJ) speed for a number of reasons such as
incompletely mixed fuel and oxidizer, the influence of the expanding flow within the reaction zone
immediately behind the wave and the upstream velocity of the unburned gas. In our models, wewill assume that the detonation propagates with the CJ speed relative to the incoming flow.
Second, we assume that there is a single detonation front within the chamber and the flow will
oscillate with frequency f = ω/2π. The flow properties will be periodic with temporal period T= 1/ f . This greatly simplifies the analysis because in a frame co-rotating with the detonation, theflow will be steady. This corresponds to transforming to a coordinate θ 0 = θ − ωt and eliminating
1Due to the large total pressure loss associated with the arrays of small slots or holes used to inject the fuel andoxidizer, there are a number of challenges in incorporating these combustion chambers into the core of a gas turbine orsimilar high-mass-flow engine.

2

Shepherd and Kasahara - September 23, 2017

(a)

(b)

Figure 1: a) A rotating detonation engine test fixture (Ishihara et al., 2015). b) Schematic of“unwrapped” flow field projected onto the θ - Z plane.
time as a variable. For example, the pressure can be written as P = P(r, θ, z, t) = P(r, θ − ωt, z) =P(r, θ 0, z). The restriction to a single wave turns out to be artificial and we will show that the resultsfor one wave can be easily extended to multiple waves.
Third, we will assume that radial variations in the flow can be neglected and the flow can beapproximated as two-dimensional in space, i.e., P(θ 0, z). Developing (“unwrapping”) the surfaceof the cylinder onto a plane and considering the flow to be uniform in the radial direction reduces
the geometry to the X-Z plane with the flow field shown in schematic of Fig. 1b. In the X-Z plane,
the flow is spatially periodic in the X direction, f (x + 2πR, z) = f (x, z) with period 2πR. This
representation is commonly employed in many numerical simulations of RDEs. Eude et al. (2011)
have compared pressure distributions computed from two and three-dimensional simulations anddemonstrated that integral (average) properties are insensitive to three-dimensional effects and twodimensional models are adequate for predicting properties such as thrust. The transformation to a
two-dimensional system that is co-rotating with the detonation wave is presented in Appendix A
By averaging over a period, the azimuthal and time dependence of the flow can be eliminatedfrom consideration to facilitate the computation of average thrust. For the purposes of developing
simple models, the variations in the radial direction are also averaged so that the flow in the
combustion chamber can be considered quasi-one dimensional with average properties dependent
only on the axial distance z. For example, the pressure can be averaged to obtain
1
P(z) =
T

∫ t+T

1

∫ Ro

t

2πRW

Ri

1
P(r, θ, z, t) 2πr dr dt =
2π

∫ 2π

1

∫ Ro

0

RW

Ri

P(r, θ 0, z) r dr dθ 0 ,(1)

where W = Ro − Ri is the radial width of the combustion chamber, assumed to be constant for the
purposes of analysis and R = (Ro + Ri )/2 is the mean radius of the combustion chamber annulus.

3

Shepherd and Kasahara - September 23, 2017

In the co-rotating reference frame, the time averaging is equivalent to averaging over the azimuthal
angle θ, as shown in the second version of the integral in (1).

2

Control Volume Analysis of Thrust

The thrust for an engine without an exit nozzle, i.e., just a straight combustion chamber, can beevaluated using the standard approach of evaluating the integral momentum balance over a control
volume surface A surrounding the engine,∫∫∫
ρuu · n dS +Pn dS =
τ · n dS + F .
(2)
A

A

A

Two possible control volumes are shown in Fig 2. These control volumes and the idealized
combustion chamber are axially symmetric. The force F is the reaction force needed to keep theengine stationary; it is equal in magnitude and opposite in direction to the thrust. For both control
volumes, we assume that the pressure on the exterior of the engine is everywhere equal to theatmospheric pressure Pa . These control volumes and the simplifying assumptions are identical tothose used in the PDE models of Wintenberger et al. (2003) and Endo et al. (2004).

(a)

(b)
Figure 2: Two approaches to formulating a control volume for a rotating detonation engine tocompute thrust base on a) chamber exit conditions. b) combustion chamber inlet conditions.
The control volume of Fig. 2a completely encloses the entire engine and crosses the exit planeof the combustor. This enables the computation of the thrust based on exit conditions. The balanceof forces in the axial (z) direction for this control volume, assuming static conditions and uniform

4

Shepherd and Kasahara - September 23, 2017pressure everywhere except at the exit, is:
∫
∫
2
Fz =
ρe uz,e dS +(Pe − Pa ) dS .
Ae

(3)

Ae

To evaluate the integrals, we apply the averaging method described in the introduction to computethe average magnitude of the thrust T = F z . The momentum of the fuel/oxidizer supply does notcontribute to the force balance because the supply is assumed to enter perpendicular to the axial
direction. In order to use this version of the thrust computation, the entire flow field within the
combustor must be analyzed in order to define the exit conditions. This is the approach taken by
many previous researchers Fievisohn and Yu (2017); Schwer and Kailasanath (2013, 2011b). If theflow at the exit of the combustor is uniform, then the thrust expression is the same as that used in
idealized rocket motor models
Û e + Ae (Pe − Pa ) .
(4)
T = MuThis expression is used to evaluate thrust for the mean axial flow model described in Section 4 andis also appropriate for combustion chambers with an ideal exit nozzle.
As pointed out by Schwer and Kailasanath (2011b), as long as the pressure ratio (plenum/exit)
is sufficiently large, the detonation propagation is unaffected by the exit conditions as the flow
throughout the combustion chamber remains supersonic. This implies that it is possible to only
consider the flow at the inlet of the combustion chamber and an alternative control volume construction analogous to that used by Wintenberger et al. (2003) and Endo et al. (2004) for PDEs todevelop an approximate analytical model for RDE thrust. This is the basis of the pressure history
model of RDE thrust developed in Section 3.
To develop this alternative method of thrust evaluation, consider the control volume of Fig. 2bwhich encloses the engine but does not cross the exit plane. Instead, the control volume is placedalong the walls of the combustion chamber (in the gas) and crosses the inlet plane of the combustor.
This is the same type of control volume that is used in the PDE models of Wintenberger et al.
(2003) and Endo et al. (2004). This enables the computation of the thrust based on the conditions
at the porous thrust surface with a correction for the drag on the walls of the chamber. The balanceof forces in the axial (z) direction for this control volume is
∫
∫
∫
2
Fz =
(Pc − Pc,1 ) dS +[ρc uz,c + (Pc,1 − Pa )] dS −
τw dS ,
(5)
Ac

|

Ac

{z
I

} |

Aw

{z
II

}

|

{z
III

}

where Ac is the area of the frontal surface of the porous injector at the upstream end of the combustor(labeled thrust wall in Fig. 1a) and Aw is the area of the annular walls of the combustor that are
subjected to the shear stress τw created by the axial flow. The pressure Pc,1 is the value of the pressurein the reactant layer just upstream of the detonation wave. In order to make quantitative estimatesof the thrust using this control volume analysis, we now consider appropriate approximations foreach of the terms in (5).
Term I represents the integrated effects of the detonation-generated pressure force on the injector
surface at the upstream end of the combustion chamber. We will assume that radial variation of thepressure is small so that the integral can be evaluated by multiplication of the area, Ac = π(Ro2 − Ri2 )= 2πRW, times the average value of the integrand,
FI = Ac (P c − Pc,1 ) .

(6)

5

Shepherd and Kasahara - September 23, 2017

This is exactly analogous to the result used to compute impulse I for a single-cycle of PDE operation
using the temporal average of pressure on the thrust surface rather than the spatial average. To
compute the integral, we must compute P c using the operating conditions and geometry of theengine.
Term I I represents the force associated with the momentum flux of the gas at the inlet to the
combustion chamber. Generally, the flow into the combustion chamber is low-speed and the term
ρu2 does not contribute substantially to the thrust although the pressure difference Pc,1 − Pa acting
on the face of the injector may be significantly different than ambient in some cases. The magnitudeof term I I will depend on the engineering details of the injection system. If we assume uniform
flow over the inlet surface, the approximate value of the integral will be
FI I ≈ Ac (Pc,1 − Pa ) + uc MÛ .

(7)

If the average density ρc ahead of the detonation wave is known instead of the velocity uc , then theestimate is≈ Ac (Pc,1 − Pa ) +

MÛ 2
.
ρc Ac

(8)

Term I I I represents friction due to the viscous forces on the side walls due to flow within the
combustion chamber. In addition to the assumptions made above, although not explicitly includedin the momentum equation, we assume that the combustion chamber is adiabatic, i.e., heat losses
to the chamber walls are neglected. The effect of heat loss will be to reduce the mean effective
pressure and thrust from the ideal values computed by the present analysis; just as in the PDE caseas discussed by Radulescu and Hanson (2005) and Kawane et al. (2011).
The total thrust will be given by the sum of the three terms
T = FI + FI I + FI I I .

(9)Using engineering estimates for each of these three terms we observe that the majority of thecontribution to the thrust will be term I, with term I I being a factor of 5-10 smaller depending on
the fuel and engine specifics and term I I I is at least one to two orders of magnitude smaller than thesum of terms I and I I. The present analysis focuses on term I, estimates term I I in order to obtainrealistic comparisons with computations and measurements of other researchers, and neglects term
I I I.
This analysis also neglects the contribution of the pressure on the base of the engine located at
the exit of the combustor. Interaction between the exit plume and the base may result in the averagepressure P b on the base being sufficiently different than the atmospheric pressure Pa that there will
be a fourth contribution(10)
FIV = Ab (P b − Pa )to the thrust. This approach could also be used to account for the presence of a converging or
converging-diverging nozzle, or a conical aero-spike attached to the inner portion of the base. In
order to calculate the average base pressure or account for a nozzle, a more detailed flow analysis
will be required.

Shepherd and Kasahara - September 23, 2017

3

6

Pressure History Model

The pressure history model is based on evaluating the thrust associated with the pressure history
on the thrust wall due to the propagation of the detonation and expansion flow of the detonation
products. In order to find the average pressure on the thrust wall, we need to model the gas
dynamics of the flow field created by the rotating detonation. A simplified model is to consider the
propagation of the detonation in a layer of height H over a solid surface, setting aside the effect ofporosity of the thrust surface and the inflow of reactants that occurs in an RDE. This is equivalent
to assuming that the flow ahead of the detonation is parallel to the thrust surface and the detonation
wave is normal to the surface.
This model is best analyzed in the X 0-Z 0 frame that is co-moving with the detonation as shownin Fig. 3. In this frame, the flow is steady and consists of a planar detonation in the reactive
layer which creates an oblique shock wave in the bounding gas and an expansion wave behindthe detonation. The idealized flow is assumed to be self-similar near the triple point where thewaves intersect, this region is shown as the detail in Fig. 3. The postshock state 5, and the state 3
following the expansion wave must be in mechanical balance which means that pressure and theflow deflections are continuous across the contact surface: P5 = P3 ; θ 5 = θ 3 . There must also bemechanical balance along the contact surface ahead of the triple point: P1 = P4 ; θ 1 = θ 4 . The
inflow from the plenum will result in θ 1 being nonzero and a tilted detonation wave as shown in
Fig. 1b, however this does not change the triple point balance considerations. More significantly,
the periodic nature of the flow will result in a nonzero azimuthal flow uθ upstream of the oblique
shock in the stationary (laboratory) reference frame so that u4 , u1 in the co-moving frame where
u4 = UC J + uθ and u1 = UC J . In this situation, the contact surface 1-4 will be a shear layer with avelocity difference of up to 1000 m/s (See Eude et al., 2011, Fig. 5).

Figure 3: Simplified model for detonation of a uniform layer of thickness H and detail of triple
point region (red box).
There are two key issues in applying this model to an RDE. First, a model of the injection
process is needed to determine the pressure just upstream of the detonation. Fievisohn and Yu

Shepherd and Kasahara - September 23, 20177(2017) discuss a simple control volume model and implementation as boundary condition in their
method of characteristics computations. Schwer and Kailasanath (2011a) compare a detailed
fluid dynamic model of discrete injectors with two approaches that model injection as a boundary
condition. Second, the state of the bounding gas is not known ahead of time but is created by
the expansion of the detonation products. In general, this expansion creates a non-uniform state
upstream of the shock wave in the bounding gas; this state must be determined in a self-consistentsolution considering the periodic flow associated with the annular geometry. Assuming that thepressure upstream of the detonation is given, solutions to the triple-point configuration can be foundby assuming isentropic expansion of the detonation products to the upstream pressure in order todefine state 4. The iterative solution of the 3-5 contact surface matching problem (P3 = P5 and θ 3= θ 5 ) is equivalent to graphically finding the intersection of the bounding layer shock polar with
the Prandtl-Meyer (PM) expansion solution between states 2 and 3 as shown in Fig. 4 for threecases. There are two possible solutions for states 3 and 5 but the attached oblique shock case thatis observed in simulations and experiments is the lower pressure solution located on the supersonicbranch of the oblique shock polar. These solutions were computed using realistic thermochemistry
and equilibrium properties using the Shock and Detonation Toolbox Browne et al. (2017). The
numerical values associated with states 1–5 that were determined by these computations are given
in Appendix B for three fuels, hydrogen (H2 ), ethylene (C2 H4 ) and propane (C3 H8 ) with an initialstate of P1 = 100 kPa and T1 = 300 K..
Note that all three ideal RDE cases shown in Fig. 4 admit solutions with attached oblique
shocks. However, in the most general case of arbitrary conditions in the bounding gas layer, asolution with an attached oblique shock may not exist. This was discussed by Fievisohn and Yu(2017) using ideal gas models of the polars and considering parametric variations of the bounding
atmosphere temperature. Our computations with realistic thermochemistry and parameters relevant
to RDE operation for a wide range of reactant mixtures give results similar to Fig. 4 and we findthat solutions are obtained for states 3 and 5 as long as state 4 corresponds to expanded detonation
products. Houim and Fievishon (2017) examined the case of a hot, stationary inert gas atmospherebounding a cold detonation layer and found that if the bounding layer is sufficiently hot, attachedoblique shock solutions are not observed; see Appendix C for further discussion.
The model of Fig. 3 has been examined previously in the context of fuel-air explosions by Sichel
and Foster (1979) who considered the impulse delivered by a fuel-air layer bounded by the groundat the bottom and the atmosphere at the top. Sichel and Foster used analytical expressions basedon a constant ratio of specific heats (but with different values for shocked air and the combustion
products) to determine the solutions for states 5 and 3 for MAPP2 gas-air and methane-air mixtures
over a range of equivalence ratios. The flow field following the detonation was determined using
the method of characteristics to compute the interaction of the expansion wave with the ground.
The method of Sichel and Foster was recently adapted by Fievisohn and Yu (2017) to the RDE
flowfield and extended to include inflow and periodic boundary conditions. Fievisohn and Yu
describe their method-of-characteristics solution technique in detail and showed that their solution
compared favorably with the finite-difference solution of Schwer and Kailasanath (2011b).
Sichel and Foster (1979) pointed out that if the detonation reaction zone is sufficiently thin,then for purposes of computing impulse from an infinitely long layer, the only relevant dimensionalparameters are the layer height H and the detonation overpressure ∆PC J . Further, it is reasonable
2MAPP gas is a mixture of Methyl Acetylene, Propadiene, Propane and other C3 and C4 hydrocarbons

8Shepherd and Kasahara - September 23, 2017

C2H4-air15

pressure P/P1

H2-air

10

5

0
0

20

40
60
flow deflection (deg)

80

100

35

C2H4-O230

pressure P/P1

25
20
15
10
5
0
0

20

40
60
flow deflection (deg)

80

100

Figure 4: Three examples of pressure-deflection diagrams for 3-5 interface mechanical matchingsolution.

9

Shepherd and Kasahara - September 23, 2017

to suppose that the main features of the flow are governed by inviscid gas dynamics and the conditions behind CJ detonations correspond to sonic flow in wave-fixed coordinates for all mixtures.
Considering these factors, they concluded that the pressure along the ground behind the detonation
front can then be written in terms of a non-dimensional function ψ and distance ξ = x/H as
P(x) − Pc,1 = (PC J − Pc,1 )ψ(x/H) = ∆PC J ψ(ξ; γe ) ,(11)where Pc,1 is the pressure in the layer of explosive gas upstream of the detonation. The function
ψ depends on the properties of the combustion products, indicated here as a dependence on theratio of specific heats γe . Using (11) for the distribution of pressure along the thrust surface inan RDE, we evaluate the average pressure (1) due to the detonation in terms of an integral of thenondimensional pressure
P c − Pc,1 = ∆PC J∫ 2πR/HH
2πR

ψ(ξ) dξ .

(12)0

The integral for a given value of 2πR/H and γe is a nondimensional function K, defined as
K=∫ 2πR/Hψ(ξ; γe ) dξ .

(13)

0

As discussed by Bykovskii et al. (2006), there is another essential scaling function of the height
H and the annulus width W. As is well-established from fundamental studies (Lee, 2008) on
detonation behavior, the magnitude of H and W relative to the detonation wave cellular structure
size (cell width or length) or equivalently reaction zone length, is a key consideration in determining
if detonations can be initiated and propagate continuously. If the layer height or annulus width
is too small (Bykovskii et al., 2006) relative to the detonation cellular structure width or length,
it will be impossible to initiate and sustain continuous detonation. These are very importantpractical considerations and will determine the range of fuel-oxidizer-diluent compositions as wellas injection conditions for successful operation of a given combustor fixture in the continuous
detonation mode. In the present analysis we are assuming that these conditions are met and theonly role of the layer height is as a gas dynamic scaling parameter.
Sichel and Foster (1979) found the function ψ to be relatively independent of the gas composition
so they proposed a universal relationship that could be used to compute the pressure and impulse
for a range of fuel-air detonations. The RDE computations by Schwer and Kailasanath (2013) for
injection face pressure also support this notion, with consistent pressure distributions for H2 anda variety of C2 and C3 hydrocarbons in air and distinct but also consistent distributions for thefuel-O2 cases; two examples are shown in Fig. 5; the lines are data digitized from Fig. 7 and 9 of
Schwer and Kailasanath (2013).
Motivated by the success of simple analytic fits to the expansion wave behind propagating
detonations (Beltman and Shepherd, 2002), we propose the following model for the pressure decay(along the injection surface) behind the detonation front,ψ(ξ) = exp(−αξ) ,(14)

where α is a constant that is determined from fitting the data. Note that this model does notdescribe the narrow pressure spike and peak at the detonation front that is shown in Fig. 5. This

10

Shepherd and Kasahara - September 23, 2017
50

80

45

70

35

60

30

50

P (atm)P (atm)

40

25
20
15

40
30
20

10

10

5

0

0
0

5

10

15
X (cm)

(a)

20

25

30

0

5

10

15
X (cm)

20

25

30

(b)

Figure 5: Pressure at injector face (lines) Schwer and Kailasanath (2013) compared to the exponential decay model (points) described in text. a) Stoichiometric C2 H4 -air and b) StoichiometricC2 H4 -O2 .
spike represents the high pressure behind the leading shock front, decaying rapidly to the CJ
pressure with the first cm behind (to the left on this graph) the leading shock front located at 14
cm. This spike contributes a negligible amount to the impulse and as in the previous model of PDE
performance, will be neglected. Integration of the normalized pressure signal yields a value of K =1.02. For comparison, analysis of the results of Sichel and Foster yields K =1.3, the higher value isconsistent with the larger confining effect3 associated with a cold air bounding layer as comparedto hot detonation products. Using the parameters given by Schwer and Kailasanath (2013) for aC2 H4 -air detonation with initial conditions of 150 kPa and 255 K, and a layer height of 3.2 cm, we
obtain from our model equations, (11) and (14), the points shown in Fig. 5a for a value of α = 1/K= 0.98. Although the model results do not exactly match the pressure profile near the wave front,
the results are both qualitatively and quantitatively reasonable.
The layer heights are substantially different for fuel-air and fuel-oxygen cases analyzed by
Schwer and Kailasanath (2013); in order to compare the present model with their results it isnecessary to use the mass flux and conditions ahead of the detonation to compute H using massconservation
MÛ
.
(15)
H=
ρc WUC J
Using this estimate, the fuel-oxygen cases of Schwer and Kailasanath (2013) were analyzed. Direct
integration of the pressure trace shown in Fig. 5b yields K = 1.54; the model pressure history
using α = 1/K = 0.65 yields the points shown in Fig. 5b. The differences between simulations andmodel pressure histories are larger than in the C2 H4 -air case due to the presence of a shock wave
within the expansion flow; apart from this obvious discrepancy, the model is a reasonable empiricalapproximation to the simulation results. The difference in α between fuel-air and fuel-oxygencases can be attributed to the differences in the combustion product properties due the the larger
energy content and substantially higher detonation product temperatures in fuel-oxygen mixtures
as compared to fuel-air mixtures. As a consequence, flow deflection angles θ 3 are larger for the
3The extent of this effect can judged from the difference in values of θ 3 ≈ 20◦ for Sichel and Foster compared to θ 3≈ 30◦ for the RDE cases (Appendix B) using fuel-air mixtures.

11

Shepherd and Kasahara - September 23, 2017

fuel-oxygen cases than for fuel-air cases but the pressures P3 are also higher. Systematic studies
using numerical simulations for a range of mixtures will be needed to further clarify the dependenceof K on mixture properties.
The value of K can be obtained by analytic evaluation which can be simplified in most cases
because the 2πR/H is between 10–20. This means that the upper limit to the integration in (13)can be taken to be +∞ without any loss in fidelity for the purpose of the present model so that
∫ ∞
1
(16)
K≈
exp(−αξ) dξ = .
α
0

3.1

Evaluating the thrust

From (12) and (13) , the average thrust due to the detonation is
FI = ∆PC J

Ac H
2πR

K = ∆PC J W HK ,(17)

and based on the analysis presented in the previous section, we use K = 1.02 for fuel-air mixturesand K = 1.54 for fuel-oxygen mixtures. The specific thrust due to the detonation is defined asÛ where MÛ is the total (fuel and oxidizer combined) mass flow rate. The total mass flow rate
FI / M,
can either be computed from the inflow velocity, MÛ = Ac ρc uc or from the rate of mass consumption
by the detonation, MÛ = HW ρc UC J . Using the second approach, the specific thrust is
∆PC J
FI
K.
=
ρc UC JMÛ

(18)

This expression can be further simplified by using the following approximate expression (Browne
et al., 2017) for ideal gas CJ pressure,ρc UC2 J.
γe + 1

(19)

FI
UC J
≈K
.
γe + 1MÛ(20)∆PC J ≈
Substituting this expression into (18), we obtain

Another typical figure of merit for RDEs and PDEs is the specific impulse per unit mass of fuelconsumed divided by the acceleration of earth’s gravity g . Defining the fuel mass fraction as Y f =Û the specific impulse per unit mass of fuel consumed in one revolution of the detonation
MÛ f uel / M,
isFI
Isp f ,I =,(21)
Y f g MÛwhich can be approximated as
Isp f ,I ≈K UC J.
Y f g γe + 1(22)Values of specific thrust and the fuel-based specific impulse due to detonation alone for a rangeof equivalence ratios and either air or pure oxygen as the oxidizer are shown in Fig. 6 for three12

Shepherd and Kasahara - September 23, 2017

fuel-based speciﬁc impulse Ispf (s)representative fuels. These computations are based on computed CJ properties with the state just
upstream of the detonation at 295 K and 100 kPa; however, we expect that these specific performancemeasures will be relatively independent of conditions upstream of the detonation. This is based on
the functional form of (20) and (22) and the known modest dependence of UC J and γe on initialconditions.
The variations in specific thrust with equivalence ratio shown in Fig. 6a are a consequenceof the dependence of the CJ speed on the energy content of the mixture and the sound speed of
mixture upstream of the detonation. This is discussed at length in Wintenberger et al. (2003) and aconsequence of the thermodynamic relations that determine the CJ state as discussed in Appendix
B of Browne et al. (2017). The trends with equivalence ratio shown in Fig. 6b are primarily due
to the inverse dependence on fuel mass fraction of the model (21) and are the same trends thatare predicted by the PDE model of Wintenberger et al. (2003). The fuel-based specific impulse ofC2 H4 and C3 H8 are indistinguishable throughout the equivalence ratio range because the energy
content per unit mass of these fuels as well as other alkanes is essentially identical.
5000
H2 -O2C2 H4 -O2
C3 H8 -O2
H2 -Air
C2 H4 -Air
C3 H8 -Air

4000
3000
2000
1000
0
0.5

1

1.5

2

2.5

3

equivalence ratio

(a)

(b)

Figure 6: a) Specific thrust T / MÛ and b) fuel-based specific impulse Isp f of RDEs as predictedusing the pressure history model.
So far, we have considered a single detonation wave but multiple waves are often observed(Kato et al., 2016; George et al., 2017). Consider a number N of waves within the chamber, allrotating in the same direction (counter-rotating waves introduce additional complications that wehave not considered). For a total mass flow rate MÛ through the combustor, each wave combusts
mass at the rateÛMÛ i = M/N.
(23)
The height Hi of each wave will be given by mass conservation
Hi =

MÛ i
= H/N ,
ρc WUC J(24)

where H is the height of a single wave at the same mass flow rate. The nondimensional integral of
the pressure has an upper bound of(2π/N)R 2πR
=
(25)
H/N
H

13

Shepherd and Kasahara - September 23, 2017

which is independent of N because the extent of the “tail” of the pressure wave behind detonation
wave will be 1/N of the total circumference of the combustion chamber. This implies that the ideal
value of the nondimensional integral is independent of the number of waves∫ 2πR/H
Ki = K =exp(−αξ) dξ(26)
0or in dimensional terms
1
K=∆PC J H∫ 2πR

Pc (x) − Pc,1 dx

(27)

0

and if there are N waves each of height Hi , the effective height of the equivalent single wave is
given by
H = N Hi .

(28)

The net thrust associated with the detonation pressure on the combustor inlet from all N waves willbe the sum of the thrust of each individual wave
FI =

N
Õ

∆PC J W Hi Ki

(29)

H
K
N

(30)

i=1

=

N
Õ
i=1

∆PC J WN
H Õ
1= ∆PC J W K
N i=1

(31)

= ∆PC J W HK ,(32)

which is identical to the result obtained for a single wave. This is a consequence of the scaling of
all spatial dimensions with the height of the detonation wave. The insensitivity of the thrust to thenumber of detonation waves is confirmed by numerical simulations (Yi et al., 2011; Tsuboi et al.,2017).

4

Axial Flow Model

Although the detonation wave rotates, this does not result in a net (average) azimuthal rotational
flow at the exit of the combustion chamber if, as is usually the case, the flow entering the chamber
does not have a net azimuthal rotation. This does not mean that the velocity in the azimuthal
direction vanishes everywhere within the flow field, but at any axial location z, there are positive
and negative values of azimuthal velocity which average out so that the total angular momentum is
zero.
To show this, consider the angular momentum about the z-axis of the coordinate system shownin Fig. 7. The angular momentum of an element of volume dV isr × ρu dV .

(33)

Shepherd and Kasahara - September 23, 201714

The conservation of angular momentum for a stationary control volume Ω with control surface ∂Ωis∫∫∂(r × ρu) n · u dS ,r × ρu dV +
(34)
M=
∂t Ω
∂Ω
where M is the total torque exerted on the control volume. For an annular control volume thatencompasses the RDE shown in Fig. 7, the imposed net torque is zero as the engine is not rotating,there are no rotating components within the flow path and we are neglecting the effect of friction
on the combustion chamber walls. Because the flow is periodic, averaging over one period of thecycle will result in the first term on the right-hand side averaging to zero, leaving∫ ∫1 T(r × ρu) n · u dS dt .
0=(35)
T 0 ∂Ω
For the geometry of Fig. 7, the integral can be written in terms of velocity components and simplified
by transforming to the wave-fixed (steady flow) frame of reference∫0=Rρuθ uz dS ,(36)
A

where R is the perpendicular distance from the point in the flow to the z axis and A is the area of
the annulus perpendicular to the z axis, the direction of the mean axial flow. If we further supposethat the flow is uniform in the radial direction, this reduces to
∫ 2π
1
(37)
ρ(θ 0, z)uθ (θ 0, z)uz (θ 0, z) dθ 0 = ρuθ uz = 0 .
2π 0
Axial and azimuthal velocities are correlated as shown by Fig. 5 of Schwer and Kailasanath (2011b)with fluctuations of uz about a mean positive value and both positive and negative excursions of uθ
that support the notion that (37) will be satisfied. The fluctuations in velocity and density are theconsequence of the flow field, particularly the oblique shock waves, which simulations show decay
in strength with distance from the detonation layer. This suggests that sufficiently far from thedetonation layer the flow will primarily be axial and relatively uniform in the azimuthal direction.
Using the nomenclature of Appendix A, uθ = v and uz = w, the vanishing of the mean angularmomentum suggests that the mean azimuthal flow speed will be small in comparison to the mean
axial flow speedvw(38)
far downstream from the detonation layer. The validity of this assumption will depend on the axial
length of the combustor and represents an idealization, one that is important in formulating theaxial flow model presented in this section.

4.1

Modeling the expansion process

The net result of the two-dimensional expansion process within the annular combustion chamber is
to convert the flow from the almost entirely azimuthal flow immediately following the detonation
to primarily axial flow at the chamber exit. The expansion simultaneously decelerates the flow
induced by the detonation and turns the azimuthal flow of the detonation products so that it becomesShepherd and Kasahara - September 23, 201715

Figure 7: Geometry of rotating detonation combustion chamber and associated axisymmetriccoordinate system.
predominantly axial. Interaction of the flow with the oblique shock wave results in further reductionof the azimuthal flow components. We anticipate that if the combustion chamber is sufficiently
long, then the azimuthal flow can be neglected, and the performance upper bound can be basedon the same notions as traditional rocket motor performance based on quasi-one-dimensional flow.
The extent to which the flow approaches a purely axial state will depend on the length of the
combustor compared to the height of the detonation layer. Examples of the diminished azimuthal
flow at the exit of the combustor can be found in the simulations of Yi et al. (2011); Schwer andKailasanath (2011b); Davidenko et al. (2007).
The basis of ideal steady-flow analyses in propulsion systems with shock and detonation waves
is the conservation of total enthalpy (Wintenberger and Shepherd, 2006) along streamlines. In
the wave-fixed (steady flow) frame in the unwrapped 2-D coordinates, this is expressed as (seeAppendix A)
1
ht = h + (u02 + w02 )2(39)which in the laboratory (unsteady flow) frame for a detonation rotating with azimuthal speed UDtransforms to1
h1 = h + (v 2 + w 2 ) − vUD .
2

(40)

16

Shepherd and Kasahara - September 23, 2017

Considering a fixed axial location, the velocity and enthalpy will be a function of time or equivalently
azimuthal location. Using the averaging technique of (1), the average of (40) can be written
w 2 = 2(h1 − h) − v 2 + 2vUD .

(41)

As discussed in the previous section, the conservation of angular momentum (37) indicatesthat if the axial flow is sufficiently uniform, then the average azimuthal speed at the end of the
combustor will be negligible, v ≈ 0. However, simulations demonstrate that fluctuations in v are notnecessarily negligible, i.e., v 2 > 0, at the end of annular combustors with dimensions commonlyused in testing. The effect of fluctuations in v on the r.m.s. value of axial flow speed w at the exit
can be estimated from (41)q"qqw 2 = 2(h1 − h) − v 2 = 2(h1 − h) 1 −v2
2(h1 − h)

# 1/2
,

(42)

by assuming
v = 0. Numerical evaluation of (42) reveals that if the r.m.s. average of the azimuthalqq2
speed, v , is 30% of the ideal axial speed, 2(h1 − h), this will result in an approximately 5%
fluctuation (r.m.s.) in the axial velocity contribution to the average specific thrust. Given this modest
effect on the thrust we will neglect the effect of azimuthal and axial speed fluctuations, taking theeffective axial speed at the exit to be the ideal value computed from the mean thermodynamic stateat the engine exitqw≈2(h1 − h) .

(43)

The enthalpy can be considered a function of composition, pressure P and entropy s. In an
expansion wave following a detonation, the flow time scales are sufficiently long compared to thechemical time scales (see Wintenberger et al. 2004) that the flow is to a reasonable approximationin chemical equilibrium and so that h = h(P, s) only. In the ideal version of the flow that we areconsidering, entropy changes are due only to the flow crossing shock waves. However, the shockwaves are relatively weak, resulting in a modest loss in total pressure, for example about 7% in
the simulations by Yi et al. (2011). As a consequence we can consider the axial velocity to be afunction of pressure only and the entropy is taken to be constant and equal to the value of the CJ
state; for this ideal model we assume that detonation propagates at the CJ speed, UD = UC J . With
these assumptions we obtain the following simple expression for the axial flow velocity, whichappears identical in form to the expression used in ideal quasi-one-dimensional flow models:p(44)
w = 2(h1 − h(P, s2 )) ,with the important distinction that h(P, s2 ) is the enthalpy on the isentrope associated with thetwo-dimensional flow created by the products of the rotating detonation. The result (44) and theunderlying assumptions have been used by previous researchers, particularly by Davidenko et al.
(2011), who carried out both multi-dimensional numerical simulations and modeling of an ideal
RDE based on axial exit flow. Davidenko et al. (2011) carried out an extensive parametric study for
H2 -O2 mixtures comparing an ideal conventional rocket engine with ideal RDE performance basedon the axial flow approximation. Using axial flow models for RDE analysis has a long history andsome elements of the present model can be found in the pioneering studies by the University of17

Shepherd and Kasahara - September 23, 2017

Michigan group (Shen and Adamson, 1972; Adamson and Olsson, 1967; Nicholls et al., 1966) on
the RDE configurations of the type shown in Fig. 1.
The enthalpy h(P, s2 ) of the combustion products downstream of the detonation as a functionof pressure is shown in Fig. 8. At the CJ state (2), the conservation of total enthalpy across thedetonation wave implies h2 = h1 + 21 (UC2 J − w22 ) and there are no solutions (the dashed portion of
the curve) to (44) for PC J ≥ P ≥ Pm , where the limiting pressure Pm is defined by h(Pm, s2 ) = h1 .
In the case shown in Fig. 8a, the limiting pressure ratio Pm /P1 = 4.18 (P1 = 100 kPa) and axial flow
solutions for Pm > P > P∗ are subsonic and for P∗ > P are supersonic. The limiting pressure Pm
corresponds to the effective one-dimensional (axial flow) stagnation pressure since the axial flow
speed is zero at this point. Computations for a wide range of compositions indicates that Pm ≈0.22-0.25 P2 .
For given values of the exit pressure P and ambient pressure Pa , the specific thrust can be
estimated by (45) as a function of exit pressure assuming uniform axial flow and the control volumeof Fig. 2a.
P − Pa
T
=w+
.
(45)
ρw
MÛ
The density ρ and axial speed w are computed assuming isentropic expansion from the CJ stateand (44) so that the right-hand side is only a function of pressure. An example of specific thrust for
Pa = 100 kPa is shown in Fig. 8b as a function of exit pressure. Excluding the sharp changes as P
→ Pm or 0, the results are insensitive to pressure over a wide range of values. This suggests that asingle value of specific thrust can be estimated from the axial model without specifying the precise
value of the exit pressure. The model is valid not only for constant-area combustors but also for
combustors with variable area; for example, those with attached ideal nozzles which can be treatedas quasi-one-dimensional flow.
3000

specific thrust (m/s)

2500
2000
1500
sonic axial flow state

1000
500
0
0

(a)

1

2
3
pressure P/P1

4

5

(b)

Figure 8: Axial flow RDE model for C2 H4 -air: a) Pressure-enthalpy states; b) Specific thrust.
Using the value of specific impulse predicted at the sonic point, the results for three fuels areshown in Fig. 9. The results are remarkably similar to those obtained using the detonation pressure
history method shown in Fig. 6. The quantitative values for specific thrust and impulse predicted18

Shepherd and Kasahara - September 23, 2017by the axial flow model can be up to 30% higher than for the pressure history model for fuel-air
cases but there is less than 10% difference for the fuel-oxygen cases. The values of the specificthrust computed with the axial flow model depend not only on the composition of the reactants
but also the thermodynamic conditions in the reactants just ahead of the detonation. Parametriccomputations for the case of C2 H4 -Air mixtures are shown in Fig. 10.
fuel-based speciﬁc impulse Ispf (s)speciﬁc thrust (m/s)

2500
2000
1500
H2 -Air
C2 H4 -Air
C3 H8 -Air
H2 -O2C2 H4 -O2
C3 H8 -O2

1000
500
0
0.5

1

1.5

2

2.5

3

6000
H2 -Air
C2 H4 -Air
C3 H8 -Air
H2 -O2C2 H4 -O2
C3 H8 -O2

5000
4000
3000
2000
1000
0
0.5

1

1.5

2

equivalence ratio

equivalence ratio

(a)

(b)

2.5

3

1600
1500
1400
1300
1200
1100
1000
900
800
700
600

250 K
300 K
350 K

0.5

1.0

1.5
2.0
equivalence ratio

(a)

2.5

3.0

specific trhust (m/s)

specific thrust (m/s)Û b) fuel-based specific impulse
Figure 9: Sonic axial flow RDE model: a) Specific thrust T / M;
Isp f .
1600
1500
1400
1300
1200
1100
1000
900
800
700
600

400 kPa
200 kPa
100 kPa
50 kPa

0.5

1.0

1.5
2.0
equivalence ratio

2.5

3.0

(b)

Figure 10: Sonic axial flow RDE model: a) Effect of initial temperature T1 for P1 = 100 kPa; b)
Effect of initial pressure P1 for T1 = 300 K.

4.2

Analytical Expressions

The one-γ analytical model of the CJ state can be used with perfect gas thermodynamic expressions
to estimate the specific thrust. The starting point is the expression (Browne et al., 2017) for CJ

19Shepherd and Kasahara - September 23, 2017
Mach number in terms of energy release per unit masss12q(γ 2 − 1)=.
MC J −
MC J
γRT1

(46)

The CJ temperature and pressure are
TC J = T1(γMC2 J + 1)2(47)

MC2 J (γ + 1)2

and
PC J = P1

γMC2 J + 1.
(γ + 1)(48)

The enthalpy of the reactants is
h1 = CPT1 + q(49)

h = CPT .

(50)

and the enthalpy of the products is

The product temperature can be computed from the change in pressure and entropy along thestreamline emanating from state 2 (CJ)s − s2
T = T2 exp
Cp



P
P2

 R/Cp.

(51)

The change in entropy can be estimated from the change in total pressure along the streamline
assuming the total temperature remains constantPts − s2 = −R ln
Pt,2



.

(52)

These expressions can be used to derive an explicit expression for specific thrust for the particularcase of exit pressure equal to the atmospheric pressure (P = Pa ). As shown in Fig. 9, the specificthrust is relatively insensitive to the exit pressure so we choose the pressure-matched condition to
make the computation of specific thrust particularly simple. Assuming4 isentropic expansion along
the streamlines, s = s2 , we obtain"  (γ−1)/γ   (γ−1)/γ # 1/2
p
T
q
Pa
P1
T2
= w = 2CpT1 1 +
−
.
(53)
CPT1
P1
P2
T1
MÛ Pa
4Evaluating the entropy change for a total pressure loss of 7% results in correction to (51) of less than 1% to theisentropic temperature.

Shepherd and Kasahara - September 23, 201720

Substituting the analytical expressions for T2 = TC J and P2 = PC J , we obtain an explicit expression! −(γ+1)/γ  1/2
 2   (γ−1)/γ

2 
1
1
1
γ+1
Pa

1
+
M
−
−
.
C
J


2
2 +1
γ
−
1
2(γ
+
1)
M
P
M
γM
C
J
1


a
CJ
CJ


(54)
For the range 5 < MC J < 8 and 1.1 < γ < 1.15 which are of interest for RDEs, numerical evaluation
gives reasonable agreement (within 10%) to the value computed using detailed thermochemistry
for P1 = Pa . Examining the range of values for MC J for fixed γ, the specific thrust is observed tobe linearly proportional to CJ velocity, just as in the RDE analysis and consistent with the pressure
history model. We conclude that the scaling relationship
T
= a1
MÛ P

s

T
∝ UC JMÛ

(55)

is a feature of both PDE (64) and RDE (20) models; this result could be anticipated from dimensional
analysis although there are a number of ways in which a characteristic velocity could be defined.
For example, (54) can be interpreted as an example of a relationship of the formT
= a1 f (MC J, γ, Pa /P1 )
MÛ

(56)

as the consequence of dimensional analysis and without further examination of the functional formof f , it is not obvious that (55) holds.
Expression (54) is useful in discussing the dependence of specific thrust on temperature andpressure that is shown in Fig. 10. The modest increase in specific thrust with a decrease in initial
temperature (Fig. 10a) is due to the dependence of MC J on T1 as given by (46). The more dramaticincrease in specific thrust with increasing P1 for a given Pa is due to the decrease in enthalpy
with decreasing pressure; the CJ temperature (T2 /T1 ) and pressure ratios (P2 /P1 ) are essentiallyindependent of initial pressure so that a decrease in enthalpy and the corresponding increase in
specific thrust is entirely due to the greater isentropic expansion resulting from smaller pressure
ratios, Pa /P2 = Pa /P1 · P1 /P2 .

4.3

Axial Flow Model and Nozzles

Sufficiently far downstream from the detonation where the axial flow approximation is valid,
the usual quasi-one dimensional flow methodology employed for ideal rocket motor performance
estimates (Sutton, 1992, Chap. 3) can be used to examine the effect of a nozzle on the performanceof a RDE. Relationship (44) is identical to the conservation of energy relation for a nonreacting
quasi-one dimensional flow if we recognize that the effective total temperature Tt is defined by
h1 = CpT1 + q = CpTt and h = CPT. This is equivalent to defining the total temperature as theconstant pressure combustion temperature as suggested by Kato et al. (2016). The constant pressure
combustion temperature has a modest dependence on pressure so that a single value of Tt can beused for a given mixture composition and a wide range of combustor pressures. The conservationof energy (44) applies along a streamline in the flow and is equally valid within the constant area
combustor and an attached nozzle. As a consequence, the results of the axial flow model and thespecific thrust estimates (45) apply to configurations with nozzles downstream of the combustor

21

Shepherd and Kasahara - September 23, 2017and as examined below, the usual ideal rocket motor analysis can be applied using the appropriate
stagnation conditions.
To relate the axial flow model to the usual perfect gas methodology, the axial speed relationpw = 2(h1 − h(P, s2 ))(57)
can be approximated using the perfect gas model asqw ≈ 2(CPTt − CpT) =

s


 1/2
T
2γRTt
1−
γ−1
Tt

(58)

for T ≤ Tt (P ≤ Pt ). The standard perfect gas isentropic relationships can be used to define
temperature and density in terms of pressure. The density, axial speed and flow area A are related
through the conservation of mass flow.
MÛ = ρw A .

(59)

These relationships can be used to evaluate the specific thrust for a given nozzle exit area A using(4). For a given mass flow rate, an effective sonic throat area A∗ can be defined based on the sonicpoint in the expansion. At the sonic point w∗ = a∗, the mass flow will be
ρ a ∗
A .
MÛ ∗ = ρ∗ a∗ A∗ = ρt at
ρt at
∗

∗

(60)

Apply the perfect gas thermodynamic relationships, we obtain the standard relationship for mass
flow through a choked nozzle in terms of the stagnation conditions"MÛ ∗ = A∗ Pt

 γ+1/γ−1 # 1/2γ
2.
RTt γ + 1

(61)

Consider the situation in which the mass flow rate through the combustor is determined by the
injector system and a nozzle is attached to the end of the combustor. For a given total temperature
Tt , the chamber stagnation pressure Pt will be determined by the constancy of mass flow for steady
operation and the value of the nozzle throat area (Kato et al., 2016).
"
 γ+1/γ−1 # 1/2
MÛ RTt γ + 1
Pt = ∗
.
A
γ
2

(62)

As the nozzle area decreases, the stagnation pressure increases as observed in the tests of Katoet al. (2016). The axial flow model computations show that Pt ∝ P2 = C · P2 /P1 · P1 where C ≈0.24 so that increases in Pt will result in proportional increases in P1 because P2 /P1 is nominallya constant for a given mixture. Likewise, for fixed values of A∗ , the pressure Pt or P1 will increase
proportional to increases in mass flow rate as observed in the tests of Kato et al. (2016).

Shepherd and Kasahara - September 23, 2017522

Comparisons with simulations

The simulation parameters for four cases are given in Table 1 together with a comparison ofsimulation results to model estimates for the fuel-based specific impulse. The specific impulse
values estimated by the present models are within 5-10% of the values computed by Schwer andKailasanath (2013); the fuel-air estimates are slightly higher than the Schwer and Kailasanath(2013) values; the fuel-oxygen values are slightly lower except for the C2 H4 -O2 axial flow model.
Term I I amounts to 20% of the thrust for the H2 -air case and 15% for the C3 H8 -O2 case. The
pressure history model results shown in Table 1 were made using (7) and the value of uc = 300
m/s given by Schwer and Kailasanath (2013). The axial flow model results used the same pressure
and temperature (1.5 atm and 255 K) ahead of the detonation as the pressure history model. The
PDE model results are consistently lower (by up to 11%) than the axial sonic flow model results
but comparable to both the pressure history model predictions and the simulation values.
Table 1: Parameters for the cases reported in Schwer and Kailasanath (2013) and comparison ofnumerical simulations with models of fuel-based specific impulse.
case

H2 -air C2 H4 -air C2 H4 -O2 C3 H8 -O2
Model Input and CJ ParametersÛM (kg/s)
0.950
1.12
1.91
0.884
ρc (kg/m3 )
1.48
2.04
2.19
2.41
H (mm)
32.4
29.9
15.1
15.4
α
0.98
0.98
0.65
0.65
UC J (m/s)
1982
1836
2402
2383
PC J (MPa)
2.75
3.24
5.97
6.46
γe
1.169
1.165
1.142
1.137
Yf
0.02852 0.06375
0.22600.2160Fuel-Based Specific Impulse
Unsteady 2D simulations, Schwer and Kailasanath (2013)
Isp f S&K (s)
4860
1990
700
1070
Pressure History Model
Isp f (I + I I) (s)
4706
1975
704
1016
Axial Flow Model
Isp f (s)
5383
2280
911
952
PDE Model, Wintenberger et al. (2003)
Isp f (s)
4761
1985
750
781

6

Comparison to Test Results

Kato et al. (2014, 2016) performed experiments using a C2 H4 -O2 fueled RDE to directly measure
the average thrust with a load cell and instrumented the plenum and combustion chamber withpressure transducers. The annular RDE chamber had an inner radius of 30.25 mm, an outer radiusShepherd and Kasahara - September 23, 201723

of 33.45 mm and was approximately 100 mm long. The fuel and oxygen were injected from separate
plenums using impingement to promote mixing with the fuel injection being axial and the oxygen
injection oriented radially outward. The plenum pressures were approximately 5 to 10 times theaverage combustion chamber pressure. The flow through the injectors was choked and mass flow
rates were measured by using an electronic scale to determine the change in mass of the supply
tanks. The thrust was measured over an operational period of 0.2 to 0.8 s depending on the mass
flow rate and equivalence ratio.
Comparisons between measured specific thrust (assumed to be the load cell output) and model
predictions for selected tests from Kato et al. (2016) are shown in Fig. 11 using both the pressure
history and axial sonic flow models. Analyses of these tests require making assumptions about
the state of the reactants just upstream of the detonation front as there are no time-resolved
measurements available for the pressure or temperature. The measured values of mass flow rateand static pressure in the downstream portion of the combustor were used to compute the totalpressure at that point assuming axial flow as discussed in Kato et al. (2016). The pressure justahead of the detonation is inferred by using the observation from the axial flow model results thatthe total pressure is approximately 0.24 times the Chapman-Jouguet pressure and performing an
inverse CJ computation to determine P1 corresponding to PC J = Pt /0.24. The static temperature
was assumed to be equal to the value in the plenum. Some tests in Kato et al. (2016) used aconvergent or converging-diverging nozzle so that strictly speaking, the pressure history model is
not applicable and should only be compared to those tests with just the constant area combustion
channel.
Experimental values of specific thrust are 20-60% less than the pressure history estimates asshown in Fig. 11a, which may indicate a highly non-ideal detonation process with significant losses,
resulting in speeds much less than the CJ value. The influence of mixing on detonation propagation
is discussed by Fujii et al. (2017) and Gaillard et al. (2015) who describe the challenges obtaining
uniform mixtures with impinging jets. The filled data in Fig. 11a are for combustors with a constantarea and the open data points are for cases with nozzles. The data in Fig. 11b are for both constantas well as variable area cases. Representative ranges of values for the pressure P1 inferred fromthe experimental conditions are 35 to 40 kPa for the low mass flow cases and 60 to 70 kPa forthe high mass flow model. The test data trends in Fig. 11b are consistent with sonic axial flow
model for the high-mass-flow cases but show significant deviations for the fuel-rich, low-mass-flow
cases. Just as in the pressure history model, the test data are lower than the axial sonic flow model
prediction values in all cases. A key source of uncertainty that may be responsible for the differences
between model predictions and measurements is the lack of experimental data for the pressure andtemperature of the reactant gas layer just upstream of the detonation. Other possible causes forthe differences between model predictions and test data include nonideal mixing and detonation
processes as well as deviations from ideal axial flow at the exit of the test fixture. Resolving theseissues will require further tests and analyses.

7Comparison of RDE and PDE analysis

The RDE performance model can be compared with the PDE model of Wintenberger et al. (2003)
if we recognize that PDE impulse per cycle is the analog of RDE thrust and PDE mass per cycleis the analog of RDE mass flow rate. PDE specific impulse I/M is therefore the analog of RDE

24Shepherd and Kasahara - September 23, 2017
2000

2500

60 kPa

1600

2000
1500
1000

Type I - no nozzle
Type I - nozzle
Type II - no nozzle
Type II - nozzle

500

specific thrust (m/s)

pressure history model (m/s)

1800
1400

40 kPa1200

30 kPa1000
20 kPa800
600

Kato et al 2016 Type I - low mdot

400

Kato et al 2016 - Type I, high mdot

200

Kato et al 2016 - Type II

0

0
0

500
1000
1500
experimental specific thrust (m/s)

2000

0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 2.2 2.4 2.6
equivalence ratio

(a)

(b)

Figure 11: a) Measured and model estimates of specific thrust T / MÛ for the C2 H4 -O2 RDE tests ofKato et al. (2016). (a) Pressure history model and (b) sonic axial flow model; lines are predictions,points are test data.
specific thrust T / MÛ and PDE fuel-based specific impulse Isp f = I/(gY f M) is the analog of RDE
Û Specific thrust and specific impulse not only havefuel-based specific impulse Isp f = T /(gY f M).
the same units for both PDEs and RDEs but the analytical models have exactly the same dependence
on the detonation thermodynamic parameters. From equation (27) of Wintenberger et al. (2003),
the analytical model for impulse per PDE cycle is
IPDE = KPDE V

∆P3
.
UC J

(63)

To compare this to the RDE impulse, the plateau pressure difference ∆P3 can be expressed interms of the CJ pressure, ∆P3 ≈ 0.35∆PC J , and we use the nominal value of KPDE = 4.3 given
by Wintenberger et al. (2003). The mass of gas that participates in each cycle of the detonation
process is M = ρcV. Substituting these relations into (63), we obtain the following expression for
the PDE analog of specific thrust due to the detonation
 
∆PC J
I
≈ 1.5
,
(64)
M PDE
ρc UC J
which is identical in form to the RDE pressure-history specific thrust expression (18) but with aconstant K = 1.5 in place of KRDE = 1.02–1.54 (fuel-air vs. fuel-oxygen). This is the reason forthe surprisingly good qualitative correlation (Fig. 12) of measured RDE performance with the PDEmodel, as originally observed by Rankin et al. (2017). As shown in Fig. 12, the pressure history
and sonic axial flow RDE models bound the data for both RDE and PDE tests. The systematic
discrepancy between the RDE pressure history model and the RDE data can be explained by notingthat only term I of the model could be estimated from the information available in Rankin et al.
(2017). The sonic axial flow model appears to be about 10% higher than the experimental data for
H2 -air cases and less than 5% higher for the C2 H4 -air cases.

25

Shepherd and Kasahara - September 23, 2017
3000

H2-Air PDE (Schauer et al)
H2-Air RDE (Rankin et al)

6000

fuel-based specific impulse (s)fuel-based specific impulse (s)7000

H2-Air PDE Model (Wintenberger et al)
RDE Model (Term I)

5000

RDE Model (axial sonic flow)

4000
3000
2000
1000
0

C3H8-air PDE (Schauer et al)
C2H4-air RDE (Rankin et al)

2500

C2H4-Air PDE Model (Wintenberger et al)
RDE Model (Term I)
RDE Model (axial sonic flow)

2000
1500
1000
500
0

0.0

0.5

1.0
1.5
2.0
equivalence ratio

(a)

2.5

3.0

0.0

0.5

1.0
1.5
2.0
equivalence ratio

2.5

3.0

(b)

Figure 12: Comparison of predicted RDE and PDE performance with measured performance using
data from Fig. 10 of Rankin et al. (2017). a) H2 -air. b) C2 H4 -air.

8

Conclusion

The analyses given in this paper provide simple mechanistic explanations and relationships topredict RDE performance using either elementary thermochemical computations or approximate
analytical models based on perfect gas relationships. The analyses also demonstrate why RDE and
PDE models give very similar trends and identical qualitative dependence on detonation properties.
The pressure history model indicates how quantitative differences in performance occur due to thedifferences in the pressure-time history at the thrust surface and the possible importance of thethrust associated with the filling process in the RDE case. The axial flow model provides a way
to account for both pressure and mass flow contributions to thrust without explicitly evaluating thecontributions to the integrals required for the pressure history model. Both the pressure history andaxial flow models require knowledge of the conditions just upstream of the detonation; reinforcing
the importance of the fuel-oxidizer injection and mixing processes to the performance of RDEs. The
axial flow model provides a simple upper bound on the effective stagnation conditions which can beused in the traditional quasi-one dimensional rocket motor methodology to estimate performanceof RDEs used in a rocket motor configuration.

Acknowledgment
The authors thank the Japan Society for the Promotion of Science for an Invitation Fellowship for
Research that supported a visit to Nagoya University in the fall of 2015.

Shepherd and Kasahara - September 23, 201726

References
T. C. Adamson, Jr. and G. R. Olsson. Performance Analysis of a Rotating Detonation Wave Rocket
Engine. Astronautica Acta, 13:405–415, 1967.
W. M. Beltman and J. E. Shepherd. Linear Elastic Response of Tubes to Internal Detonation
Loading. Journal of Sound and Vibration, 252(4):617–655, 2002. doi: 10.1006/jsvi.2001.4039.
S. Browne, J. Zeigler, and J. E. Shepherd. Numerical Solution Methods for Shock and Detonation Jump Conditions. GALCIT FM2006-R3, California Institute of Technology, Pasadena,
CA, January 2017. URL http://shepherd.caltech.edu/EDL/publications/reprints/
ShockDetonation.pdf.
F. A. Bykovskii, S. A. Zhdan, and E. G. Vedernikov. Continuous Spin Detonations. Journal of
Propulsion and Power, 22(6):1204–1216, 2006. doi: 10.2514/1.17656.
D. M. Davidenko, I. Gokalp, and A. N. Kudryavstev. Numerical Simulation of the Continuous
Rotating Hydrogen-Oxygen Detonation with a Detailed Chemical Mechanism. In West-East
High Speed Flow Field Conference, Moscow, Russia, November 2007.
D. M. Davidenko, Y. Eude, I. Golkap, and F. Falempin. Theoretical and Numerical Studies on
Continuous Detonation Wave Engines. In 17th AIAA International Space Planes and Hypersonic
Systems and Technologies Conference, San Francisco, California, April 2011. doi: 10.2514/6.
2011-2334.
T. Endo, J. Kasahara, A. Matsuo, K. Inaba, S. Sato, and T. Fujiwara. Pressure History at the Thrust
Wall of a Simplified Pulse Detonation Engine. AIAA Journal, 42(9):1921–1930, 2004. doi:
10.2514/1.976.
Y. Eude, D. M. Davidenko, I. Gokalp, and F. Falempin. Use of the Adaptive Mesh Refinement for
3d Simulations of a CDWRE (Continuous Detonation Wave Rocket Engine). In 17th AIAA International Space Planes and Hypersonic Systems and Technologies Conference, San Francisco,California, April 2011. doi: 10.2514/6.2011-2236.
R. T. Fievisohn and K. H. Yu. Steady-State Analysis of Rotating Detonation Flowfields with
the Method of Characteristics. Journal of Propulsion and Power, 33(1):89–99, 2017. doi:
10.2514/1.B36103.
M. L. Fotia, J. L. Hoke, and F. R. Schauer. Experimental Performance Scaling of Rotating
Detonation Engines Operated on Gaseous Fuels. Journal of Propulsion and Power, 33(5):1187–1196, 2017. doi: 10.2514/1.B36213.
J. Fujii, Y. Kumazawa, A. Matsuo, S. Nakagami, K. Matsuoka, and J. Kasahara. NumericalInvestigation on Detonation Velocity in Rotating Detonation Engine Chamber. Proceedings of
the Combustion Institute, 36:2665–2672, 2017. doi: 10.1016/j.proci.2016.06.155.
T. Gaillard, D. Davidenko, and F. Dupoirieux. Numerical Optimisation in Non Reacting Conditionsof the Injector Geometry for a Continuous Detonation Wave Rocket Engine. Acta Astronautica,
111:334–344, 2015. doi: 10.1016/j.actaastro.2015.02.006.

Shepherd and Kasahara - September 23, 201727

A. St. George, R. Driscoll, V. Anand, and E. Gutmark. On the Existence and Multiplicity ofRotating Detonations. Proceedings of the Combustion Institute, 36:2691–2698, 2017. doi:
10.1016/j.proci.2016.06.132.
R. W. Houim and R. T. Fievishon. The Influence of Acoustic Impedance on Gaseous Layered
Detonations Bounded by an Inert Gas. Combustion and Flame, 179:185–198, 2017. doi:
10.1016/j.combustflame.2017.02.001.
K. Ishihara, Y. Kato, K. Matsuoka, J. Kasahara, A. Matsuo, and I. Funaki. Performance Evaluationof Rotating Detonation Engine with Conical-Shape Tail. In 53rd AIAA Aerospace Sciences
Conference, Kissimmee, Florida, January 2015. doi: 10.2514/6.2015-0630.
K. Kailasanath. Recent Developments in the Research on Rotating-Detonation-Wave Engines. In
55th AIAA Aerospace Sciences Conference, Grapevine, TX, January 2017. doi: 10.2514/6.
2017-0784.
Y. Kato, K. Ishihara, K. Gawahara, K. Matsuoka, J. Kasahara, A. Matsuo, and I. Funaki. Thrust
Measurement of Rotating Detonation Engine by Sled Test. In 50th AIAA/ASME/SAE/ASEE Joint
Propulsion Conference, Cleveland, OH, July 2014. doi: 10.2514/6.2014-4034.
Y. Kato, K. Ishihara, K. Matsuoka, J. Kasahara, A. Matsuo, and I. Funaki. Study of Combustion Chamber Characteristic Length in Rotating Detonation Engine with Convergent-Divergent
Nozzle. In 54th AIAA Aerospace Sciences Conference, San Diego, CA, January 2016. doi:
10.2514/6.2016-1406.
K. Kawane, S. Shimada, J. Kasahara, and A. Matsuo. The Influence of Heat Transfer and Friction
on the Impulse of a Detonation Tube. Combustion and Flame, 158:2023–2036, 2011. doi:
10.1016/j.combustflame.2011.02.017.
J. H. S. Lee. The Detonation Phenomenon. Cambridge University Press, New York, NY USA,2008.
F. A. Lyman. On the Conservation of Rothalpy in Turbomachines. Journal of Turbomachinery,
115:521–526, July 1993. doi: 10.1115/1.2929282.
J. A. Nicholls, R. E. Cullen, and K. W. Ragland. Feasibility Studies of a Rotating Detonation Wave
Rocket Motor. Journal of Spacecraft and Rockets, 3(6):893–898, 1966. doi: 10.2514/3.28557.
M. I. Radulescu and R. K. Hanson. Effect of Heat Loss on Pulse-Detonation-Engine Flow Fields andPerformance. Journal of Propulsion and Power, 21(2):274–285, 2005. doi: 10.2514/1.10286.
B. A. Rankin, M. L. Fotia, A. G. Naples, C. A. Stevens, J. L. Hoke, T. A. Kaemming, S. W.
Theuerkauf, and F. R. Schauer. Overview of Performance, Application, and Analysis of Rotating
Detonation Engine Technologies. Journal of Propulsion and Power, 33(1):131–143, January
2017. doi: 10.2514/1.B36303.
D. Schwer and K. Kailasanath. Effect of Inlet on Fill Region and Performance of Rotating Detonation
Engines. In 47th AIAA/ASME/SAE/ASEE Joint Propulsion Conference and Exhibit, San Diego,CA, 2011a. doi: 10.2514/6.2011-6044.

Shepherd and Kasahara - September 23, 201728

D. Schwer and K. Kailasanath. Numerical Investigation of the Physics of Rotating-DetonationEngines. Proceedings of the Combustion Institute, 33:2195–2302, 2011b. doi: 10.1016/j.proci.
2010.07.050.
D. Schwer and K. Kailasanath. Fluid Dynamics of Rotating-Detonation-Engines with Hydrogenand Hydrocarbon Fuels. Proceedings of the Combustion Institute, 34:1991–1998, 2013. doi:
10.1016/j.proci.2012.05.046.
P. I-W Shen and T.C. Adamson, Jr. Theoretical Analysis of a Rotating Two-Phase Detonation in
Liquid Rocket Motors. Astonautica Acta, 17(4-5):715–728, 1972. ISSN 0004-6205.
M. Sichel and J. C. Foster. The Ground Impulse Generated by a Plane Fuel-Air Explosion withSide Relief. Acta Astronautica, 6(3-4):243–256, 1979. doi: 10.1016/0094-5765(79)90096-1.
G. P. Smith, D. M. Golden, M. Frenklach, M. W. Moriarty, B. Eiteneer, M. Goldenberg, C. T.
Bowman, R. K. Hanson, S. Song, W. C. Gardiner, V. V. Lissianski, and Z. Qin. GRI-Mech, July
1999. URL http://www.me.berkeley.edu/gri_mech/.
G. P. Sutton. Rocket Propulsion Elements. Wiley, New York, NY USA, 6th edition, 1992.
N. Tsuboi, E. Seiichiro, A. K. Hayashi, and T. Kojima. Front Cellular Structure and Thrust
Performance on Hydrogen-Oxygen Rotating Detonation Engine. Journal of Propulsion andPower, 33(1):100–111, 2017. doi: 10.2514/1.B36095.
E. Wintenberger and J. E. Shepherd. The Stagnation Hugoniot Analysis for Steady Combustion
Waves in Propulsion Systems. Journal of Propulsion and Power, 22(4):835–844, 2006. doi:
10.2514/1.12779.
E. Wintenberger, J. M. Austin, M. Cooper, S. Jackson, and J. E. Shepherd. An Analytical Modelfor the Impulse of a Single-Cycle Pulse Detonation Tube. Journal of Propulsion and Power, 19(1):22–38, 2003. doi: 10.2514/2.6099.
E. Wintenberger, M. Cooper, F. Pintgen, and J. E. Shepherd. Reply to Comment on "Analytical
Model for the Impulse of Single-Cycle Pulse Detonation Tube" by M. I. Radulescu and R. K.
Hanson. Journal of Propulsion and Power, 20(5):957–959, 2004. doi: 10.2514/1.9441.
T.-H. Yi, J. Lou, C. Turangan, J.-Y. Choi, and P. Wolanski. Propulsive Performance of a Continuously
Rotating Detonation Engine. Journal of Propulsion and Power, 27(1):171–181, 2011. doi:
10.2514/1.46686.
R. Zhou, D. Wu, and J.-P. Wang. Progress of Continuously Rotating Detonation Engines. Chinese
Journal of Aeronautics, 29(1):15–29, 2016. doi: 10.1016/j.cja.2015.12.006.

Shepherd and Kasahara - September 23, 2017

A

29

Equations of Motion

The equations of motion for inviscid compressible flow can be written as:∂ρ+ ∇ · (ρu) = 0 ;∂t∂u+ ρ(u · ∇)u = −∇P ;
ρ
∂t 



∂
1 2
1 2
∂P
ρ
h + |u| + ρ(u · ∇) h + |u| = −
.
∂t
2
2
∂t

(65)
(66)
(67)

In terms of the components of velocity (u, v, w) and the cylindrical coordinate system (r, θ, z) shownin Fig. 7, these equations are:
∂
1 ∂
∂
∂ρ
+ (ρu) +(r ρv) + (ρw) = 0 ;
∂t ∂r
r ∂θ
∂z
∂u
∂u v ∂u∂u v 2
1 ∂P
+u
+
+w
−
=−
;
∂t
∂r r ∂θ
∂z
r
ρ ∂r
∂v
∂v v ∂v
∂w uv
1 ∂P
+u
+
+w
+
=−
;
∂t
∂r r ∂θ
∂z
r
ρr ∂θ
∂w
∂w v ∂w
∂w
1 ∂P
+u
+
+w
=−
;
∂t
∂r
r ∂θ
∂z
ρ ∂z
 



∂
∂
v ∂
∂
1 ∂P
1 2
1 2
+w
.
h + |u| + u +h + |u| = −∂t2∂r r ∂θ
∂z
2
ρ ∂t

(68)
(69)
(70)
(71)
(72)

The relation to the 2D, “unwrapped” coordinate is through the approximations of neglecting motionin the r direction by setting u = 0 and taking ∂/∂r = 0. This results in an inconsistency for theradial momentum which reduces in the limit u = 0 tov 2 1 ∂P
=
r
ρ ∂r

(73)

and because v does not vanish, there must be a radial pressure gradient. However, if the radiusof the annulus R is sufficiently large, the pressure gradient in the radial direction will be smallin comparison to the azimuthal and axial gradients and can be neglected. Formally this can beaccomplished by taking the limit r → ∞ in the 3D equation set.
To transform to wave-fixed coordinates in 2D, define the following coordinate and velocity
transformations using the mean radius R of the combustion chamber.
x 0 = UD t − Rθ ;z0 = z ;t0 = t ;

(74)
(75)
(76)

u0 = UD − v ;w0 = w .

(77)
(78)

and

Shepherd and Kasahara - September 23, 2017

30

This transformation corresponds to going from the laboratory frame (θ-Z) representation in Fig. 1b
to the unwrapped wave-fixed frame (X 0-Z 0)in Fig. 3. The velocity components in the unwrappedframe are (u0, w0) corresponding to the coordinates (x 0, z0). The detonation is located at x 0 = 0 andthe thrust surface is located at z0 = 0 in the transformed variables. The quantities ρ, P and h are
invariant, i.e., the same in both coordinate systems. The equations of motion transform exactly to
the 2D steady-flow compressible (Euler) equations in this limit:
∂
∂
0
(ρu
)
+
(ρw0) = 0 ;
0
0
∂x
∂z
0
∂u0
1 ∂P
∂u
;
u0 0 + w0 0 = −∂x∂z
ρ ∂ x01 ∂P∂w0∂w0
u0 0 + w0 0 = −
;
∂x
∂z
ρ ∂z0



∂
1
∂
u0 0 + w0 0 h + (u02 + w02 ) = 0 .
∂x
∂z
2

(79)
(80)
(81)
(82)

One implication of the energy equation is that the total enthalpy is constant along streamlines in
the wave-fixed frame,1
ht = h + (u02 + w02 ) = constant .
2(83)

In this frame of reference, the flow is steady and total enthalpy is conserved across both shock
and detonation waves enabling evaluation of the total enthalpy at any convenient point along a
streamline. In the flow just ahead of the detonation h = h1 , w0 = 0 and u0 = UD so that1ht = h1 + UD2 .
2(84)

If the flow ahead of the detonation has uniform properties, then this is the unique value of stagnation
enthalpy for the entire flow field. Transforming (83) to the laboratory reference frame using (77),the total enthalpy is1ht = h + (v 2 + w 2 + UD2 ) − vUD = constant .
2(85)

Subsituting for the value of ht from (84) we obtain the alternative expression for total enthalpy
conservation, which is used for the axial flow model discussed in Section 4,1
h1 = h + (v 2 + w 2 ) − vUD .
2(86)

In the limit as the azimuthal velocity vanishes, v → 0, the axial flow velocity is given by
w| axial =p2(h1 − h) .

(87)Shepherd and Kasahara - September 23, 201731

Although at first glance this appears to be inconsistent with the same limit for (83), in the co-rotating
frame, conservation of angular momentum requires that u0 → UD in the limit of purely axial flowin the laboratory frame. Substituting (84) in (83), we find that1h1 = h + (w02 + u02 − UD2 ) ,2(88)which gives a result identical to (87) in the limit of purely axial flow, recognizing that w0 = w.
In the context of turbomachinery, the expression on the right-hand side of (86) is referred to as
the rothalpy and UD is tangential speed of the rotor. A detailed discussion of this concept is given
by Lyman (1993).

Shepherd and Kasahara - September 23, 2017

B

32

Solutions to Triple-Point Configuration

The pressure-deflection states shown in Fig. 4 were computed using equilibrium thermochemistry to
determine the detonation state, oblique shock state and Prandtl-Meyer expansion with the methods
described in Browne et al. (2017). Selected numerical values associated with each state are given in
Table 2. Initial conditions for all cases are stoichiometric mixtures with P1 = 0.1 MPa and T1 = 300
K; velocities are all given in the wave-fixed frame. State 4 was determined by isentropic expansionof the detonation state 2 to the initial pressure P1 assuming one-dimensional flow parallel to thex axis. The solutions for states 3 and 5 were determined by numerical solution of the pressureand deflection matching conditions at the 3-5 contact surface. The reported values of enthalpy
are mixture averages including the heats of formation for each species with the concentrations
determined by the equilibrium solutions. The 53 species found in GRI-Mech 3.0 (Smith et al., 1999)
were used for these computations with the thermodynamic data fits corrected for discontinuities
and the temperature range extended to 6000 K.

33Shepherd and Kasahara - September 23, 2017

Table 2: Triple-point and axial-flow solution states for representative fuel-air and fuel-oxygen
mixtures.
H2 -air C2 H4 -air H2 -O2 C2 H4 -O2
State 1 - Initial Conditions
ρ1(kg/m3 ) 0.838
1.15
0.482
1.24
h1
(kJ/kg)
2.61
122.
4.47
425
u1 (UC J )
(m/s)
1968
1823
2836
2373
a1
(m/s)
409
346
539
328
State 2 - CJ
T2
(K)
2943
2923
3675
3930
P2
(MPa)
1.55
1.82
1.86
3.31
3
ρ2
(kg/m ) 1.51
2.09
0.884
2.30
h2
(MJ/kg) 1.34
1.28
2.83
2.42
a2
(m/s)
1091
1005
1542
1281
γ
1.163
1.161
1.129
1.139
State 3 - PM Expansion
T3
(K)
2479
2455
3288
3423
P3
(MPa)
0.429
0.470
0.461
0.626
ρ3
(kg/m3 ) 0.506
0.655
0.255
0.526
h3
(MJ/kg) 0.145
0.204
0.104
0.234
θ3
(deg)
28.6
30.3
32.5
38.1
u3
(m/s)
1895
1777
2800
2451
a3
(m/s)
1002
918
1422
1156
State 4 - Upstream of Oblique Shock
T4
(K)
1948
1922
2935
2993
3
ρ4
(kg/m ) 0.151
0.180
0.0646
0.101
h4(MJ/kg) -0.952-0.780
-2.45
-1.75
u4
(m/s)
2405
2264
3598
3160
a4
(m/s)
902
825
1310
1045
State 5 - Downstream of Oblique Shock
T5
(K)
2602
2593
3346
3528
3
ρ5
(kg/m ) 0.480
0.616
0.245
0.489
h5
(MJ/kg) 0.479
0.548
1.07
1.38
u5
(m/s)
1710
1572
2429
1929
a5
(m/s)
1024
942
1451
1200
Axial Sonic Flow
∗
w
(m/s)
951
877
1376
1130
P∗
(MPa)
0.203
0.236
0.253
0.426
T / MÛ
(m/s)
1347
1297
2117
1901

C3 H8 -O2
1.36
-506.
2356
308
3822
3.59
2.53
1.46
1269
1.134
3324
0.653
0.557
-0.736
39.0
2452
1145
2898
0.103
-2.74
3164
1036
3439
0.515
0.482
1891
1193
1123
0.460
1908

34

Shepherd and Kasahara - September 23, 2017C

Existence of Attached Oblique Solutions

Houim and Fievishon (2017) call attention to the role of acoustic impedance ratio Z across theinterface between the detonation layer and inert bounding gas, using two-dimensional unsteadysimulations with a single-step model of detonation. They define the acoustic impedance ratio Z
in terms of the ratio of the impedance ρa of the reactants to that of the bounding gas layer, both
evaluated upstream of the wave complex. Using the notation of Fig. 3, this is
Z14 =

(ρa)4
.
(ρa)1

(89)

p
For constant property model used in Houim and Fievishon (2017), this reduces to Z = T1 /T4 andvalues of Z = 1.73, 1, 0.71, 0.55 and 0.29 were used in their simulations (corresponding to T4 =100, 300, 600, 1000 and 3500 K for T1 = 298 K. They found that in all cases except for T4 = 3500 K,
a solution was obtained with an attached oblique shock wave, corresponding to the configurations
shown in Fig. 4. For case 5 (T4 = 3500 K) of Houim and Fievishon (2017), the solution consisted ofa detached precursor shock (oriented normal to the interface and propagating in the inert bounding
layer) followed by an oblique shock in the cold reactant layer that intersects a nominally normal
detonation wave in the reactant layer. A secondary oblique shock is needed to match the conditions
associated with the PM expansion within the detonation products. Although these cases generally
still result in an expansion wave following the detonation, the scaled variation of pressure with
distance along the thrust surface may be significantly different for detached and attached oblique
shock solution. Solutions of this type (normal precursor shock in the inert bounding layer) will be
unsteady because the flow behind the detached normal shock is subsonic relative to the shock wave,
resulting in acoustic feedback coupling the precursor shock propagation with the flow downstreamof the triple point. This is in contrast to the attached oblique shock cases discussed in Section 3;
these all have supersonic flow behind the oblique shock and disturbances generated by the flow
downstream of the triple point do not propagate upstream.
Simulations by Houim and Fievishon (2017) raise the interesting possibility that more complex
flowfields may result than shown in Fig. 1 or 3. Analyses of triple-point configuration for cases
similar to those examined by Houim and Fievishon (2017) are shown in Fig. 13. In Fig. 13a,
all cases except for 3500 K (comparable to Case 5 of Houim and Fievishon (2017)) display anintersection between the oblique shock polar and PM expansion behind the detonation. These cases
admit a mean flow solution with an attached oblique shock, which is consistent with the unsteadysimulations of Houim and Fievishon (2017). The polars do not intersect for the 3500 K case,consistent with the detached shock solution observed in their unsteady simulations.
There are several key differences between the simulations of Houim and Fievishon (2017) and
the RDE situation. One key difference is that the expanded detonation products are moving relative
to the reactant layer, resulting in shear flow between the inert gas layer and reactants. Parametriccomputations with four layer speeds are shown in Fig. 13 for detonation of a stoichiometric C2 H4 O2 layer (T1 = 300 K, P1 = 1 atm) bounded by a layer of air at various initial temperatures andvelocities. These solutions show that attached oblique shock solutions are possible if the bounding
gas speed is sufficiently high; the speed of the flow of the bounding inert gas layer is given relative
to the reactant layer in the direction opposed to the detonation motion in the laboratory frame. The
velocity required to obtain attached shock solutions are comparable to those of state 4 in the model
RDE cases we have examined, see Table 2. Other significant factors include the bounding gas layer

35

Shepherd and Kasahara - September 23, 2017

temperature, which is much lower for fuel-air than fuel-oxygen cases, and significant differences inthe value of γ, the ratio of specific heats compared to what was used in the Houim and Fievishon(2017) simulations. The consequences of these differences is that for RDEs, attached oblique shock
solutions such as shown in Fig. 4 appear to be possible in all the cases we have examined.

(a)

(b)

Figure 13: Pressure-flow deflection analyses. (a) Effect of air temperature (b) Effect of shear flowfor 3500 K air temperature.
All other factors being equal, the existence of attached shock solutions can be correlated withthe ratio Z14 but is it clear from examining the behavior of the shock polars in Fig. 13 that this iscaused by decrease in the inert layer Mach number with increasing sound speed for fixed U1 andT1 , resulting in the shift of the shock polar relative to the PM expansion. An even larger shift in
the relative location of the PM and oblique shock polar occurs due to the change in the reactant
composition from fuel-oxygen to fuel-air, as shown in Fig. 4. In general, the existence of attachedshock solutions will depend on M1 , γ1 , M4 and γ4 (using one-γ approximations for each layer) andis not correlated with Z14 alone.
Houim and Fievishon (2017) unsteady simulations indicate that transmission of the transverse
waves across the boundary between detonation products (state 3) and shocked inert gas (state 5)
does appear to be a significant factor is determining if detonations will successfully propagate
within the reactant layer. If we treat the interactions of the transverse waves with the 3–5 interface
in the acoustic limit, then the ratio
(ρa)5
(90)
Z35 =
(ρa)3
will be relevant to the contact surface interaction of transverse waves associated with detonation
instability rather than Z14 . From Table 2, values of Z35 are 1; whereas the values for the contact
surface upstream of the detonation Z14 are between 0.25 and 0.40. This difference between upstream
and downstream values of Z can explain the confinement of transverse waves to the detonation
layer in case 5 of Houim and Fievishon (2017) and stability of detonation propagation despite thelack of a stable attached oblique shock solution.