A Mixing Process Influenced by Wall Jet-Induced Shock Waves in Supersonic Flow

Featured Application

The present proposed mixing enhancement strategy can be a good candidate for application in Scramjet and RBCC combustors of hypersonic propulsion systems.

Abstract

With the development of hypersonic air-breathing propulsion systems, such as the supersonic combustion ramjet (Scramjet) and rocket-based combined cycle (RBCC) engines, the mixing process of supersonic airstream with fuel in the engine combustor has been drawing more and more attention. Due to the compressibility effects, the mixing process in a supersonic condition is significantly inhibited. In the present paper, the novel strategy of wall-jet induced shock waves (WJISW) is put forward to realize mixing enhancement. The interaction process between WJISW and the supersonic mixing layer is researched and the enhanced-mixing mechanism is revealed, employing large eddy simulation (LES) methods. The fine vortex structures of the flow field are well captured and presented, utilizing the numerical schlieren technique. Detailed visualization results indicate that WJISW in a low frequency condition can result in the ‘region action mode’ (RAM) never reported before. The drastic dynamic behaviors including growth, deformation, and distortion in the interaction region can undoubtedly promote the mixing of upper and lower streams. The Reynolds stress distributions along the streamwise x-direction suggest that more intense fluctuations can be achieved with a low frequency WJISW. Moreover, a sharp increase in mixing layer thickness can be realized in the interaction region. The dynamic mode decomposition (DMD) analysis results show that the mixing layer evolution process is dominated by the mode induced by WJISW, which leads to the coexistence of both large- and small-scale structures in the flow field. The entrainment process corresponding to large-scale vortices and the nibbling process corresponding to small-scale vortices can obviously promote mixing enhancement. It is suggested that the present proposed strategy is a good candidate for enhanced-mixing with application to Scramjet and RBCC combustors.

1. Introduction

In recent years, with the rise of research on rocket based combined cycle (RBCC) engines and scramjet, efficient mixing of air and fuel under high-speed compressible conditions has become a hot issue. For example, in the current popular RBCC combustion chamber, rich fuel and air mix and diffuse in the flow channel, and sufficient mixing is of great importance to improve engine combustor performance. Whereas, under the limited size of the combustor, the residence time of supersonic flow in the chamber is on the order of milliseconds [1], and this time provided for the mixing of fuel and air is extremely short. Apart from that, although in a low-speed condition the diffusion process of the mixing layer is very fast, the existence of a compressible effect significantly inhibits the growth process of the mixed layer under supersonic conditions. Therefore, considering the limited mixing time and the influence of a compression effect [2,3,4], it is urgent to adopt mixing enhancement strategies to improve mixing efficiency.

In order to realize mixing enhancement of supersonic flow, researchers worldwide have put forward a variety of mixing enhancement techniques, which can be divided into passive active methods. The passive mixing method perturbs the mixing layer by changing the configuration of the middle partition plate of the mixing layer or adding a disturbance device to speed up the mixing. Our previous work comprehensively reviewed passive strategies proposed by investigators to promote mixing [5]. Typical passive mixing enhancement measures include: lobe mixer [6], slope device [7], cavity [8,9], strut [10], and so on. A fact that cannot be ignored is that due to the introduction of a passive device, the inherent defect, including non-ignorable flow loss [11], limits the wide usage of these techniques. Active mixing enhancement methods are to t external energy into the mixing layer by introducing external excitation, and mixing enhancement is achieved by the regulation of the mixing layer. Typical active mixing enhancement strategies include acoustic excitation [12], plasma magnetic fluid excitation [13], blowing and suction [14], reflux device [15], vibration device [16,17], etc. The work conducted by Suzuki [18] et al. indicated that in the low-speed shear layer, sinusoidal excitation can promote earlier instability of the flow, and the rolled-up K–H structures can notably increase the mixing level. Whereas, in supersonic conditions, the effect of an electromagnetic flap is negligible.

Therefore, it can be concluded that the defect of both active and passive strategies calls for newer methods to realize mixing enhancement. It is worth noting that as the inherent flow structures in the confined supersonic mixing layer flow, shock waves [19] and their dynamic behaviors can significantly influence the flow evolution process. Specifically, the interaction between shock waves induced by a wedge shock generator and the mixing layer can notably affect the K–H vortex development. Menon [20] obtained the first experimental visualization results and found that shock wave induced by a wedge generator contributed to mixing enhancement. Lu and Wu [21] used direct numerical simulation (DNS) to research the effect of shock waves on the mixing layer and found that the influence of a shock wave on the mixing layer is very local and the installation should be set as early in the upstream as possible to short the initial mixing distance. Shau [22] claimed that the growth rate of the mixing layer increased slightly under the impingement of shock, but it returned to the original value within a short distance downstream. Gènin and Menon [23] develop a hybrid computation methodology to investigate turbulent flows containing shock waves, concluding that shock can promote high-level mixing in the non-similar region by amplifying the local turbulence, but this trend is limited in local space.

It can be found that a shock wave can achieve mixing enhancement, but with the introduction of a wedge generator, the promoted area is mainly near the interaction point. It is worth noting that as the most typical and simple means of fuel injection, a wall transverse jet and its flow field structure [24,25] and noise control [26] have been studied in detail by many scholars. Due to the jet’s obstruction of the incoming flow, a bow shock wave will be formed in the upstream area of the jet. In addition, because the position of the jet hole is controllable, the position of the shock wave generated by the jet is also controllable, which provides a good idea for us to use the controllable shock wave induced by the jet to achieve mixing enhancement. The previous research [27,28] related to wall jets is limited to the wall jet injecting fuel, whereas our research attempts to use the wall jet to induce a shock wave to act on the plane mixing layer. Based on this, we propose a new strategy called: wall-jet induced shock wave, which retains the advantages of shocks to enhance the growth rate of the downstream mixing layer and overcomes the disadvantage of the local range of shock interaction. Meanwhile, introducing shock waves through the wall-jet can conquer the drawbacks of the passive method and it does not require the participation of fixed devices in the flow field, which reduces flow loss. Furthermore, this scheme holds the advantage of an active method and the position of the shock wave interacting on the mixing layer can be controlled. By introducing shock wave induced through controlled wall jet to promote the mixing enhancement of supersonic mixing layers, it can expand the application of the classical problem of Shock wave/mixing layer interaction (SWMLI) and provide a new method to support the realization of efficient mixing under supersonic conditions.

In the present work, the large eddy simulation (LES) method is employed to analyze the effect of shock waves induced by a wall jet on the supersonic planar mixing layer, and the flow and growth properties of the mixing layer is studied by means of turbulence characteristics and evaluation of the thickness of the mixing layer. Besides, the modal decomposition is performed by the DMD method to further explore the deep-seated mechanism of mixing enhancement influenced by a wall jet. Through extensive reading of the literature, we found that DMD has been applied to the decomposition of the flow field where the moving shock wave [29] is located, which shows the feasibility of DMD. The remainder of this paper is organized as follows: Section 2 gives a brief introduction to the numerical method and the DMD process. Section 3 shows the simulation results in terms of flow field visualization and growth characteristics and analyzes the flow mechanism by mode decomposition. Finally, Section 4 gives the conclusions.

2. Methodology Details

2.1. Governing Equations and Numerical Methods

LES was conducted to solve the flow process by filtering the turbulent motion and decomposing the turbulent flow into large scale pulsations and small-scale pulsations. The large-scale pulsations are solved by direct numerical simulation (DNS) and the mass, momentum, and energy transport of the small-scale pulsations on the large-scale motion are modeled to close the equations, where the model is called a subgrid model. The Navier–Stokes equations, filtered by Deardroff box filtering [30], are as follows:

$$\frac{\partial \overline{\rho }}{\partial t}+\frac{\partial }{\partial x_i}\left(\overline{\rho }{\tilde{u}}_i\right)=0$$

(1)

$$\begin{array}{c}\frac{\partial }{\partial t}\left(\overline{\rho }{\tilde{u}}_i\right)+\frac{\partial }{\partial x_j}\left(\overline{\rho }{\tilde{u}}_i{\tilde{u}}_j+\overline{p}{\delta }_{ij}-{\tilde{\tau }}_{ij}+{\tau }_{ij}^{sgs}\right)=0\end{array}$$

(2)

$$\frac{\partial \overline{\rho }\tilde{E}}{\partial t}+\frac{\partial }{\partial x_i}\left[\left(\overline{\rho }\tilde{E}+\overline{p}\right){\tilde{u}}_i+{\overline{q}}_i-{\tilde{u}}_j{\tilde{\tau }}_{ij}+H_i^{sgs}+{\sigma }_i^{sgs}\right]=0$$

(3)

In the above equations, the mark ‘¯’ denotes spatial filtering and ‘~’ denotes Favre filtering. Superscript ‘sgs’ is variances associated with subgrid model. Here, ${\tau }_{ij}$ is the viscous stress tensor, $q_i$ is heat flux vector, and ${\delta }_{ij}$ is unit tensor. $E$ is energy, $H_i$ is the enthalpy flux vector, and ${\sigma }_i$ is the energy flux vector.

The kinetic-energy transport sub-grid model is adopted to measure the influence of small-scale transport on large-scale pulsation. The sub-grid viscosity coefficient is defined as

$${\mu }_t=C_k\rho k_{sgs}^{1/2}\mathsf{\Delta }f$$

(4)

The subgrid viscous stress tensor meets the equation as following:

$${\tau }_{ij}-\frac{2}{3}\rho k_{sgs}{\delta }_{ij}=-2C_k\rho k_{sgs}^{1/2}\mathsf{\Delta }f{\overline{S}}_{ij}$$

(5)

In this equation,

$$k_{sgs}=\frac{1}{2}\left(\overline{u_K^2}-{\overline{u}}_k^2\right)$$

(6)

$$\mathsf{\Delta }f\equiv V^{1/3}$$

(7)

$V$ is the volume of the grid. Moreover, subgrid kinetic energy can be obtained by the following equation:

$$\rho \left(\frac{\partial {\overline{k}}_{\mathrm{sgs}}}{\partial t}+\frac{\partial {\overline{u}}_j\overline{k_{\mathrm{sgs}}}}{\partial x_j}\right)=-{\tau }_{ij}\frac{\partial {\overline{u}}_i}{\partial x_j}-C_{\epsilon }\rho \frac{k_{sgs}^{3/2}}{\mathsf{\Delta }f}+\frac{\partial }{\partial x_j}\left(\frac{{\mu }_t}{{\sigma }_k}\frac{\partial k_{\mathrm{sgs}}}{\partial x_j}\right)$$

(8)

It should be noted that $C_k$, $C_{\epsilon }$ and ${\sigma }_k$ are constants in above formula. In present work, the implicit method based on density is applied to numerical simulation. The filter equation is discretized by finite volume method and a second-order upwind scheme is utilized to calculate the convection terms. Besides, a second-order central difference scheme is applied to calculating the viscous term, and a second-order implicit scheme is employed to discretize the time. Such numerical schemes [17] are widely used in previous studies on mixing layers.

2.2. Computational Model and Boundary Conditions

As a preliminary work concerning the proposed strategy and relevant enhanced-mixing mechanisms, the two-dimensional compressible N–S equations are solved. The two-dimensional simulation strategy is acceptable in the present study since the computational model is based on the inner flow of a scramjet engine, and a great deal of research with two-dimensional models has already been performed to reveal the relevant flow mechanisms in inner engine flows [31,32]. Moreover, the compressibility effects on mean velocity field and turbulent statistics are researched, employing the two-dimensional method since it is the first step towards understanding the realistic three-dimensional flow, which is much more complex.

The simulation in this paper is based on the model of Yee [33], where the length of the computational domain is 200 mm and the width is 40 mm, corresponding to the height of the inlet on both the upper and the lower sides is 20 mm. The low-speed inlet on the lower side is defined as inlet 1 and the high-speed inlet on the upper side is defined as inlet 2, as shown Figure 1. To enhance the mixing effect of the incoming flow and to simulate the random disturbance of the actual combustor, white noise disturbance is added to both the upstream and the downstream incoming flows. Unlike the previous approach of generating a shock wave through a wedge shock generator, a wall jet (16 < x < 17) at the upper boundary is used in the present work to generate a controllable shock wave with a jet nozzle, and the throat width is 1 mm.

The incoming flow parameters for this calculation are listed in

Table 1. Inflow parameters.

	Ma	V/(m·s−1)	P/kPa	T/K	Mc
inlet1	1.7	658.4	101	346.6	0.6
inlet2	5.6	986.3	101	76.7

.

The convective Mach number is defined as

$$M_c=\frac{u_1-u_2}{a_1+a_2}$$

(9)

and the subscripts 1 and 2 represent the lower and the upper layers of the flow field, respectively. The local sound velocity can be obtained from the relation:

$$a=\sqrt{\gamma RT}$$

(10)

As for the selection of the boundary condition, the pressure inlet condition is set in the inlet, outlet, and jet nozzle; the slip wall is applied to upper and lower wall. The jet nozzle is set just for the choked state, which is marked by the superscript ‘*’, and the nozzle flow will reach the maximum at this state. The choked flow calculation formula is as follows.

$$\frac{p^{\ast }}{p_0{}^{\ast }}=0.528$$

(11)

$p^{\ast }$ and $p_0{}^{\ast }$ are the static and the total pressure of the choked flow, respectively. The total pressure of the jet nozzle is set to fluctuate in the form of a sine curve around 0.02 times the total pressure of the upper incoming flow.

$$p_0^{\ast }=p_{02} \left[A\sin \left(2\pi ft\right)+\zeta \right]$$

(12)

$p_{02}$ is the total pressure of the upper incoming flow. Here, we set $A$ = 0.01 and $\zeta$ = 0.02, and $f$ takes different values in different cases to control the frequency of the jet. This paper sets a case of a mixing layer and two cases of jet shock interference, where $f$ is taken as 2$f_0$ and 0.2$f_0$, respectively, and $f_0$ is the characteristic frequency of an unstable mode of the free mixing layer, the value of which is $3\times {10}^4$ Hz. The frequencies of different jet shock waves are shown in

Table 2. Simulation conditions for three cases.

	Frequency of Jet	State of Jet Shock
Case 1	None	shock free
Case 2	2f0	high-frequency shock
Case 3	0.2f0	low-frequency shock

.

Because we focus on the interaction of shock wave and mixing layer, the grid is only densified in the the mixing area, and the grid of another area is evenly distributed, which can be seen in Figure 2. Furthermore, the boundary condition does not influence the effect of the shock wave on the mixing layer, so the upper and lower boundaries of the computational domain are set as a slip wall. The advantage of this slip condition over a non-slip condition is that the slip condition demands fewer grids at the boundary, and the shock wave is not smeared in the lower grid resolution. This method is used by many scholars in similar research [34].

In the present work, three sets of grids are used for numerical calculations to verify grid independence, and the number of grids is 1000 × 300, 2000 × 300 and 2500 × 300, respectively. The grids are uniformly distributed in the streamwise x-direction and the transverse grids are refined in the core area. Figure 3 compares the distribution of the mean Mach number along the stream direction at the centerline of the flow field under different grid conditions. It can be found that the value of the time-averaged Mach number along the streamwise direction at the centerline of the flow field under the grid conditions of 2000 × 300 and 2500 × 300 is approximately the same, while it is significantly underestimated under the poor grid condition (1000 × 300) between 0.05 m and 0.07 m in the streamwise x-direction. In order to economize on computational resources, a grid condition of 2000 × 300 is employed for calculation in this paper.

2.3. Code Validations

In order to verify the correctness of the numerical method used in this calculation, the velocity profiles at different streamwise locations are compared between the calculations in this paper and the Geobel–Dutton [35] experimental results, with all conditions consistent. The schematic diagram of the experimental device of a is shown in Figure 4.

The incoming flow conditions of the Geobel–Dutton experimental are listed in

Table 3. Incoming condition of Geobel–Dutton experiment.

Stream	Ma	V/(m·s−1)	T/K	P/kPa	Mc
1	1.9	700	334	49	0.45
2	1.3	400	215	49

.

Figure 5 gives the comparison between the calculated and experimental velocity profiles at x = 50 mm, x = 150 mm, and x = 200 mm, respectively: $\mathrm{H}$ is the distance between the upper and lower boundaries of the flow field, and ${\mathrm{y}}_0$ is the coordinate value of the centerline of the mixing layer; and ${\mathrm{U}}_2$ represents the time-averaged velocity value of the high-speed side, and ${\mathrm{U}}_1$ represents time-averaged velocity value of the low-speed side. Furthermore, we estimated the absolute error of the edge of the mixing layer with the most severe shear, which is 0.035, 0.061, and 0.018 at the streamwise direction of 50 mm, 150 mm, and 200 mm, respectively. The results demonstrate that the LES model used in our research can capture the velocity distribution characteristics of the flow field to a large extent, and the error is small, which proves that the algorithm of this paper is feasible.

Besides, we employed the method used by Smirnov et al. [36] to evaluate the error accumulation of numerical simulation, which is related to grid resolution and physical time. The results are shown in

Table 4. Error estimates.

Allowable Error (%)	Grid Resolution	Number of Time Steps	Accumulated Error	Allowable Number of Time Steps	Reliability
5	3000 × 200	72,000	9.97 × 10−6	1.81 × 1012	2.51 × 107

. As shown in the table, the higher the reliability, the smaller the error. We can conclude that the grid distribution and the numerical procedures can provide high reliability for the research.

Furthermore, we adopt the method used by Celic et al. [37] to estimate the discretization error in CFD calculation and the detailed description of this method can be found in Ref. [37].

Table 5. Sample calculation of discretization error.

N1, N2, N3	${\mathit{r}}_{21}$	${\mathit{r}}_{32}$	${\mathsf{\varphi }}_1$	${\mathsf{\varphi }}_2$	${\mathsf{\varphi }}_3$	p	${\mathit{e}}_{\mathit{a}}^{32}$	$\mathit{G}\mathit{C}{\mathit{I}}_{\mathit{f}\mathit{i}\mathit{n}\mathit{e}}^{32}$
750,000, 600,000, 300,000	1.118	1.414	984.422	984.3101	983.8613	1.09	0.046%	0.12%

illustrates the calculation procedure for selected grids: N1, N2, and N3 respectively represent the amount of the three sets of selected grids; $r_{21}$ is the ratio of the average length of the second sets of grids to the average length of the first set of grids, and $r_{32}$ is the ratio of the average length of the third sets of grids to the average length of the second set of grids; ${\varphi }_1$, ${\varphi }_2$, and ${\varphi }_3$ respectively represent the statistical velocity magnitude of the reference location, and the reference location is chosen at (100 mm, 20 mm) at the flow field of the three sets of selected grids; and P, $e_a^{21}$ and $GC{I}_{fine}^{21}$ are the key parameters for the error estimate, and the detailed calculation process can be found in Ref [37]. According to Table 1, the numerical uncertainty in the medium-grid solution for the velocity magnitude should be reported as 0.12%, which meets the demand of CFD.

2.4. DMD Method

The DMD method put forward by Rowley [38] can extract valuable information from the flow field snapshots obtained from computation or experimentation. The snapshot sequence is represented by the matrix $V_0^N$.

$$V_0^N=\left\{v_0,v_1,v_2,v_3,\dots ,v_N\right\}$$

(13)

The snapshot sequence is separated by a sampling time interval ∆t, and $v_i$ indicates the i-th data, while subscript 0 and $N$ represent the first and last element in the sequence. If the sampling interval is small enough, the data sequence can be regarded as linear dynamical systems.

$$v_{i+1}=A{v}_i$$

(14)

Here, A is a constant matrix. Moreover, the following equation can be deduced.

$$V_1^N=A{V}_0^{N-1}$$

(15)

Next, an eigenvalue decomposition of A is performed with $\mathsf{\Lambda }$ as its eigenvalue and $S$ as its eigenvector.

$$A=S^{-1}\mathsf{\Lambda }S$$

(16)

Since the dimension of the columns of matrix A are much larger than that of the rows, it is very complicated to directly perform singular value decomposition (SVD), so the similarity transformation of A is applied to obtain the similarity matrix of A, which simplifies the solution of eigenvalues and eigenvectors.

$$V_0^{N-1}=U\sum D^{\ast }$$

(17)

$$A=V_1^{N}D{\sum }^{-1}{U}^{\ast }$$

(18)

$$\tilde{A}=U^{\ast }AU=U^{\ast }{V}_1^{N}D{\sum }^{-1}$$

(19)

The eigenvalues of A are the same as those of $\tilde{A}$, and the eigenvectors of A have a transformation relationship with those of $\tilde{A}$. Therefore, the eigenvalues and eigenvectors of A can be obtained by performing SVD on $\tilde{A}$. After the above solution process, $S^{-1}$ and $\Lambda$ are expressed as follows.

$$S^{-1}=\left[y_1,y_2,\dots ,y_n\right]$$

(20)

$$\mathsf{\Lambda }=\left[\begin{array}{ccc}\begin{array}{cc}\begin{array}{c}{\lambda }_1\\ 0\end{array}& \begin{array}{c}0\\ {\lambda }_2\end{array}\end{array}& \begin{array}{c}\cdots \\ \ddots \end{array}& \begin{array}{c}0\\ ⋮\end{array}\\ \begin{array}{cc}⋮& \ddots \end{array}& \ddots & 0\\ \begin{array}{cc}0& 0\end{array}& \dots & {\lambda }_r\end{array}\right]$$

(21)

Matrix $W_i$ is introduced here.

$$W_i=S{v}_i$$

(22)

The following expression can be derived from the above formula.

$$v_i=S^{-1}{W}_i=S^{-1}{\wedge }^i{W}_0$$

(23)

$$W_0={(S^{-1})}^{-1}{v}_0$$

(24)

Thus, the flow field at any time can be obtained by:

$$v_i={\displaystyle \sum }_{k=1}^r{y}_k{\lambda }_k^i{W}_k^0$$

(25)

Here, $W_k^0$ is the k-th element of vector W₀. For an approximately linear system, there is a certain error in using Equation (24) to calculate $W_0$, so the least squares method is used to calculate $W_0$.

$$\underset{F}{minimize }J\left(W_0\right):=\Vert V_1^N-U\tilde{A}\Sigma D^{\ast }\Vert {}_F^2$$

(26)

Next, the data sequence can be described by the following expression.

$$\left[v_0,v_1,v_2,\dots ,v_n\right]=\left[y_1,y_2,\dots ,y_r\right]\left[\begin{array}{ccc}{W}_1^0& & \\ & W_2^0& \\ & \ddots & \\ & & W_r^0\end{array}\right]\left[\begin{array}{ccc}\begin{array}{cc}1& {\lambda }_1\end{array}& \cdots & {\lambda }_1^n\\ \begin{array}{cc}\begin{array}{c}1\\ ⋮\end{array}& \begin{array}{c}{\lambda }_2\\ \cdots \end{array}\end{array}& \ddots & ⋮\\ \begin{array}{cc}1& {\lambda }_r\end{array}& \cdots & {\lambda }_r^n\end{array}\right]$$

(27)

In this way, the DMD method completes the decomposition of the flow field. The left side of the equation is a data sequence of the original flow field, the first term on the right side of the equation is a spatial distribution of physical quantity, the second term is the amplitude, and the third term is the evolution of the time coefficient in each mode, the frequency $f$, and the growth rate $\xi$, which can be obtained from the three terms.

$$f_k=\frac{1}{\mathsf{\Delta }t}ln\left|{\lambda }_k\right|$$

(28)

$${\xi }_k=\frac{1}{\mathsf{\Delta }t}arg\left({\lambda }_k\right)$$

(29)

3. Results

3.1. Flow Structure Visualization

The numerical schlieren technique proposed by Hadjadj [39] is employed here to obtain the instantaneous flow structures; it displays the flow field image based on density gradient. The strategy can be defined as

$$S\left(x,y\right)=\beta exp\left(-\frac{\kappa \left|\nabla \rho \right|}{{\left|\nabla \rho \right|}_{max}}\right)$$

(30)

$$\left|\nabla \rho \right|=\sqrt{{(\frac{\partial \rho }{\partial x})}^2+{(\frac{\partial \rho }{\partial y})}^2}$$

(31)

Here, the parameter $\beta$ determines the color corresponding to the zero gradient, the parameter $\kappa$ governs the amplification degree of the small gradients, and the flow field is usually better displayed when the $\beta$ and $\kappa$ are taken as 0.8 and 15, respectively.

Figure 6 shows the numerical schlieren contour for present three cases. Due to the Kelvin–Helmholtz(K–H) instability generated by different upper- and lower-layer speeds, the typical coherent structure called K–H vortices in the mixing layer rolling up can be seen in all three cases. The K–H vortices gradually grow, and the phenomenon of pairing and merging between vortices appears. Browand [40] confirmed that the dynamic of pairing and merging between adjacent vortex structures is an important mechanism of free mixing layer growth. Apart from that, the shocklets in three cases are clearly visible. The shock waves generate by the local flow being impeded by the large-scale structure in the turbulent mixing layer when the local Mach number relative to the large-scale structure is greater than 1.0 [41]. Generally, shocklets appear in the high Mach number region of the flow field, corresponding to the upper flow area in our cases. It can be found that shocklets become more obvious due to the compression effect induced by jet shock in Case 2 and in Case 3.

In Case 1, the mixing layer develops around the central region of the flow field. However, in Case 2 and Case 3, due to the introduction of jet shock, the mixing layer is deflected downward after the incident shock, then deflected upward after the reflection shock. Meanwhile, the jet shock is first refracted through the mixing layer, then reflected upward on the lower wall, and then again refracted on the mixing layer. Refraction shock interacts with shocklets and finally becomes segmented state. Relative to Case 1, a new phenomenon of vortices rolling up closer to the inlet can be seen in Case 2 and Case 3, indicating that there are subsonic regions at the shock/mixing interaction area. This leads disturbance to transmit backward and can promote the vortices to roll up. Compared Case 2 with Case 3, it can be found that incident shock will be segmented in Case 2, while relatively straight in Case 3. Fast jet frequency in Case 2 makes the shock wave propagation speed lower than the jet disturbance speed, eventually it causes the ‘disconnection’ phenomenon. Overall, the vortex structures belonging to the interaction region possess characteristics of coexistence of both large- and small- scale vortices, and this dynamic is of great importance for promote mixing. Thus, detailed analysis concerning the structure behaviors is essential.

To better analyze the interaction process between jet shock and mixing layer in case 3, the detailed visualization around the interaction position is presented in Figure 7. One important phenomenon is that compared to the interaction process in case 2, in case 3 the interaction position is not fixed. In Figure 7a the interaction position located at x = 83 mm approximately, while in Figure 7b at about x = 62 mm. Indeed, by analyzing our calculated results it can be found that these two positions correspond to the farthest and nearest interaction locations in the flow field. Through reviewing previous literatures focusing on shock induced by wedge shock generator and its interaction process with mixing layer [23], it can be concluded that the interaction position was fixed and can be called the ‘point action mode’ (PAM). However, in present work new mode never reported before is found here and we defined it as ‘region action mode’ (RAM). This newly found mode plays important roles in the mixing enhancement process. Since as mentioned above, at different time the interaction point moves between the Farthest interaction point (FIP) and the nearest interaction point (NIP).

Numerous researches concerning mixing augmentation affected by wedge inducing oblique shock concluded that at the interaction point, the vortices distort and breakdown into small-scale structures [42]. This indicates that in present work, different interaction points induced by WJISW in a low frequency condition can result in the unique dynamical behaviors of vortex structures located between NIP and FIP. These dynamical behaviors including breakdown and deformation of vortex structure. The breakdown of vortex structure is related to the turbulent fluctuation of the flow field. Next, we clarify the relationship between the unique dynamic behavior of the vortex structure in RAM and the fluctuation characteristics of the flow field by quantitatively analyzing the distribution of the turbulent intensity of the flow field.

3.2. Turbulence Intensity Analysis

The turbulence intensity reflects the momentum exchange characteristics on both sides of the interface caused by the motion of fluid particle. Larger turbulence intensity means more significant three-dimensional characteristics of the flow, and large-scale structures are easily broken into small-scale structures [43]. In order to further explore the relationship between the unique dynamic behavior of the vortex and the turbulent fluctuation characteristics of the flow field under RAM, the distribution of maximum turbulence intensity at intersecting surface versus the streamwise x-direction and its integration over the entire streamwise x-direction are compared and analyzed for the three cases, as shown in Figure 8 and Figure 9, respectively.

The turbulence intensity characterized by four indexes, including the dimensionless streamwise Reynolds normal stress ($RS11$), transverse Reynolds normal stress ($RS22$), Reynolds shear stress ($RS12$), and turbulent kinetic energy ($TKE$) are calculated as follows,

$$RS11=\frac{\overline{\rho }\widetilde{u\prime u\prime }}{{\rho }_1\mathsf{\Delta }{u}^2}$$

(32)

$$RS22=\frac{\overline{\rho }\widetilde{v\prime v\prime }}{{\rho }_1\mathsf{\Delta }{u}^2}$$

(33)

$$RS12=\frac{\overline{\rho }\widetilde{u\prime v\prime }}{{\rho }_1\mathsf{\Delta }{u}^2}$$

(34)

$$TKE=\frac{1}{2}\left(\frac{\widetilde{u\prime u\prime }+\widetilde{v\prime v\prime }}{\mathsf{\Delta }{u}^2}\right)$$

(35)

Here, $\overline{\rho }$ is the average density, ${\rho }_1$ is the density of lower flow, $u\prime$ is streamwise pulsation velocity, $v\prime$ is transverse the pulsation velocity and $\mathsf{\Delta }u$ is the velocity difference between upper and lower flow. It can be seen that all the turbulence intensity characteristic values of the benchmark case have ‘one peak’. The reason of ‘one peak’ is the roll-up of vortices in the flow field, and the strong shearing effect leads to large pulsations in the mixing layer, resulting in a peak in the turbulence intensity term at this location. Whereas, in both Case 2 and Case 3 two higher peaks are detected along the streamwise x direction due to the effect of incident and reflected shock waves. Previous researches demonstrated that this phenomenon of ‘multiple peaks’ [44] can significantly promote the mixing process, so shock wave is beneficial to mixing augmentation. By comparing Case 2 and Case 3, it can be found that the growth rate of turbulence intensity is almost the same in both cases. However, under the influence of the RAM in case 3, the area where the shock wave impacts on the mixing layer is larger, which makes the turbulence intensity peak higher and the range of large values of turbulence intensity wider. This distribution property indicates that the RAM has a considerable effect on the large-scale enhancement of turbulence intensity in the mixing layer. Comparing the distribution of different turbulence intensity characteristics in the streamwise x-direction, it can be seen that the turbulence intensity of $RS22$ reaches the peak earlier than that of $RS11$ in Case 3, making the enhancement of the transverse Reynolds stress peak in this condition is not as significant as that of the streamwise Reynolds stress. This is because the shock oscillation is in streamwise direction in the RAM. Thus, the effect on the transverse Reynolds stress is not significant. The decline after $RS12$ reaches its peak is gentler than that of $RS22$. This is because the shear Reynolds stress is a combined measure of turbulent pulsations in both the streamwise and transverse directions, and the superposition of the streamwise delayed peak and the transverse trough makes the tangential Reynolds stress decrease more slowly. The development of turbulent kinetic energy along the streamwise direction is similar to the streamwise Reynolds stress, because the turbulent kinetic energy is the superposition of the streamwise and transverse Reynolds stress. Furthermore, the shock wave oscillated in streamwise x-direction causes the streamwise Reynolds stress to dominate the flow fluctuation.

In summary, the turbulence intensity in the area of shock interaction is significantly increased under the RAM, and that in other regions of the flow field can also be effectively increased. Combined with the numerical schlieren observed in Figure 6, it can be seen that the vortex structure is more fragmented in the region of large values of turbulent intensity in the flow field, which proves that the large-scale growth of turbulent intensity in the flow field under the RAM is the principal reason for the fragmentation of vortex structure in Case 3.

Integrating the turbulence intensity maximum values along the streamwise direction can reflect the turbulence characteristics of the whole flow field more clearly [44]. Figure 9 shows the comparison of the integrated values of $RS11$, $RS22$, $RS12$, and $TKE$ along the streamwise x-direction for the three cases. It can be found that the value in Case 2 is slightly higher than that in Case 1, while the value in Case 3 is much higher than that in Case 2. In the comparison of the three cases in $RS22$, the difference between Case 3 and Case 2 is not much higher than the difference between Case 2 and Case 1. This is also consistent with the conclusion obtained from the above analysis of the distribution of the turbulence intensity maximum value along the streamwise x-direction. Therefore, considering the effect of enhancing the whole flow-field turbulence intensity, WJISW in a low frequency condition is much superior to WJISW in high frequency condition.

3.3. Mixing Layer Growth Process

To quantitatively evaluate the effect of WJISW on the mixing characteristics of the mixing layer, a velocity thickness index ${\delta }_b$ is adopted, and it can be defined as the transverse distance between two positions with normalized mean velocity of 0.1 and 0.9,

$${\delta }_b=y_{\overline{u}=0.9}-y_{\overline{u}=0.1}$$

(36)

where the normalized time-averaged velocity is defined as

$$\overline{u}=\frac{U\left(y\right)-U_1}{U_2-U_1}$$

(37)

The closer the normalized time-averaged velocity is to 1, the flow velocity is closer to the high-speed side and this position is more influenced by the high-speed side. The closer the normalized velocity is to 0, the closer the flow velocity is to the low-speed side, and this position is more influenced by the low-speed side. The distribution of velocity thickness versus streamwise x-direction for the three cases are depicted in Figure 10. The velocity thickness of the free mixing layer shows a linear growth characteristic, which is consistent with the results of previous work [42]. The velocity thickness in Case 2 decreases after two shock wave interaction points. Due to the compressibility of the shock wave, the shock wave causes a density increase, making the velocity thickness over a large area of the flow field in Case 2 smaller than that in Case 1. Nevertheless, the velocity thickness growth rate increases downstream of the shock interaction point, and it finally catches up with the velocity thickness of Case 1 at x = 200 mm. For Case 3, due to the unique advantage of the RAM, velocity thickness decreases after a long shock interaction zone. Particularly, velocity thickness increases steeply in the shock interaction zone owing to the dramatic increase in turbulence intensity at the shock interaction position, resulting in a sharp oscillation of the flow field. The increase in velocity thickness of the mixing layer in Case 2 at the shock interaction position is hardly noticed because the shock wave only acts on a single point in the mixing layer, which has little effect on the turbulence intensity of the flow field. Downstream of the shock interaction position, Case 3 has the same growth rate of velocity thickness as Case 2, and the final velocity thickness at x = 200 mm in Case 3 has increased by 62% compared to Case 1, which greatly improves the mixing level.

Because the interaction between shock wave and the mixing layer will cause the total pressure loss of the flow field, and the total pressure loss coefficient is an important index to evaluate the actual operation of the engine, we perform a detailed analysis of the total pressure loss coefficient, which is defined as,

$$\begin{array}{c}\zeta =1-\frac{\int {\mathrm{P}}_{\mathrm{x}}\mathsf{\rho }\mathrm{udA}}{\int {\mathrm{P}}_0\mathsf{\rho }\mathrm{udA}}\end{array}$$

(38)

where ${\mathrm{P}}_{\mathrm{x}}$ is the time-averaged total pressure of a certain streamwise position and ${\mathrm{P}}_0$ is the time-averaged total pressure at the inlet position; $\mathsf{\rho }$ and $\mathrm{u}$ represent the time-averaged density and velocity, respectively. The total pressure loss distribution along the streamwise x direction for three cases is shown in Figure 11. It can be seen from the figure that when jet shock interacts with the mixing layer, although the total pressure loss coefficient increases to a large extent at the jet ejection position, the growth rate of the total pressure coefficient obviously slows down at x = 75 mm. Eventually, at x = 200 mm, the total pressure loss of Case 3 is larger than Case 2, and the total pressure loss of Case 2 is larger than Case 1. Notably, due to the introduction of shock waves, many techniques [45,46] for mixing enhancement cause a large loss in total pressure. In comparison, the flow loss of our scheme is acceptable. Besides, the total pressure loss of Case 3 is only slightly larger than that of Case 2, whereas Case 3 can achieve a higher mixing enhancement than Case 2. Therefore, in terms of comprehensive measurement of total pressure loss and mixing enhancement capacity, WJISW in a low frequency condition is a better choice.

3.4. DMD Analysis

As a newly proposed mixing technique, it can be found that the RAM generated by WJISW at a low frequency provides a significant degree of mixing enhancement, and the reveal of its mechanism contributes to the design optimization of the mixing augment device. In this paper, our cases are nonlinear dynamic systems, and the original flow field includes a large amount of information. Thus, it is not easy to grasp the focus of the data. It is worth noting that as an emerging reduced-order method applied in the flow field, DMD can realize the decomposition of the original data to obtain a specific flow field structure with a dominant frequency and growth rate [38]. This method can provide a convenient means for flow field dimensionality reduction and information acquisition, and it has achieved some application in the field of compressible flow [47,48,49,50], especially in revealing the typical structural characteristics of the flow [51]. In this paper, 1200 snapshots of the flow field with a time interval of 1.5μs are extracted from the transient flow field of three cases close to 8 circulation cycles for DMD analysis. The penalty terms Υ of the three cases in our paper are chosen as 1 × 10⁴, 1 × 10⁵, and 1 × 10⁴, respectively. After optimizing the amplitude of the modes by the DMD method, the number of non-zero modes of the three cases reduced from the original 1199 to 385, 423, and 401, respectively, which reflects the advantage of the DMD method in dimensionality reduction.

To confirm the reliability of the DMD method, the original flow field and the reconstructed flow field employing DMD method are compared and analyzed, as shown in Figure 12. The flow fields are characterized utilizing the dimensionless vorticity ${\omega }_{nor}$ and it can be defined as

$${\omega }_{nor}=\left|\frac{{\omega }_z{\delta }_{\left(0\right)}}{u_2-u_1}\right|$$

(39)

$${\omega }_z=\frac{\partial u}{\partial y}-\frac{\partial v}{\partial x}$$

(40)

It can be seen that the DMD reconstructed flow field has basically restored the vortex behavior in the flow field, and the position and size also match with the original flow field. The roll-up, pairing, and merging of vortex structures can be clearly seen from the DMD reconstructed flow field, indicating that the vortex structure dynamics, which are important for the development of the mixing layer, can be captured by the DMD method.

Furthermore, a quantitative comparison is performed by extracting the dimensionless vorticity on the centerline of the reconstructed flow field and the original flow field, as shown in Figure 13. From the figure, it can be more intuitively seen that the overall trend of the reconstructed values remains consistent with the original values, which also coincides with the results analyzed above. It can also be seen from the figure that the dimensionless vorticity on the centerline of the flow field of Case 1 fluctuates more, while Case 2 and Case 3 fluctuate less after x = 70 mm, which is due to the fact that the mixing layer impacted by WJISW will be deflected by the shock wave interaction, so the region of the most intense turbulent pulsation is not on the centerline of the flow field.

Figure 14 shows the frequencies and energy proportion distribution of DMD modes in three cases. Except for the first mode, all the modes are conjugate modes. The conjugate mode refers to modes with opposite frequencies, but with the same growth rate and modal energy proportion. Thus, here only frequencies greater than 0 are shown. The black points represent the values with eigenvalues less than or equal to 1 or a modal growth rate less than or equal to 0 (the two concepts can be found equivalent by simple logarithmic calculation), which correspond to stable or periodic modes, while the red points represent the values with eigenvalues greater than 1 or a modal growth rate greater than 0, which corresponds to unstable modes. Moreover, blue triangles correspond to the first eight modes with a higher energy proportion, which play dominant roles in the flow field and will be mainly analyzed below. It can be found that the frequency distribution of the first eight modes in Case 1 is concentrated around 100 kHz, while the frequency distribution of the first eight modes of Case 2 and Case 3 is more scattered, ranging from 50 kHz to 240 kHz and 6 kHz to 180 kHz, respectively. Meanwhile, it can also be seen that the frequencies of the first few modes in Case 3 are much smaller than those in Case 1 and Case 2.

Figure 15, Figure 16 and Figure 17 show the first eight orders of spatial mode contours in three cases, respectively. The 1st modes show the time-averaged solution of dimensionless vorticity [52], while the rest of the modes reflect the oscillation characteristics generated by the flow field pulsation. The pulsation of the flow field in Case 1 mainly comes from the mixing of the fluid above and below the mixing layer, that is to say, the roll-up of K–H vortices. The upper and lower boundaries of the topology extracted from each mode can be regarded as the upper and lower boundaries of the mixing layer, and its width increases with the distance in streamwise x-direction, indicating that the thickness of the mixing layer grows in the streamwise x-direction. The analysis of the 2nd to 8th order modes shows that the topological structure, with alternating positive and negative values of dimensionless vorticity, appears only in the core region of the flow field. This structure can be considered as the spatial form of the wave, and the distance between the centers of positive and negative values can be regarded as half a wavelength [53]. The values in other regions of the flow field are zero, indicating that pulsations occur mainly in the mixing layer. Since the original flow field can be seen as a superposition of the DMD modes obtained from its decomposition, the vortices’ behavior can be seen in the superposition image. However, in the decomposed modes of different frequencies, only the wave structure can be seen in the flow field image. The above findings together converge to a conclusion that the vortex in the flow field is essentially a superposition of waves of different frequencies.

It can be seen from Figure 16 that the first few order modes in Case 2 cover a larger wavelength range than Case 1, and the 2nd, 3rd, 4th, and 5th order spatial modes of the flow field have significant shock effects, while the 6th, 7th, and 8th order spatial modes do not. This indicates that the 2nd, 3rd, 4th, and 5th order spatial modes of the flow field are shock-induced modes, and the energy of these modes accounts for a higher proportion of the total flow field energy compared to the later modes that are not shock-induced, i.e., K–H instability-induced modes. Therefore, the jet shock-induced modes play a dominant role in the flow field.

From Figure 17, we can see that the topology of the first few orders modes in Case 3 are very huge, especially the wavelength of the 2nd order mode has exceeded half of the flow field area. The 2nd, 3rd, 4th, and 6th order spatial modes can clearly see the effect of the shock wave, indicating that they are modes induced by jet shock wave, and the topological stratification can be seen more clearly than the other modes, induced only by K–H instability. Moreover, it can be clearly seen that the wavelengths of the modes induced by the shock are longer than those of the first two cases, but the wavelengths of these modes are varying in length. The modes with long wavelengths are favorable to induce large vortex structures, while the modes with short wavelengths are favorable to induce small vortex structures, resulting in the coexistence of large and small scales of vortex structures in the flow field. Entrainment and engulfment of large-scale structures allow surrounding fluids to enter the mixing region, and the nibbling effect of the small-scale structure facilitates the expansion of the surface area of scalar contact in the mixing layer. These dynamic behaviors accelerate the exchange of mass, momentum, and energy between the upper and lower layers of the flow [54], which is the reason why Case 3 has the best degree of mixing enhancement in the three cases. Comparing these three cases, Case 1 has a single vortex size, and Case 2 can induce different scales of vortex structures, but the largest vortex size is still small, and the large-scale entrainment effect is insufficient. Whereas, the interaction of large- and small-scale vortices in Case 3 leads to the fragmentation of the vortex structure in the flow field, which is consistent with the structure of the flow field in the numerical schlieren of Case 3.

Figure 18 compares the evolution of the 2nd–8th order mode amplitudes with time in three cases. In the figure, the horizontal ordinate is the dimensionless time $t_{nor}$.

$$t_{nor}=\frac{t}{0.0000015}$$

(41)

From Figure 18a, it can be seen that there is not much difference between the modal amplitude and frequency of the first few order modes of Case 1 and the 2nd mode, which has a lower initial amplitude and is ranked ahead of the 3rd mode. This is because the DMD mode ranking criterion in this paper not only considers the initial amplitude, but also focuses on its evolution in the entire time history of available snapshots. The 3rd mode shows a weak decay characteristic with a negative growth rate, and the amplitude decreases slowly with time, while the 2nd mode is an unstable mode with a positive growth rate, and the amplitude increases gently with time. Therefore, the overall contribution of the 3rd mode to the flow field is greater than that of the 2nd mode, making the 2nd mode ranked before the 3rd mode, which also reflects the advantage of the modal ranking criterion in this paper. A large difference in the modal frequencies of the first few order modes of the flow field in Case 2 is presented in Figure 18b. Comparing with Figure 16, we can find that the mode with a small frequency corresponds to a long spatial mode wavelength and vice versa. The modes with the smallest frequencies and the longest spatial mode wavelengths are the 2nd mode and the 7th mode, respectively, and these structures are conducive to the induction of larger vortices in the flow field, which are found to be the modes induced by the jet shock and K–H instability, respectively, corresponding to Figure 16. However, the amplitude of the 2nd mode is more than twice that of the 7th mode, which also indicates that the modes induced by the jet shock are more dominant compared to those induced by K–H instability. It can be seen from Figure 18c that the frequencies of the first few order modes in the flow field in Case 3 are very different, corresponding to the larger wavelength range covered by the spatial modes in Figure 17. In addition, the 2nd mode amplitude is nearly twice as large as the remaining other modes, indicating that it accounts for a higher proportion of energy in the flow field.

The energy proportion, frequency, and growth rate of the first eight modes for three cases are analyzed to better research the vortex dynamics and mixing process.

Table 6. Modes’ characters of the first eight modes in Case 1.

Order	Proportion	Frequency/Hz	Growth Rate
1	0.110951	5 × 10−12	6
2	0.005012	113,335	255
3	0.004626	106,678	−1074
4	0.004621	90,923	575
5	0.004261	88,690	393
6	0.004004	97,742	−428
7	0.003983	80,179	−265
8	0.003970	106,430	89

,

Table 7. Modes’ characters of the first eight modes in Case 2.

Order	Proportion	Frequency/Hz	Growth Rate	Frequency/fj
1	0.10638	−8 × 10−11	2.65	0
2	0.01716	60,000	−40	1
3	0.00802	119,998	−50	2
4	0.00495	179,999	79	3
5	0.00421	239,974	92	4
6	0.00385	97,445	−25	1.62
7	0.00375	52,148	−631	0.87
8	0.00359	118,078	149	1.96

and

Table 8. Modes’ characters of the first eight modes in Case 3.

Order	Proportion	Frequency/Hz	Growth Rate	Frequency/fj
1	0.05807	−2 × 10−10	−7.69	0
2	0.01402	6018	−94	1
3	0.00703	12,028	10	2
4	0.00413	18,055	67	3
5	0.00243	73,954	−239	12.32
6	0.00239	24,025	−300	4
7	0.00234	67,617	−511	11.26
8	0.00229	186,742	−225	32.12

show the first eight orders of modal characteristics for the three cases, respectively. The frequencies of the first few modes of the flow field are high, except the 1st mode in Table 6. In addition, there is little difference in the 2nd- to 8th-order mode energy, suggesting that the energy in the flow field is dispersed, and no dominant modes control the development of the mixing layer.

From Table 7, we can see that the energy proportion of the first-order mode in Case 2 is 0.10638, which is lower than that in Case 1, which indicates that the perturbation energy proportion of Case 2 in the flow field is higher. Furthermore, the energy of the 2nd order and the 3rd order mode is much larger than the later modes, which indicates that these two modes of Case 2 are dominant in the flow field. Moreover, it can also be found that the energy of the 2nd, 3rd, 4th, and 5th order modes of the flow field are one, two, three, and four times of the wall jet frequency, respectively, which further confirms that the 2nd to 5th order modes are multiple frequency modes induced by the jet shock wave as mentioned above. It is noteworthy that the growth rate of the 2nd to 3rd order modes, which account for a relatively large amount of energy in the flow field, is negative and shows a weak decay, suggesting that the shock-induced modes will eventually decay completely if the flow field has a long streamwise length. In addition, although the 7th order mode is also induced by K–H instability, it has a smaller frequency and corresponds to a longer wavelength than the first few order modes in Case 1, perhaps due to the effect of the shock wave on the K–H instability of the flow field.

As can be seen from Table 8, the energy proportion of first-order mode in Case 3 is 0.05807, which is lower than that of both Case 1 and Case 2, indicating that the perturbations energy proportion of Case 3 in the flow field is higher. Similarly, the mode with multiple frequency of the jet ranks high in the energy proportion of each mode, which is the same as the conclusion reached in Case 2: the mode induced by the shock wave dominates the flow field. Moreover, the low frequency jet can excite high multiple frequency modes, suggesting that high frequency flow control can be achieved through low-frequency excitation. The decay rate of the 6th mode, which is a high-frequency mode with a mode frequency four times the excitation frequency, is found to be fast. As the decay rate of the higher multiple frequency mode is larger, the distance acting in the flow field is shorter, which provides a new idea to guide the design of the jet position. If the jet excitation is continuously introduced in front of the position where the higher multiple frequency mode decays more, the jet shock can interact on this position, which can make the low-frequency excitation mode and the high multiple frequency mode together dominate in the flow field, obviously enhancing the mixing degree.

4. Conclusions

In this paper, we research the mixing enhancement characteristics induced by the strategy of WJISW in a supersonic flow field through large eddy simulations. The instantaneous flow structures, turbulence pulsation characteristics, mixing layer growth process, and modal decomposition results are obtained and analyzed. In order to accentuate the mixing enhancement effect of this strategy, the free mixing layer under the same conditions is compared with it. The following conclusions can be drawn:

(1)

WJISW makes the mixing layer deflected. In particular, WJISW in a low frequency condition can result in the movement of the interaction position between NIP and FIP in the flow field. This new mode never reported before can be defined as ‘region action mode’ (RAM). In this mode, the vortex structures exhibit unique dynamic behavior, including breakdown and deformation. This is significantly different from the flow field of the mixing layer induced by a wedge shock generator;

(2)

The turbulence intensity of the mixing layer under WJISW appears multi-peak in the streamwise x-direction, which can significantly promote the mixing process. In addition, WJISW in a low frequency condition can make the turbulence intensity peaks higher in the streamwise x-direction and the range of large values wider, indicating that the ‘region action mode’ found in this paper has a considerable influence on the increase of turbulence intensity in the mixing layer, which significantly enhances the turbulent pulsation in the flow field;

(3)

The growth characteristics of the mixing layer under WJISW in a high frequency condition are similar to that of the mixing layer induced by wedge studied previously. The mixing layer thickness decreases after the shock interaction point, while the growth rate increases downstream, and eventually the mixing layer thickness will catch up with the free mixing layer. In contrast, under WJISW in a low frequency condition, the velocity thickness of the region where the shock wave acts with the mixing layer increases sharply, making the final velocity thickness increase by 62%, compared with the free mixing layer, and it greatly improves the mixing level;

(4)

The modal decomposition results show that the shock-induced modes dominate the flow field under WJISW. The coexistence of large- and small-scale structures in the mixing layer under WJISW in a low frequency condition is due to the dominance of the shock-induced modes. The entrainment and engulfment of large-scale structures promote surrounding fluids into the mixing region, and the nibbling effect of the small-scale structure facilitates the expansion of the surface area of scalar contact in the mixing layer, accelerating the exchange of mass, momentum, and energy between the upper and lower layers of the flow, both of which can effectively promote the mixing process;

(5)

The proposed strategy well shows the potential application in Scramjet and RBCC combustors for hypersonic propulsion systems. The newly observed ‘region interaction mode’ in the flow field influenced by wall-jet induced shock waves can obviously accelerate the mixing process. Future work concerning the optimization is essential to obtain a deeper understanding of the proposed strategy.