4.1. Governing equations

In this study, neglecting the effects of heat conduction, molecular diffusion and viscosity, the governing equations are the two-dimensional reactive Euler equations coupled with a simple two-step chain-branching reaction model (Short & Sharpe Reference Short and Sharpe2003; Ng, Radulescu & Higgins Reference Ng, Radulescu, Higgins, Nikiforakis and Lee2005), which can be written as

(4.1) \begin{align} \frac {\partial U}{\partial t} + \frac {\partial F}{\partial x} + \frac {\partial G}{\partial y} = S, \end{align}

where $U$ is the vector of conserved variables, $F$ and $G$ are both the fluxes, and $S$ is the reaction source term. These can be written as

(4.2) \begin{align} U &= [\rho , \rho u , \rho v , E , \rho y_{\textit{I}} , \rho y_{\textit{R}} ], \nonumber\\ F &= \big[\rho u , \rho u^{2} + p , \rho u v , u \left (E + p\right )\! , \rho u y_{\textit{I}} , \rho u y_{\textit{R}}\big], \nonumber\\ G &= \big[\rho v , \rho u v , \rho v^{2} + p , v \left (E + p\right )\! , \rho v y_{\textit{I}} , \rho v y_{\textit{R}}\big], \nonumber\\ S &= \left [0,0 , 0,0 , \omega _{\textit{I}} , \omega _{\textit{R}}\right ]\!, \end{align}

where, $p$ , $\rho$ are pressure and density, respectively, $E$ is energy per unit volume, and $u$ and $v$ are the velocities in the $x$ and $y$ directions, respectively. The present two-step chemical reaction model has a thermally neutral chain-branching induction zone and an exothermic heat-releasing zone. Here, $y_{\textit{I}}$ and $y_{\textit{R}}$ represent the progress variables characterising the induction and exothermic processes in the two-step reaction model, and $\omega _{\textit{I}}$ and $\omega _{\textit{R}}$ represent the reaction rates of the induction and reaction processes, respectively.

For an ideal gas, the total energy and state equation can be written as

(4.3) \begin{align} &E = \frac {p}{\gamma - 1} + \frac {\rho }{2} \big(u^{2} + v^{2}\big) - \rho y_{R} Q, \nonumber \\ &p= \rho R T, \end{align}

where $Q$ is heat release of the chemical reaction, $T$ is the temperature, and $R$ is the universal constant of gas.

The extent of the induction zone is controlled by a reaction rate of Arrhenius form. At the end of the induction zone, the conversion of reactant into radical is finished. The reaction rate of the induction stage can be expressed as

(4.4) \begin{align} \omega _{\textit{I}} = - H \!\left (1 - y_{\textit{I}}\right ) k_{\textit{I}} \rho \exp \left [E_{\textit{I}}\left (\frac {1}{R T_{\textit{S}}} - \frac {1}{R T}\right )\right ]\!, \end{align}

where $E_{\textit{I}}$ is the activation energy of the induction stage, $T_{\textit{S}}$ is the temperature behind the leading shock of the ZND detonation with CJ speed, $k_{\textit{I}}$ is a pre-exponential factor to control the induction length $\varDelta _I$ , and $H(1-y_{\textit{I}})$ is a step function, written as

(4.5) \begin{align} H(1 - y_{\textit{I}}) = \begin{cases} 1, & {\textrm{if }} y_{\textit{I}} \lt 1, \\ 0, & {\textrm{if }} y_{\textit{I}} \geqslant 1. \end{cases} \end{align}

The rate of the exothermic process where the radical converts to the product can be expressed as

(4.6) \begin{align} \omega _{\textit{R}} = \left [ 1 - H(1 - y_{\textit{I}}) \right ] k_{\textit{R}} \rho (1 - y_{\textit{R}}) \exp \left (\frac {E_{\textit{R}}}{R T}\right )\!, \end{align}

where $E_{\textit{R}}$ is the activation energy in the exothermic stage, and $k_{\textit{R}}$ is a pre-exponential factor to control the reaction length $\varDelta _R$ . The change in the ratio of the induction length $\varDelta _I$ to the reaction length $\varDelta _R$ can affect the instability of the detonation.

To compare with the existing literature, a scaled activation energy parameter is given here as

(4.7) \begin{align} \epsilon _{\textit{I}} = \frac {E_{\textit{I}}}{R T_{\textit{S}}}, \quad \epsilon _{\textit{R}} = \frac {E_{\textit{R}}}{R T_{\textit{S}}}, \end{align}

where $T_{\textit{S}}$ can be given as

(4.8) \begin{align} T_{\textit{S}} = \frac {\left (2 \gamma M_{\textit{CJ}}^{2} - \gamma + 1\right ) \left [2 + \left (\gamma - 1\right ) M_{\textit{CJ}}^{2}\right ]}{\left (\gamma + 1\right )^{2} M_{\textit{CJ}}^{2}}. \end{align}

The activation energies were defined to represent a typical hydrocarbon–oxygen mixture with weakly unstable stability, as discussed in the work of Ng et al. (Reference Ng, Radulescu, Higgins, Nikiforakis and Lee2005). The Arrhenius rate coefficients were calculated to make the induction length 1 mm, which is quite convenient for scaling the geometry with the detonation characteristic length. Calculating the detonation properties of stoichiometric ${\textrm{H}}_2{-}{\textrm{O}}_2{-}{\textrm{Ar}}$ with dilution from 0 $\%$ to 20 $\%$ using a detailed chemical reaction mechanism reveals that the activation energies range from 4 to 6, and the Mach number of detonations is approximately 5. Thus the present model and parameters mimic the stoichiometric ${\textrm{H}}_2{-}{\textrm{O}}_2{-}{\textrm{Ar}}$ with dilution from 0 $\%$ to 20 $\%$ with initial pressure below 50 kPa. The reaction and thermodynamic parameters are calculated and listed in table 1, corresponding to a weakly unstable detonation.

Table 1.

Reaction and thermodynamic parameters.

4.2. Numerical methods

For the two-dimensional reactive Euler equation, a second-order Strang operator splitting method effectively decouples the convection and chemical reaction source terms. The convection term is discretised using the explicit second-order MUSCL–Hancock method based on the Roe-HLL solver (Deiterding Reference Deiterding2003), while the reaction source term is numerically integrated through a semi-implicit fourth-order general Runge–Kutta stiff ordinary differential equation solver (GRK4A) (Kaps & Rentrop Reference Kaps and Rentrop1979). A sufficiently small grid size is necessary to capture the fine structure of the detonation wavefront. In this study, adaptive mesh refinement technology, based on the AMROC framework (Deiterding Reference Deiterding2009, Reference Deiterding2011a ,Reference Deiterding b ), is employed to reduce both computing time and memory requirements.

Figure 6.

Numerical smoked foils with four different grid resolutions: (a) 8 pts/ $\varDelta _{\textit{I}}$ , (b) 16 pts/ $\varDelta _{\textit{I}}$ , (c) 32 pts/ $\varDelta _{\textit{I}}$ , (d) 48 pts/ $\varDelta _{\textit{I}}$ .

Figure 7.

Temperature contours of the detonation front in the non-uniform gas mixture with different grid resolutions: (a) 8 pts/ $\varDelta _{\textit{I}}$ , (b) 16 pts/ $\varDelta _{\textit{I}}$ , (c) 32 pts/ $\varDelta _{\textit{I}}$ , (d) 48 pts/ $\varDelta _{\textit{I}}$ .

4.3. Grid convergence test

A grid convergence test for detonation propagation in a non-uniform gas with $\Delta T = 55\ \text{K}$ , as shown in figure 1, using four different grid sizes – 8 pts/ $\varDelta _{\textit{I}}$ , 16 pts/ $\varDelta _{\textit{I}}$ , 32 pts/ $\varDelta _{\textit{I}}$ and 48 pts/ $\varDelta _{\textit{I}}$ – was carried out to determine the adequate grid size. The resulting numerical smoked foils and flow fields at a specific moment ( $t = 200$ $\unicode{x03BC}$ s) are illustrated in figures 6 and 7, respectively. A comparison of these flow fields revealed that grid resolutions 32 pts/ $\varDelta _{\textit{I}}$ and 48 pts/ $\varDelta _{\textit{I}}$ exhibited more intricate details compared to 8 pts/ $\varDelta _{\textit{I}}$ and 16 pts/ $\varDelta _{\textit{I}}$ . Notably, the flow fields from the 32 pts/ $\varDelta _{\textit{I}}$ and 48 pts/ $\varDelta _{\textit{I}}$ simulations demonstrated considerable similarity, particularly regarding the fluid structures. In addition, the maximum pressure history curves along $y = 0.015$ m from figure 6 are plotted in figure 8 for comparison. It was observed that the maximum pressure history curves for grid sizes 32 pts/ $\varDelta _{\textit{I}}$ and 48 pts/ $\varDelta _{\textit{I}}$ almost completely overlapped. These results indicate that a grid size of 32 pts/ $\varDelta _{\textit{I}}$ provides an accurate flow field resolution, justifying its use in the present study.

Figure 8.

Maximum pressure curves along $y=0.015$ m with different grid resolutions.

Figure 9.

Numerical smoked foils of shocks in non-uniform gas with different temperature amplitudes: (a) $\varDelta T=35$  K and (b) $\Delta T=105$ K. Line AB is the temperature interface; line AC is the trajectory of the first triple-point; line AD is the trajectory of the second triple-point; line CD is the trajectory of the first triple-point after Mach reflection.

5.1. Shock propagation in non-uniform mixtures

Figure 9 shows numerical smoked foils of detonation propagation in the non-uniform gas for a case with low-temperature amplitude ( $\Delta T=35$ K) and a case with high-temperature amplitude ( $\Delta T=105$ K). The corresponding transient structures of the shock front are demonstrated in figure 10. The observed difference in the trajectories’ inclinations between the two smoked foils in figure 9 indicates that the curvature of the wavefront varies with the temperature amplitude of the non-uniformity. The strength of OS1 can be assessed by comparing the inclination angles of line AC. As shown in figure 9, the inclination angles for the two cases are $22.42^\circ$ and $28.15^\circ$ , corresponding to temperature amplitudes 55 and 105 K, respectively. Therefore, OS1 is weaker in the case with temperature amplitude 55 K compared to the case with temperature amplitude 105 K.

Figure 10.

Initial transition process of shocks in the non-uniform mixtures: (a) $\Delta T=35$ K, 7.5  $\unicode{x03BC}$ s; (b) $\Delta T=35$ K, 45  $\unicode{x03BC}$ s; (c) $\Delta T = 105$ K, 7.5 ${\unicode{x03BC} }$ s; (d) $\Delta T = 105$ K, 40 ${\unicode{x03BC} }$ s.

Figure 11.

Temperature variations along the upper wall and central lines in the channel: (a) $\Delta T=35$ K, wall temperature; (b) $\Delta T=35$ K, midline temperature; (c) $\Delta T=105$ K, wall temperature; and (d) $\Delta T=105$ K, midline temperature.

In accordance with the extended GSD results in § 3, two oblique shocks, namely OS1 and OS2, are initially formed in the cold zone due to the presence of the temperature interface, as shown in figure 10. Additionally, lines AC and AD represent the trajectories of the two triple-points at the shock front. The analysis further reveals that the evolution of the shock front involves the reflection of the first oblique shock OS1 on the top wall, resulting in the formation of a Mach stem and a new triple-point trajectory line CD. Subsequently, the collision and merging of the two triple-points at point D lead to the development of a single triple-point, leaving a trajectory line DE. The horizontal trajectory DE signifies the asymptotic approaching of a steady state by the shock front following the initial transit process. An examination of the small-temperature amplitude scenario reveals a correspondingly small angle $\angle$ BAD, indicating the weakness of the oblique shock OS2 adjacent to the interface. Consequently, the weak oblique shock OS2 effortlessly merges when encountering the reflected triple-point, forming a weak von Neumann Mach reflection structure characterised by a smooth connection between the Mach stem and the oblique shock OS2. Conversely, under conditions with large-temperature amplitude, the interaction between the strong oblique shock OS2 and the Mach stem (shock–shock interaction) deflects the trajectory AD to DE. This culminates in the eventual establishment of a stable state of coexistence featuring the single Mach reflection adjacent to the oblique shock OS2. These two steady-state structures are consistent with the strong and weak reflections of the triple-wave structure of the detonation front proposed by Strehlow (Reference Strehlow1964), henceforth referred to as weak and strong structures in the present study.

Figure 11 shows the temperature curves at different locations along the wall and midline of the channel. In figure 11(c,d) for the case with temperature amplitude $\Delta T = 105$ K, upon entry into the non-uniform gas, it is observed that the temperature of the shock front decreases to 1750 K in the cold layer, and increases to 2700 K in the hot layer. The formation of a Mach stem on the wall results in an increase in the temperature from 1750 K to 1950 K. This is followed by a transition process where the wall temperature increases from 1950 K to 2300 K. Eventually, a temperature 2300 K is reached, slightly exceeding the post-shock temperature 2241 K in the uniform gas. Subsequently, the temperature and pressure decay due to the rarefaction wave from the interface in the hot layer along the midline. In the case with $\Delta T=35$ K, as shown in figure 11(a,b), the transition process is observed to be more gentle, resulting in smaller changes in temperature and pressure compared to the case with $\Delta T = 105$ K.

Figure 12.

Steady structures of shock waves: (a) schlieren plot with streamlines, $\Delta T = 35$ K; (b) schlieren plot with streamlines, $\Delta T = 105$ K; (c) temperature, $\Delta T = 35$ K; (d) temperature, $\Delta T = 105$ K; (e) pressure and temperature curves, $\Delta T = 35$ K; ( f) pressure and temperature curves, $\Delta T = 105$ K.

Figure 12(a,b) illustrate the steady-state configuration of shock waves for two distinct temperature amplitudes, $\Delta T = 35$  K and $\Delta T = 105$ K, representing weak and strong shock structures. In the case of a weak structure, the shock front in the hot layer exhibits slight convexity, while it displays a slight concavity in the cold layer. A smooth and continuous transition links the convex front in the hot layer to the concave front in the cold layer. Conversely, the strong structure requires an oblique shock to connect the front in the cold layer with the front in the hot layer, giving rise to a Mach reflection structure. The magnitude of the oblique shock OS2 is heavily contingent on the temperature amplitude across the horizontal interface. In the strong structure, the normal velocity of the fluid particle passing through OS2 is significantly lower than that passing through the Mach stem. Consequently, the temperature behind OS2 is considerably lower than behind the Mach stem, and much lower than behind the hot layer front, as observed in figure 12(c,d). Behind OS2, the temperature interface transforms into a contact discontinuity that couples with the shear layer of the Mach reflection structure, forming a new cold zone that segregates the gases in the cold and hot layers. Figure 12(a,b) also present the streamline in the frame of the leading shock. The flow field is partitioned into five regions in the $y$ -direction by the front in the hot layer, two oblique shocks, and two Mach stems. The streamlines behind the front in the hot layer and oblique shocks deflect towards the cold layer, affecting the continuous compression of the gas within the cold layers. The temperature and pressure profiles along the designated dashed line in the temperature contours are shown in figures 12(e) and 12( f). The temperature of the cold zone is observed to be 1000–1200 K lower than that of the gas on either side of the interface.

The transition from a weak to a strong structure occurs as the temperature amplitude increases from a small value (35 K) to a large value (105 K). A distinctive feature between the two structures is the presence of oblique shocks linked to the interface. Figure 13 shows the variations in oblique shock height at different temperature amplitudes. We also compare the variation in the post-shock temperature and pressure behind the oblique shock with temperature amplitude. When the temperature amplitude is below a specific value, approximately 40 K, discernible oblique shocks are absent, resulting in a weak structure. In contrast, a higher temperature amplitude increases the oblique shock height for the strong structure. A longer oblique shock indicates a greater difference in position between the hot and cold zones, leading to a wider band of unreacted gas. Therefore, with a larger temperature difference, the oblique shock has a more significant effect on shaping the structure.

Figure 13.

Heights and post-shock states of (a) oblique shock (OS) and (b) Mach stem (MS) for cases with different temperature amplitudes $\Delta T$ .

The study of the inert shock provides insights into the propagation of ZND detonation in non-uniform gases. A key distinction between ZND detonation and the inert shock lies in the induction and chemical reaction. The rates of induction and exothermic reactions are influenced by the gas state behind the leading shock. The release of chemical energy significantly affects the leading shock in turn. Notably, the changes in fluid state have a profound effect on the induction and exothermic processes for detonations. The relationship between density, temperature, and the induction and exothermic rates is crucial for the two-step reaction model used in the present study. However, due to the exponential form of the Arrhenius law, the temperature gradient tends to dominate the induction and exothermic rates. In the hot layer, the high temperature behind the leading shock leads to a shorter induction zone, causing the front to be overdriven. Conversely, the low temperature behind the leading shock results in a lengthened induction zone in the cold layer, causing the front to be underdriven. It is important to note that the fronts in the two layers are not isolated from each other. They undergo a transition process before reaching a relatively stable state.

Figure 14.

Numerical smoked foil of a ZND detonation propagating into non-uniform mixtures: (a) $\Delta T=35$  K, (b) $\Delta T=55$ K, (c) $\Delta T=75$  K and (d) $\Delta T=105$ K. Two blue shadows and one red shadow are added to the numerical smoked foils to indicate the range of non-uniformity.

5.2. The ZND detonation propagating in non-uniform mixtures

In this subsection, we examine the propagation of a ZND detonation in the non-uniform mixture to investigate the interaction between chemical reaction and non-uniformity on the front structure and its evolution. As shown in figure 14, the numerical smoked foils for cases with different temperature amplitudes $\Delta T$ depict the presence of non-uniformity with two cold layers and one hot layer. The trajectories left by the triple-points moving between two channel walls are observable on the smoked foil. These triple-points originate from the double triple-points structure on the shock front in the cold layers, consistent with the initial transit process of the inert shock, as discussed in the previous subsection. The corresponding velocity profiles of the ZND detonation in non-uniform mixtures with different temperature amplitudes are shown in figure 15. The velocity profile is calculated by measuring the difference in wavefront position over a time interval. To capture the dynamics of high-frequency oscillations from transverse waves, data frames are collected every 5  $\unicode{x03BC}$ s. To accurately analyse the average velocity of the detonation wave, local fluctuations from these waves are filtered out, and consecutive frames are sampled at 50  $\unicode{x03BC}$ s intervals.

In figure 14(a), we can still observe that the horizontal triple-point trajectories generated by the oblique shock are consistent with cases of inert shock. However, in the present case with chemical reaction, the transverse waves severely interfere with the triple-point of the oblique shock associated with the temperature interface, causing its trajectory to be discrete and oscillated. The transverse waves dominate the propagation, and the temperature interface can only disturb the cell structure due to slight temperature differences across the interface. After the position where $x = 0.9$ m, as shown in figure 14(a), the two triple-points no longer propagate symmetrically, and the diamond-shaped trajectories begin to twist. The loss of asymmetry of the triple-point trajectories can be attributed to the intrinsic instability of the detonation overcoming the initial regularity from the symmetry of the temperature interface. To verify this conclusion, we conducted a new simulation with $k_R = 0.9\, {\textrm{s}}^{-1}$ , which corresponds to a very stable mixture, in contrast to the weakly unstable mixture used in our current study ( $k_R = 2.0\, \text{s}^{-1}$ ). The numerical smoked foils are given in figure 16 for comparison. When comparing these two mixtures, we found that symmetry is maintained over a long distance (1.6 m) for the very stable mixture. In contrast, the weakly unstable mixture loses symmetry after just 1.2 m, highlighting the impact of intrinsic instability.

Figure 15.

Velocity profiles of a ZND detonation in non-uniform mixtures with (a) $\Delta T=35$  K, (b) $\Delta T=55$  K, (c) $\Delta T=75$  K and (d) $\Delta T=105$ K.

Figure 16.

Numerical smoked foil of a ZND detonation propagating into non-uniform mixtures with $k_{\textit{R}} = 0.9 \text{ s}^{-1}$ and $\Delta T=55$ K.

The numerical smoked foil shown in figure 14(b) exhibits characteristic features a case with $\Delta T = 55$ K. Notably, ‘diamond’ cells manifest along triple-point trajectories. When the ZND detonation penetrates the non-uniform area, its front evolution is similar to that of a shock wave due to the velocity disparity on the front, resulting in the appearance of triple-wave structures (denoted as B, E, H, L, O, etc.) within the cold layer. As the triple-point of the expanding reactive Mach stem traverses the interface into the hot layer (D, G, N, etc.) and intersects with the symmetrical side’s triple-point, a ‘diamond’ pattern emerges on the numerical smoked foil. At the position where $x=1.2$ m, the cell pattern begins to distort, gradually losing its symmetrical ‘diamond’ shape. The trajectory pattern closely resembles that depicted in figure 14(a). In this case, the temperature amplitude surpasses that in figure 14(a), consequently intensifying the oblique shock wave and significantly disrupting the cellular structure. The theoretical analysis and numerical results regarding the shock have demonstrated that as the amplitude of temperature increases, the strength of the oblique shock correspondingly intensifies. When the strength of the oblique shock is relatively weak compared to that of the transverse wave, the transverse wave is less likely to be influenced by the oblique shock. Conversely, the transverse waves can significantly affect the oblique shock, which may effectively ‘swallow’ the oblique shock and modify the local structure of the detonation front, particularly as shown in the case with the smallest temperature amplitude ( $\Delta T=35$ K). In cases where the strength of the oblique shock exceeds that of the transverse wave, the oblique shock exhibits considerable resistance to disturbance from the sweeping transverse waves. In contrast, the presence of a strong oblique shock effectively obstructs the transverse waves, as evidenced in the scenario with the largest temperature amplitude ( $\Delta T=105$ K).

Figure 17.

Evolution of a ZND detonation wave immediately after entering the non-uniform region ( $\Delta T=55$  K): (a) 0  $\unicode{x03BC}$ s, (b) 15  $\unicode{x03BC}$ s, (c) 20  $\unicode{x03BC}$ s, (d) 30  $\unicode{x03BC}$ s, (e) 35  $\unicode{x03BC}$ s, ( f) 40  $\unicode{x03BC}$ s, (g) 45  $\unicode{x03BC}$ s and (h) 55  $\unicode{x03BC}$ s.

Figure 18.

Profiles of the temperature and concentration of the radical ( $\Delta T=55$ K): (a) wall and (b) midline.

Figures 17 and 18 illustrate the formation of diamond-shaped trajectories, as shown in figure 14(b). The initial stage, depicted in figure 17(b), demonstrates significant front bending, presenting a convex front in the hot layer and two oblique shocks denoted as OS1 and OS2 in the cold layer. Notably, the induction zone lengthens as it approaches the temperature interface in the hot zone due to the pronounced influence of rarefaction waves near the interface. The flow separation in the cold and hot layers is distinguished by a new strip area of unreacted gas behind OS2, confined by the temperature interface and the reflected shock associated with the Mach stem. The flow particles behind OS2 initially move towards the upper wall, then deflect back towards the hot layer due to the gas compression behind OS1. In figure 17(c), an inert Mach stem, denoted as MS1 formed by OS1, elongates, accelerating the reaction and causing a significant bulging of the reaction front. Subsequently, in figure 17(d), the reaction front becomes fully coupled to the inert Mach stem, forming a strong reactive Mach stem denoted as RMS. Profiles of temperature and the concentration of intermediate products can be observed along the top wall and midline of the channel in figure 18(a). The profiles at 15–25 ${\unicode{x03BC} }$ s reveal that the Mach stem in the cold layer strongly accelerates the reaction, thereby decreasing the induction length.

It is worth noting that as illustrated in figure 17(d–f), the expansion of RMS into the hot layer leads to the shedding of the unreacted gas strip due to the disappearance of the oblique shock OS2. Concurrently, the detachment of the detonation front in the cold layer is initiated, causing an increasing distance between the leading shock and the reaction surface. Moreover, due to the high temperature in the hot layer, the rapid reaction of the two RMS entering the hot layer accelerates the transverse waves, resulting in the propagation towards the midline of the channel and colliding to form a new strong Mach stem expanding outwards. In figure 18(b), it can be observed that the induction zone of the front in the hot layer continues to lengthen. Additionally, the gradual decrease in temperature behind the leading shock in the hot layer indicates a continuous attenuation of the front intensity, persisting until the new strong Mach stem is formed in the hot layer (45 ${\unicode{x03BC} }$ s). The above discussion indicates that a competitive relationship exists between the detonation’s transverse waves and the oblique shock associated with the temperature interface.

Figure 14(c) shows the smoked foil for the case with temperature amplitude $\Delta T = 75$  K. Before $x=0.4$ m, discernible transverse waves are not observed entering the hot layer. This is because the chemical reaction behind the Mach stem in the cold layer is inhibited due to lower initial temperature, resulting in the Mach stem not being strong enough to extend into the hot layer. After $x=0.4$ m, the triple-point trajectories gradually display cellular patterns, though the trajectories in the hot layer remain relatively weak. Compared to the cellular patterns depicted in figures 14(a) and 14(b), the cell length in figure 14(c) is more extended, indicating a weaker transverse wave in the hot layer with a lower velocity. Upon increasing the non-uniformity to $\Delta T = 105$ K, as illustrated in figure 14(d), no prominent triple-point trajectories (or cells) manifest in the channel, except for faint traces generated by extremely weak transverse shocks. The smoked foil resembles the case of inert shock, with the notable distinction of the trajectory of the triple-point attached to the oblique shock now exhibiting fluctuation. This observation suggests that the oblique shock associated with the interface is the dominant factor, with the transverse wave being highly suppressed by the substantial temperature non-uniformity. At low temperatures, the chemical reaction rate in the cold layer slows down significantly, causing the front to separate, and restricting the generation of transverse waves. During this phase, the movement of the triple-point connected to the oblique shock is only slightly affected. However, the explosion of the hotspot in the shocked zone behind the separated shock front in the cold layer causes the front to speed up and push the oblique shock back to the interface. The expansion of the front in the cold layer is not sustainable, leading to a new cycle of separation in the cold layer. In the stage where $0.4\lt x\lt1.2$ m, this process repeatedly occurs during propagation, resulting in highly disturbed triple-point trajectories. When $x\gt 1.2$ m, high-frequency weak transverse waves in the cold layer stabilise the reaction rate, causing the disturbance of the triple-point trajectories to gradually diminish and reach a steady state.

As shown in figure 14(d), the extremely low temperature in the cold layer (190 K) compared to the other cases significantly suppresses the chemical reaction behind the Mach stem. This suppression prevents the expansion of the Mach stem and its merging with the oblique shock that links to the temperature interface, thus prohibiting its entry into the hot layer. Although a higher initial temperature can promote the reaction, the shock front of the hot layer decays due to the curvature effect, which is consistent with the detonation propagation subject to boundary loss, as described by Lee (Reference Lee2008).

In the case where $\Delta T=105$ K, the propagation of the detonation wave can be analysed in four distinct stages. The first stage occurs approximately between $x= 0.05$  m and $x= 0.17$ m, during which the ZND detonation enters the non-uniform region, leading to the decoupling of the front in the cold layer. Subsequently, the weakened front in the hot layer results in the elongation of the induction zone. If no stimulation follows the leading shock, then the detonation in the hot layer will fail completely. However, as depicted in figure 19, the Mach reflection of the decoupled shock in the cold layer produces a reflected shock that strengthens the decoupled front in the hot layer. The discussion in § 2 suggests that the transmitted shock accelerates within the hot gas mixture as it passes through the temperature interface. Consequently, a reflected wave is formed behind the oblique shock OS2 in the cold layer, transmitting into the hot detonation products across the temperature interface, as demonstrated in figure 19(a,b). In the hot product, the reflected wave experiences acceleration, separating from the initial reflected wave and merging with the triple-point, and advances toward the induction zone located behind the leading front. Figure 19(c–e) illustrate that two reflected shocks within the induction zone converge, intensifying local pressure and temperature. As a result, a region of heightened temperature, known as a hotspot, emerges within the disengaged induction zone, thereby facilitating the detonation process.

The pressure profiles along the midline of the channel are shown in figure 20(a), demonstrating the compression and re-initiation process of the two reflected shocks. The reflected waves increase pressure from 0.75 MPa to 1.0 MPa, which generates a local hotspot at $x=0.145$ m, triggering a local hotspot explosion. The compression wave produced by the hotspot explosion catches up with and strengthens the leading shock, which can be observed in figure 20(b). Furthermore, the expanding gas and transverse waves stemming from the hotspot explosion propagate into the cold region, causing gas compression. However, at $x=0.15$ m in the hot zone, the hotspot explosion fails to initiate a chemical reaction in the cold zone.

Figure 19.

Stimulation of cold layer transverse waves to the hot layer front ( $\Delta T=105$ K): (a) 20  $\unicode{x03BC}$ s, (b) 25  $\unicode{x03BC}$ s, (c) 35  $\unicode{x03BC}$ s, (d) 40  $\unicode{x03BC}$ s, (e) 45  $\unicode{x03BC}$ s, ( f) 55  $\unicode{x03BC}$ s, (g) 65  $\unicode{x03BC}$ s and (h) 70  $\unicode{x03BC}$ s.

Figure 20.

Formation of a hotspot under transverse waves compression ( $\Delta T=105$ K): (a) pressure and (b) $x{-}t$ diagram.

The second stage, as shown in figure 21(a–d), occurring approximately between $x=0.17$ m and $x=0.33$ m, is characterised by the complete suppression of chemical reaction in the cold layer and significant suppression in the hot layer. The detonation is severely underdriven during this stage. The trajectory of the triple-point on the numerical smoked foil exhibits a horizontal line. In the front velocity profiles depicted in figure 15(d), the average front velocity is merely $0.8V_{\textit{CJ}}$ , and the oscillations of the front velocity curves in the two layers are minor, indicating a relatively stable state of the detonation during this stage. The chemical reaction provides weak support for the front in the hot layer, resulting in a structure similar to the strong, steady-state structure of the inert shock wave. The second stage is short-lived, as compression waves enter the induction zone in the hot zone, eventually triggering new hotspots, as shown in figure 21(a–d).

Figure 21.

Sequences of front structures for the case with temperature amplitude $\Delta T=105$ K: (a) 105  $\unicode{x03BC}$ s, (b) 120  $\unicode{x03BC}$ s, (c) 135  $\unicode{x03BC}$ s, (d) 150  $\unicode{x03BC}$ s, (e) 435  $\unicode{x03BC}$ s, ( f) 450  $\unicode{x03BC}$ s, (g) 465  $\unicode{x03BC}$ s, (h) 475  $\unicode{x03BC}$ s, (i) 750  $\unicode{x03BC}$ s, (j) 760  $\unicode{x03BC}$ s, (k) 770  $\unicode{x03BC}$ s and (l) 780  $\unicode{x03BC}$ s.

In the third stage, the chemical reaction in the cold layer is triggered by the compression effect generated by the hotspot in the hot layer. The gas behind the leading shock in the cold layer starts to react close to the position where $x=0.36$ m, as shown in figure 14(d) ( $\Delta T=105$ K). The typical evolution of the detonation in this stage is depicted in figure 21(e–h). As shown in figure 21(e), the hot layer front outpaces the cold layer front, resulting in the appearance of a double triple-point structure in the cold layer. Moving on to figure 21( f), the reactive Mach stem chases the hot layer front, and it is noteworthy that two reflected waves converge towards the hot layer front. The chemical reaction accelerates and elongates the Mach stem, driving the triple-point towards the midline of the channel, and causing the trajectory on the numerical smoked foil to deflect towards the midline. However, the low temperature in the cold layer weakens the front before the triple-point enters the hot layer, and as a result, the trajectory deflects towards the walls. As shown in figure 21(g), the triple-point comes close to the interface but does not cross it. Subsequently, in figure 21(h), the reflected wave of the triple-wave structure once again strengthens the hot layer front, expanding and creating fluctuating trajectories on the numerical smoked foil. At this stage, the velocity curves of the front in the two layers experience significant oscillations and exhibit a phase difference between them.

The alternating strengthening of hot and cold layers characterises the behaviour of detonation waves during the third stage. Once the waves have strengthened, the entire detonation wave returns to a steady state, similar to the strong structure of an inert shock. This steady state also marks the beginning of the fourth stage. As shown in figure 14(d) ( $\Delta T=105$ K), the numerical smoked foil after $x=1.3$ m is partitioned into three regions by the trajectories, resembling cases of inert shocks shown in figure 9(b). This similarity in structure suggests a resemblance between the detonation front and the shock front. Figures 21(i–l) illustrate the detonation fronts at four different moments, presenting a steady-state structure comprising a convex front in the hot layer, two oblique detonations associated with the temperature interface, and two planar fronts in the cold layers. At these moments, the higher temperatures in the hot layer and the weak transverse waves in the cold layers sustain the entire wavefront into a steady state. Further analysis, as shown in figure 14(d), reveals the presence of a local single-headed mode with weak transverse waves in the two thin cold layers near the walls. The coupling of the detonation front in the hot layer is primarily attributed to the curvature not being large enough to attenuate the detonation, and to a lesser extent, the presence of weak and disorganised transverse waves. The discrete nature of the entire detonation front, divided into three independent segments by temperature discontinuities, suggests that it is primarily the oblique shock rather than the transverse waves that governs the propagation mechanism in the non-uniform gas.

The plots in figure 15 present the velocity profiles of the detonation, offering insights for comprehending the underlying physics. Without non-uniformity, the ideal ZND detonation wave propagates at a constant CJ speed despite local fluctuation. However, with the presence of non-uniformity, the wavefront velocities exhibit distinct behaviour. Initially, the front velocities increase in the hot layer and decrease in the cold layer, consistent with prior discussions. Following a brief acceleration, the front velocity appears to decrease in both the cold and hot layers. Moreover, a higher temperature amplitude leads to a greater loss in velocity. A re-initiation process is then observed with the velocity exceeding $V_{\textit{CJ}}$ slightly, indicating that a steady state is reached after the transient phase. For scenarios with the largest temperature non-uniformity ( $\Delta T=105$ K), the amplitude and frequency of velocity oscillation are notably reduced due to the suppression of transverse waves, with oblique shocks associated with the temperature interface dominating the propagation. Once the detonation attains a relatively steady state, the average detonation speed exceeds $V_{\textit{CJ}}$ to a greater extent as the temperature amplitude ( $\Delta T$ ) increases. This is attributed to the slight change in the total energy of the gas with $\Delta T$ . Notably, the detonation speed exhibits a nonlinear relationship with the initial temperature due to its dependency on total energy (or enthalpy), which differs for non-uniform gases at the same averaged temperature. Specifically, under constant initial pressure, the enthalpy of a non-uniform gas slightly increases with rising $\Delta T$ compared to that of a uniform gas.

Figure 22 illustrates the detonation fronts at $x=1.4$ m for four different $\Delta T$ . In figure 22(a), the front structure exhibits symmetry around the midline of the channel, indicating minimal influence of the temperature interface on the detonation propagation. The transverse waves are observed to propagate easily through the temperature interfaces. However, in figures 22(b) and 22(c), the front loses its symmetry around the midline, suggesting that the influence of non-uniformity is involved. Additionally, in figure 22(d), the front displays division due to the two temperature interfaces comprising cold and hot layers fronts along with two oblique shocks. This observation suggests that the non-uniformity completely prevents the transverse waves from crossing the interface at this stage.

Figure 22.

Detonation fronts near $x=1.4$ m for cases with different temperature amplitudes: (a) $\Delta T=35$  K, (b) $\Delta T=55$ K, (c) $\Delta T=75$  K and (d) $\Delta T=105$ K.

In the present study, we used the reactive Euler equations with a simple two-step reaction model, ignoring some diffusion effects. To verify the approximations of the chemical reaction model, we conducted a new simulation with a detailed chemical reaction mechanism of a mixture $2{\textrm{H}}_2+{\textrm{O}}_2+7{\textrm{Ar}}$ with initial pressure 6670 kPa and 295 K, which is close to the case used in the present study. The results in figure 23 show that the difference in the chemical reaction model does not generate a considerable variation in the smoked foils and detonation front with the presence of the non-uniformity of temperature.

Figure 23.

Numerical smoked foils of cellular detonations in a mixture $2{\textrm{H}}_2+{\textrm{O}}_2+7{\textrm{Ar}}$ with different non-uniformities: (a) $\Delta T=55$  K and (b) $\Delta T=105$ K.

Figure 24.

Numerical smoked foils for the detonation propagation with a temperature non-uniformity ( $\Delta T=15$ K) of different wavelengths: (a) $s = 0.01$ m, (b) $s = 0.02$ m, (c) $s = 0.04$ m and (d) $s = 0.12$ m.

To investigate how the length scale of temperature non-uniformity affects detonation propagation, we conducted simulations to examine the relationship between the wavelength of non-uniformity, cell size and propagation velocity. Figure 24 displays the numerical smoked foils with varying wavelengths of non-uniformity, all maintaining the same temperature amplitude 15 K and a reaction pre-exponential factor $k_R = 2.0\text{ s}^{-1}$ . Figures 25 and 26 illustrate the changes in cell size and velocity along the direction of propagation, respectively.

Figure 25.

Cell sizes for the detonation propagation with a temperature non-uniformity ( $\Delta T=15$ K) of different wavelengths.

Figure 26.

Velocity profiles for the detonation propagation with a temperature non-uniformity ( $\Delta T=15$ K) of different wavelengths: (a) $s = 0.01$ m, (b) $s = 0.02$ m, (c) $s = 0.04$ m and (d) $s = 0.12$ m.

In figure 24(a), a temperature non-uniform region with wavelength 0.01 m ( $0.25\lambda$ ) at $x = 0.06$ m causes fine cells to form at $x = 0.1$ m, matching the wavelength. Due to the instability of the detonation wave, cell length increases at $x = 0.05$ m and stabilises near $x = 0.8$ m. In figure 24(b), increasing the wavelength to 0.02 m ( $0.5\lambda$ ) delays the enlargement of cell length, which begins at $x = 0.8$ m and stabilises at $x = 1.0$ m, indicating that broader wavelengths suppress transverse wave development. In figure 24(c), when the wavelength reaches 0.04 m ( $1.0\lambda$ ), the cell size synchronises with the inhomogeneity scale, stabilising at $x = 0.1$ m. This highlights a strong stabilising effect of temperature discontinuities at wavelengths close to the CJ cell width. In figure 24(d), with wavelength 0.12 m ( $3\lambda$ ), oversized cells form due to transverse waves, creating a ‘double cellular structure’ noted in Chen et al. (Reference Chen, Huang, Lei, Wang, Chen and Han2023). This ‘double cellular structure’ exhibits oscillatory characteristics below $x = 0.8$ m, and remains unstable until the waves merge beyond $x = 1.2$ m, after which cell length stabilises.

The velocity profiles demonstrate that the wavelength affects both the frequency and amplitude of oscillations in the detonation wave velocity during the initial propagation stage. As shown in figure 26, larger wavelengths lead to velocity oscillations characterised by lower frequencies and higher amplitudes. The reduced frequency is linked to greater spacing between transverse waves, while the higher amplitude indicates a stronger intensity of these waves. Although a larger wavelength produces fewer transverse waves, the present waves are more powerful. When the wavelength exceeds the cell width (as illustrated in figure 26), it suppresses the development of intrinsic detonation instability, causing irregular velocity oscillations to emerge later. In figure 26, the velocity signature of the ‘double cellular structure’ (highlighted within the red dashed boxes) reveals that large-scale fluctuations indicate oversized cells. In contrast, the small-scale oscillations inside the red box are attributable to the substructure cells, which emerge from the intrinsic instability of the mixture.

Next we compare the effects of different temperature amplitudes $\Delta T$ on detonation wave evolution with fixed wavelength $s = 0.02$ m. Temperature perturbations primarily influence reaction decoupling in cold zones, subsequently triggering secondary effects.

Figure 27.

Numerical smoked foils for a ZND detonation propagation with temperature non-uniformity of different amplitudes: (a) $\Delta T=15$ K, (b) $\Delta T=55$  K and (c) $\Delta T=105$ K.

In the numerical smoked foil, as shown in figure 27(a) ( $\Delta T = 15$ K), no reaction decoupling occurs in cold zones, and transverse waves develop freely. Cellular structures within $x \lt 0.8$ m are induced by temperature discontinuities, with cell widths matching perturbation wavelengths. Beyond $x \gt 0.8$ m, cells arise from detonation instabilities. Trajectories remain smooth, minimally affected by temperature discontinuities, and cell sizes approach those of unperturbed cases. In figure 27(b) with $\Delta T = 55$ K, when $x \lt 0.4$ m, reaction decoupling in cold zones suppresses transverse wave propagation, significantly increasing cell size. When $0.4\lt x \lt0.8$ m, chemical reactions recover in cold zones. Transverse wave velocities increase, reducing cell size. When $x \gt 0.8$ m, gas instabilities enlarge transverse wave spacing, increasing cell length. Temperature discontinuities obstruct transverse waves, distorting triple-point trajectories. The cell sizes are smaller than the case with $\Delta T = 15$ K. In figure 27(c) with $\Delta T = 105$ K, chemical reactions in cold zones are strongly suppressed. Oblique shock fronts and low-temperature gas bands at temperature discontinuities severely hinder transverse wave propagation. No typical cellular patterns form on the soot foil. Near $x = 0.6$ m, transient transverse waves emerge due to gas instabilities but rapidly dissipate, demonstrating the strong obstruction effect of large temperature amplitude.

As shown in figure 28, the cell size decreases as the temperature amplitude increases. The cells vanish when the temperature amplitude becomes particularly large, and the cell size can be zero. However, at this point, the detonation wave does not fail. Instead, it exists in a discrete wavefront state, where the detonation wave propagates within the individual hot and cold regions in the form of a planar detonation wave. Based on the velocity profiles in figure 29, larger temperature amplitudes cause more significant initial velocity deficits in detonation waves. Wave velocities gradually recover only after transverse waves emerge. Notably, in figure 29(c) ( $x \gt 1.0$ m), the velocity curve exhibits high-frequency oscillations due to localised transverse wave disturbances in cold zones. These disturbances occur because transverse waves on the wavefront cannot penetrate temperature discontinuities.

Figure 28.

Cell sizes for the ZND detonation propagation with temperature non-uniformity of different amplitudes.

Figure 29.

Velocity profiles for the ZND detonation propagation with temperature non-uniformity of different amplitudes: (a) $\Delta T=15$ K, (b) $\Delta T=55$  K and (c) $\Delta T=105$ K.

5.3. Propagation of cellular detonations in non-uniform mixtures

In the previous subsection, the analysis focused on the propagation of ZND detonations in non-uniform gases characterised by the presence of a temperature non-uniformity. This subsection shifts the focus towards examining the propagation of cellular detonations in similar non-uniform gases. A fully developed cellular detonation front is initially situated at $x=0.1$ m upstream of the non-uniform area started at $x=0.3$ m.

Figure 30.

Numerical smoked foils of cellular detonation propagation with different non-uniformities: (a) $\Delta T=35$ K, (b) $\Delta T=55$ K, (c) $\Delta T=75$  K and (d) $\Delta T=105$ K.

Figure 31.

Velocity profiles of cellular detonation propagation with different non-uniformities: (a) $\Delta T= 35$ K, (b) $\Delta T=55$ K, (c) $\Delta T=75$  K and (d) $\Delta T=105$ K.

Figure 32.

Formation of the new cells of the cellular detonation for the case with $\Delta T=55$ K.

In figure 30, numerical smoked foils of cellular detonations with four different non-uniformities are depicted. In the case with small $\Delta T$ , the propagation is primarily influenced by the transverse waves, with the oblique shock formed by the temperature interface having a limited effect. However, for cases with large $\Delta T$ , the oblique shock dominates the evolution of the front, with the high-frequency weak transverse waves causing only small fluctuations. Figure 31 illustrates the front velocities of the cellular detonations in non-uniform gases. A comparison with the front velocity of the ZND detonation depicted in figure 15 shows that the initial decay state disappears due to the presence of the cellular structures. However, in the case with $\Delta T=105$ K, the initial decay state persists, which implies that the cellular structures delay the effect of non-uniformity compared to the ZND detonation.

The regular and symmetrical non-uniformity is expected to create a symmetrical distribution along the detonation front. If the detonation is highly stable and characterised by regular transverse waves and regular detonation cells, then the presence of non-uniformity should not alter the symmetry of the flow field. However, the front loses its regularity and symmetry when dealing with unstable detonation. Interestingly, the presence of regular, spaced non-uniformities tends to enhance the regularity of unstable detonations. In this study, the width of the channel corresponds to the cell width. As a result, the regularly spaced temperature interfaces disturb the detonation front, leading to a symmetrical but locally distorted cell for an extended detonation travel. The type of detonation used in this study is weakly unstable, meaning that intrinsic instability becomes significant after a considerable distance as the initial regularity fades. If a highly stable detonation were employed, then it should maintain the symmetry of the cell even with the presence of regular, spaced non-uniformities.

Figure 33.

Numerical smoked foils of multi-headed cellular detonations with different non-uniformities: (a) $\Delta T=35$ K, (b) $\Delta T=55$ K, (c) $\Delta T=75$  K and (d) $\Delta T=105$ K.

Figure 32 shows the formation of new cells for cellular detonations in non-uniform gases. When the detonation reaches the non-uniform area at $x = 0.3$ m, an oblique shock wave (A and B) is observed on the wavefront. The reactive Mach stems merge into the oblique shock, leading to the simultaneous appearance of a kink (C and D) on the Mach stem. The reflection on the wall enhances the kink and transforms it into a stronger triple-point (E). Similarly, identical to the above discussion, a new oblique shock is formed (F and G), with H, I and J representing the new triple-points. The trajectories of the original and newly generated triple-points collectively give rise to the smaller cells. As mentioned earlier, when a shock wave interacts with a temperature interface, both a transmitted shock and a reflected shock (or a rarefaction wave) can be generated simultaneously. In the case of cellular detonations, the transverse shock can also transform into a transmitted shock and a reflected shock (or a rarefaction wave) upon encountering the interface. As such, the interface acts as a disturbance source that significantly impacts the detonation front and its propagation mechanism.

Furthermore, figure 33 shows the numerical smoked foils of the cellular detonation as the number of cells in the channel increases. The corresponding velocity profiles are plotted in figure 34. The findings reveal that the cell pattern remains relatively unchanged with only slight distortion with a neglectable weak non-uniformity, as shown in figure 33(a). However, significant changes in the cell pattern are observed in figure 33(b,c) with an increased non-uniformity. One change is the generation of fine cells, and the other is the high cell shape distortion, characterised by the deflection of the triple-point trajectories across the temperature interface. Following an extended transition process, a final steady state is reached at approximately $x=1.2$ m, featuring irregular cellular patterns and small subcells, which suggests that moderate non-uniformity may promote cellular instability. In the case with large non-uniformity (figure 33 d), the propagation of transverse waves is dramatically suppressed, and the detonation front appears isolated by the temperature interface. This results in dominating the oblique shock associated with the temperature interface. These observations reveal the interplay between the intrinsic instability of the detonation and the non-uniformity of flow, also highlighting the competition between the two factors.

Figure 34.

Velocity profiles of multi-headed cellular detonations with different non-uniformities: (a) $\Delta T=35$ K, (b) $\Delta T=55$ K, (c) $\Delta T=75$  K and (d) $\Delta T=105$ K.

Figure 35 shows the structures of the detonation front for various non-uniformities near $x=1.4$ m. The entire detonation front maintains regularity when $\Delta T=35$ K, as illustrated in figure 35(a). Notably, completed and extended transverse waves are observable behind the front, and the temperature interface has a negligible effect on the detonation front, merely causing a slight change in the induction length. The flat averaged velocity profiles depicted in figure 34(a) further confirm the regularity of the cellular pattern in the non-uniform area. Subsequently, in figure 35(b,c), with increasing non-uniformity of the flow field ( $\Delta T=55$ or 75 K), the influence of non-uniformity comes into play, with the emergence of smaller cells and trajectories twisting at the interface. Additionally, small-scale unreacted gas pockets are observed to develop behind the front, stemming from the oblique shocks connected to the interface. The crossing of transverse waves over the interface engulfs the oblique shocks, leading to the generation of unreacted gas pockets. The averaged velocity profiles shown in figure 34(b,c) exhibit low-frequency fluctuations, reaching amplitudes up to 10 % around $V_{\textit{CJ}}$ . Furthermore, in figure 35(d), the separation of the detonation front by the temperature interface results in the potent suppression of intrinsic instability due to the significant non-uniformity, leading to the disappearance of transverse waves. Consequently, a large mass of unreacted gas pockets forms behind the wavefront, and these are subsequently consumed slowly in the field of hot products, with no feedback on the front. As a result, the averaged propagation velocity of the detonation decreases by 10 % compared to the CJ speed. The pressure distribution along a vertical straight line near the detonation front is also plotted in figure 35. The results reveal that the significant temperature discontinuity enhances the distinctness and regularity of the detonation front by suppressing the development of transverse waves.

Figure 35.

Structures and pressure distributions along the detonation fronts for four different non-uniformities of temperature near $x=1.4$ m: (a) $\Delta T=35$ K, (b) $\Delta T=55$ K, (c) $\Delta T=75$  K and (d) $\Delta T=105$ K.