Numerical study of pressure profiles in an annular combustor and a manifold during continuous spin detonation of nonstoichiometric hydrogen-air mixtures

. A closed mathematical model of continuous spin detonation of nonstoichiometric hydrogen-air mixtures in a quasi-three-dimensional unsteady gas-dynamic formulation taking into account the reverse effect of the oscillation processes in the combustor on the mixture injection system is utilized in pressure profiles in the combustion chamber and in the manifold investigation. Numerical simulations are performed for the flow-type annular combustor 503 mm in diameter with an annular gap of width 18 mm and other geometric dimensions corresponding to experiments. Three-wave mode of continuous spin detonation is calculated and the flow structure is analyzed for the equivalence ratio of 0.6 and the specific hydrogen-air mixture flow rate of 390.6 kg/(s·m 2 ). The static, the mean static and the mean total pressure distributions in the combustor are calculated. The results are validated against experimental data.


Introduction
The results of studies of continuous spin detonation (CSD) based on the scheme proposed by Voitsekhovskii [1] are presented in [2,3].CSD of a hydrogen-air mixture was obtained experimentally in [4] and numerically investigated in [5].The calculations of CSD were performed for a stoichiometric hydrogen-air mixture with the equivalence ratio  = 1 [5].Shock waves were shown entering the manifold and producing a significant effect on the entire CSD process.The new series of experiments for nonstoichiometric hydrogenair mixtures [6] allowed to measure pressure profiles in the combustion chamber and in the manifold in detail.The objectives of this study are to utilize a closed mathematical model of CSD of hydrogen-air mixtures in a quasi-three-dimensional unsteady gas-dynamic formulation developed in [5] to numerically study the CSD dynamics of mixture with the equivalence ratio  = 0.6 in a ramjet annular combustor with a diameter of 503 mm, to obtain pressure profiles in the combustor and the manifold, and to compare the results with experiments [6].

Formulation of the problem 2.1 Mathematical model
Let us consider the problem of mathematical modeling of detonation combustion of a nonstoichiometric hydrogen-air mixture in a flow-type annular cylindrical combustor (diameter dc, length Lc, expansion of the annular combustor channel from δ to Δ) in the case of mixture exhaustion from a manifold (length Lm, contraction of the channel from Δm to δ).To compare the calculated parameters with experimental data, the combustor geometry was taken from [6].The mixture is injected from a high-pressure receiver through injectors in the end face of an annular manifold.For experimental annular combustors [4,6], the following inequalities hold: δ < Δ << dc/2, Δm << dc/2.Therefore, the problem can be considered in a quasi-three-dimensional approximation, similar to [5].The annular domain is cut and unfolded into a rectangular solution domain Ω1  Ω2 shown in Figure 1.end of the combustor through which detonation products escape to the outer space, l is the period of the problem.
Let a certain amount of energy sufficient for detonation initiation be released in the solution domain area Ω3 at a certain time instant after the hydrogen-air mixture is supplied to the combustor through the boundary Γ1.An unsteady detonation wave propagates in the domain Ω2 as a result of initiation.The obtained hydrogen-air reacting mixture flow in the computational domain Ω1  Ω2 is described by the following system of equations of non-stationary gas dynamics with chemical reactions [5]: Here t is the time,  is the density, u and v are the velocity vector components; p is the pressure; E = U + (u 2 + v 2 )/2, where U(T, ) is the total internal energy of the gas, T is the temperature ,  is the current molar mass of the mixture, and Y is the fraction of the chemical induction period.
The cross-sectional area of the manifold and combustor channel S = S(x) along the x coordinate with smooth sinusoidal convergence (L < x < 0) to the width  and subsequent smooth sinusoidal divergence (0 < x < L) to the combustor channel width  is defined in the form: As in [5] the energy release is described within the framework of the two-stage kinetic model [7]: induction stage (0  Y  1, f5 =  1/tind, f6 = 0); chemical transformation stage (Y = 0, f5 = 0, f6  0).The chemical delay of ignition of the hydrogen-air mixture in the induction region (0  Y  1) has the form: where a  = 17.15 kcal/mol, Кa = 5.3810 11 (mols/l), R is the universal gas constant, O2 is the molar mass of oxygen, z = O2/(O2 + 2H2 + in) is the mass fraction of oxygen in the hydrogen-air mixture, H2 is the molar mass of hydrogen,  is the fuel-to-oxidizer equivalence ratio,  = 3.772, in = 28.148kg/kmol.The structural formula for the hydrogen-air mixture is presented as 2Н2 + O2 + 3.7275 N2 + 0.0445 Ar.System (1) is supplemented with the equations of state: where Uth, Uch are the thermodynamic and chemical components of internal energy, respectively, which were determined in a similar manner as in [5,8].
At the chemical transformation stage (Y = 0, f5 = 0, f6  0) the rate of energy release is determined by the rate of chemical reactions [8]: where K+ is the generalized recombination rate constant, K-is the equilibrium constant, T0 is the initial temperature of the mixture,  = 1+ max/(max/min  1) , min and max are the molar masses of the gas in the ultimately dissociated, and ultimately recombined states, max is the molar fraction of triatomic molecules in the ultimately recombined state;  is the effective temperature of excitation of vibrational degrees of freedom of molecules, Ed is the mean energy of dissociation of reaction products.The quantity K-was calculated in accordance with the equation of chemical equilibrium for the corresponding composition.System (1)-( 5) is closed and completely determines the unsteady motion of the hydrogen-air reacting mixture with variable heat release in the reaction zone behind the detonation wave.

Boundary conditions
Injection of the combustible mixture into the manifold through the entrance boundary Γ0 (x = -Lm; 0 ≤ y ≤ l) is modelled by the income through a system of Laval micronozzles similar to [5].
At the combustor exit (boundary Γ2: x = Lc; 0 ≤ y ≤ l) the boundary fluxes of mass, momentum, and energy due to the interaction of the gas with the environment at pressure p = pa were determined using the method of [10].
On the left and right boundaries of the domain Ω1  Ω2, the condition of the periodicity of the solution (with period l) was used:
Constants included into the description of the thermodynamic properties of the gas and its combustion products for the hydrogen-air mixture with the value of  = 0.6 are: z = 0.2276, 0 = 23.539kg/kmol, min = 17.202 kg/kmol, max = 26.168kg/kmol,  = 3000 K,  = 1.4286,К = 2091.8kmol/m 3 .As in [5], the chemical equilibrium constant К is determined from the chemical kinetics equation ( 5) by substituting the Chapman-Jouguet detonation parameters at f6 = 0 determined using the method of [9]: DCJ = 1.707 km/s, TCJ = 2427.3K, PCJ = 13.13 atm, CJ = 26.106kg/kmol.Here DCJ is the velocity, TCJ is the temperature, PCJ is the pressure and CJ is the molecular weight of the ideal Chapman-Jouguet detonation calculated in the hydrogen-air mixture at the nonstoichiometric equivalence ratio  equal to 0.6.

Governing parameters of the model
For fixed initial constants of the model, geometric parameters of the combustor and the manifold (7), the solution of the non-stationary CSD problem ( 1)-( 6) depends on the governing parameters: p * , T * , S*, , pa, l (8) Here p * and T * are the stagnation pressure and temperature in the mixture injection system, S* is the total cross-sectional area of the micronozzle throat at the manifold entrance,  is the equivalence ratio, pa is the ambient pressure, and l is the period of the problem along the y axis.

Calculation results
Problem (1)-( 7) is solved numerically.The solution domains Ω1 and Ω2 are covered by a motionless grid with uniform cells in the y direction (the number of cells is 700) and nonuniform cells in the x direction (the number of cells is 140 and 300 correspondingly).Similar to [5], equations are integrated by the Godunov-Kolgan second-order finite difference scheme [10,11].
To compare the calculated in the computational domain Ω1  Ω2 parameters with experimental data, the following two-dimensional non-stationary periodic problem with a period l = Π/3 = π(dc -Δ)/3 = 50.74cm and the parameters in the injection system p * /p0 = 8.05, T * /T0 = 0.89, S*/Sm = 0.12,  = 0.6 with the specific flow rate of the mixture gΣ0 = 390.6 kg/(s•m 2 ) corresponding to the experiments [6] was solved.Here Π = π(dc -Δ) is the mean combustor perimeter, Sm = π(dc -Δm)•Δm = 473.26cm 2 .For comparisons of the flow pattern with the three-wave (n = 3) CSD regime [6], the parameter l is assumed to be equal to Π/3.At the combustor end the counterpressure pa/p0 = 1 was specified.The nonstoichiometric hydrogen-air mixture at rest with the fuel-to-oxidizer equivalence ratio  = 0.6 was used as initial parameters in the computational domain.
A certain amount of energy sufficient for detonation initiation is instantaneously released at the initial time (t = 0) in the domain Ω3 (see Figure 1).The calculations of the initiation of a transverse detonation wave (TDW) and the dynamic of its formation show that by the time t ≈ 18.7 ms, self-sustained CSD regime is established.

Pressure profiles
The time dependence of the dimensionless static pressure P = p/p0 at a various fixed points of the chamber and the manifold is presented in Figure 2. One can see periodic oscillations with a period Δt ≈ 0.343 ms, which corresponds to a rotation frequency of the TDW f ≈ 2.915 kHz, detonation velocity D ≈ 1.48 km/s, and the ratio D/DCJ = 0.87.In the process of the TDW propagation the pressure peak is firstly reached for P1, then for Pm,1 and then for Pm,2.Note that the value of Pm,2 is slightly higher than the value of Pm,1.This is explained by the fact of location of the point (-13.3,0) close to the boundary Γ0, where the pressure wave is reflected from the end face of the manifold.It is seen that the static pressure Pex also has periodic fluctuations corresponding to the moment when an oblique shock wave (tail) departed downstream from the TDW reaches the point at the combustor end.The peaks of Pex are 3-4 times lower than that in the upper part of the combustor (P1), as the shock wave passes over the detonation products and becomes weaker.

TDW structure
The distribution of the calculated dimensionless temperature T/T0 for the three-wave (n = 3) CSD regime at time t = 18.7 ms is presented in Figure 3. Three periodic solutions calculated for l = Π/3 are presented side by side.The upper part of Figure 3 (x < 0) illustrates the gas-dynamic flow in the manifold (Ω1), and the lower part of Figure 3 (x > 0) shows the flow in the combustor (Ω2).It is evident from the temperature field that the TDW front moves from left to right along the triangle of the incoming cold mixture.An oblique shock wave moves down to the left.Behind the TDW, the detonation products gradually expand, the pressure of the products decreases, becomes smaller than the injection pressure in the manifold, the products are displaced downward by the cold hydrogen-air mixture injected from the manifold through the boundary Γ1.For the present variant with gΣ0 = 390.6 kg/(s•m 2 ), the attempt to continue the computations with the problem period reduced to l = Π/4 (n = 4) resulted in failure and decay of the rotating TDW.According to the classification [3] the above-described numerical solution for the three-wave CSD in the combustor with parameters (7) is the only possible solution.

Comparison with experiments
In addition to static pressures, the mean static and the mean total ones were also calculated in the chamber and the manifold.The computed and the experimental data [6] are summarized in Table 1, where n is the number of TDWs accommodated over the combustor perimeter, D is the TDW velocity, h is the TDW height; <p1>, <p2>, <pex> and <pm,in> are the mean static pressures near the combustor entrance (x = 1.5 cm), (x = 6.1 cm), at the combustor exit (x = 60 cm) and the end face of the manifold (x = -13.9cm) respectively; p0,ex is the mean total pressure at the combustor exit.The results show good agreement with the experiment for the number of waves in the combustor and the TDW velocity, significant difference in the TDW height (50% lower than the experimental one).This is explained by the ideal model of mixing of the components of the mixture as in [5].*Note that the value of <p2> is measured for (x = 8.1 cm) in the experiment and for (x = 6.1 cm) in the present calculations.The point for calculations of <p2> is shifted by 2 cm towards the combustor entrance to compensate the lower TDW height than in the experiment.
The values of the mean static pressures correlate with the experiments.The difference between the calculated and experimental data is no more than 10%.The mean total pressure at the combustor exit p0,ex is 20% higher than the experimental one.This may be due to the fact that the considered mathematical model does not take into account non-stationary processes which lead to the total pressure losses.

Conclusions
A closed quasi-three-dimensional non-stationary mathematical model of continuous spin detonation for a hydrogen-air mixture in a flow-type annular combustor [5] with the chemical kinetics equation [8] is utilized in numerical simulation for geometric dimensions of the combustor and the manifold corresponding to the experiment [6].At a specific flow rate of the mixture gΣ0 = 390.6 kg/(s•m 2 ), the three-wave CSD regime with the TDW velocity D ≈ 1.48 km/s and the height of the front h ≈ 7.5 cm is obtained.The TDW dynamics and its twodimensional structure are numerically studied.The static, the mean static and the mean total pressures at the various points of the chamber and the manifold for the three-wave CSD regime are obtained.Comparison with experiments show complete agreement on the number of waves in the combustor (n = 3), good agreement on the TDW velocity, significant difference in the TDW height (50% lower).The difference between the calculated and the experimental data for the mean static pressures is no more than 10%.The mean total pressure at the combustor exit p0,ex is 20% higher than the experimental one.Thus, the considered mathematical model [5] is able to adequately describe the distribution of static, mean static and mean total pressures in the flow-type annular combustors.

Fig. 1 .
Fig. 1.Schematic of the annular cylindrical combustor model and the solution domain of the periodic problem Ω1  Ω2.