The bistability of curved compression ramp flows

The governing equations solved are dimensionless 3D Navier–Stokes equations for unsteady, compressible flow in curvilinear coordinates,

$$\frac{\partial\mathbf{Q}}{\partial t}+\frac{\partial\left(\mathbf{F}_{c}+\mathbf{F}_{v}\right)}{\partial\xi}+\frac{\partial\left(\mathbf{G}_{c}+\mathbf{G}_{v}\right)}{\partial\eta}+\frac{\partial\left(\mathbf{H}_{c}+\mathbf{H}_{v}\right)}{\partial\zeta}=0,$$

(1)

where $\mathbf{Q}$ is the conservative vector flux, $\mathbf{F}_{c}$, $\mathbf{G}_{c}$ and $\mathbf{H}_{c}$ are the invicid convection fluxes, $\mathbf{F}_{v}$, $\mathbf{G}_{v}$ and $\mathbf{H}_{v}$ are the viscous fluxes. Here, $\mathbf{Q}$, $\mathbf{F}_{c}$ and $\mathbf{F}_{v}$ are defined as

$$\mathbf{Q}=J\left(\begin{array}[]{c}\rho^{*}\\ \rho^{*}\mathbf{U}\\ \rho^{*}e^{*}\end{array}\right),\mathbf{F}_{c}=J\left(\begin{array}[]{c}\rho^{*}\mathbf{U}^{T}\mathbf{J}_{\xi}\\ \rho^{*}\mathbf{U}^{T}\mathbf{J}_{\xi}\mathbf{U}+p^{*}\mathbf{J}_{\xi}\\ \left(\rho^{*}e^{*}+p^{*}\right)\mathbf{U}^{T}\mathbf{J}_{\xi}\end{array}\right),\mathbf{F}_{v}=J\left(\begin{array}[]{c}0\\ \frac{1}{\operatorname{Re}_{\infty}}\boldsymbol{\tau}\mathbf{J}_{\xi}\\ \left(\frac{1}{\operatorname{Re}_{\infty}}\mathbf{U}^{T}\boldsymbol{\tau}-\frac{1}{(\gamma-1)M_{\infty}^{2}\operatorname{Pr}}\mathbf{q^{T}}\right)\mathbf{J}_{\xi}\end{array}\right),$$

(2)

where $\mathbf{J}_{\xi}=(\xi_{x},\xi_{y},\xi_{z})^{T}$ are the metric coefficients and $J$ is the determinant of the Jacobian matrix transforming the Cartesian coordinates $(x,y,z)$ into the computational coordinates $(\xi,\eta,\zeta)$. The velocity vector $\mathbf{U}$, the heat flux vector $\mathbf{q}$ and the stress tensor $\boldsymbol{\tau}$ are defined as

$$\mathbf{U}=\left(\begin{array}[]{c}u^{*}\\ v^{*}\\ w^{*}\end{array}\right),\quad\mathbf{q}=\left(\begin{array}[]{c}q_{\xi}^{*}\\ q_{\eta}^{*}\\ q_{\zeta}^{*}\end{array}\right),\quad\boldsymbol{\tau}=\left(\begin{array}[]{c}\tau_{\xi\xi}^{*},\tau_{\xi\eta}^{*},\tau_{\xi\zeta}^{*}\\ \tau_{\eta\xi}^{*},\tau_{\eta\eta}^{*},\tau_{\eta\zeta}^{*}\\ \tau_{\zeta\xi}^{*},\tau_{\zeta\eta}^{*},\tau_{\zeta\zeta}^{*}\end{array}\right).$$

(3)

The total energy $e^{*}$ is defined as

$$e^{*}=\frac{u^{*2}+v^{*2}+w^{*2}}{2}+\frac{1}{\gamma-1}\frac{p^{*}}{\rho^{*}}.$$

(4)

The other four flux terms have the similar forms: $\mathbf{G}_{c}$ and $\mathbf{H}_{c}$ are similar in form to $\mathbf{F}_{c}$; $\mathbf{G}_{v}$ and $\mathbf{H}_{v}$ are similar in form to $\mathbf{F}_{v}$. The viscosity $\mu$ is determined with the Sutherland’s law, and perfect gas equation, relating the pressure $p$ the density $\rho$ and the temperature $T$, is used to close the equations set. The Prandtl number $Pr=0.7$ and the specific heat ratio $\gamma=1.4$ are chosen in the simulations. The non-dimensional variable are normalized using the inflow parameters: $\rho^{*}=\rho/\rho_{\infty}$, $u^{*}=u/u_{\infty}$, $v^{*}=v/u_{\infty}$, $w^{*}=w/u_{\infty}$, $T^{*}=T/T_{\infty}$, $e^{*}=e/u_{\infty}^{2}$ and $p^{*}=\rho^{*}T^{*}/(\gamma Ma_{\infty}^{2})$. The reference length is chosen as $1$ mm.

In terms of the numerical methods, the time integration is performed by the third-order TVD-type Runge–Kutta method; the inviscid fluxes are discreted by the fifth-order WENO method (Jiang & Shu, 1996); the viscid fluxes are discreted with the sixth-order central difference scheme. The DNSs are conducted with the in-house code OPENCFD-SC, whose capability has been well examined (Hu et al., 2017; Xu et al., 2021), especially in studying compression ramp flows.

The choices of the compression ramp geometry (DCR) and the flow parameters are based on the recent shock-tunnel experiments (Roghelia et al., 2017a, b; Chuvakhov et al., 2017). As the DCR configuration shown in Fig. 2, the DCR configuration is mainly composed of two parts, the flat plate and the ramp, both of which have the same length $L=100$mm and width $W=30$mm. The ramp angle is $\alpha=15^{\circ}$. The flow parameters are shown in Tab. 1, including the Mach number $Ma_{\infty}=7.7$, the Reynolds number $Re_{\infty}=4.2\times 10^{6}$, the velocity $u_{\infty}=1726\text{m s}^{-1}$, the density $\rho_{\infty}=0.021\text{kg m}^{-3}$, the temperature $T_{\infty}=125\text{K}$ and the pressure $p_{\infty}=760\text{Pa}$. The isothermal wall condition is used with $T_{w}=293\text{K}$, since the run time of the shock tunnel is very short.

The computational domain, i.e., the mesh region, is shown in Fig. 3. The streamwise region is $-4\text{mm}\leq x\leq 196.6\text{mm}$ (the leading edge is at $x=0$ mm); the wall-normal height is $25$ mm; the spanwise width is $30$ mm. Two mesh resolutions, $1020\times 250\times 200$ (M1) and $1420\times 250\times 300$ (M2), are considered for time and grid convergence studies. Both M1 and M2 cluster the grids near the wall with the first wall-normal grid height being fixed to $\Delta y_{\text{w}}=0.008$ mm, yielding the non-dimensional height $\Delta y_{\text{w}}^{+}\approx 0.3$ mm at $x/L\approx 0.5$. Spanwise grid spaces of both M1 and M2 are uniform, as $\Delta z=0.15$ mm (M1) and $0.1$ mm (M2). The streamwise grid spaces are uniform as $\Delta x_{\xi}=0.2$ mm in M1. To investigate the pressure gradient effects, near the separation point $x_{S}/L\approx 0.55$ and the reattachment point $x_{R}/L\approx 1.28$ (also shown in Fig. 4) in M2, the $\Delta x_{\xi}$ spaces are severally clustered with $200$ points in the streamwise direction.

For the initial conditions, the initial flow fields for M1 and M2 are both extruded (or extended) spanwisely from the 2D simulation result (denoted as D${}^{\text{I}}_{0}$), in which field the separation point $x_{S}/L\approx 0.55$ and the reattachment point $x_{R}/L\approx 1.28$. For the boundary conditions, the free stream is set both at the numerical inlet boundary ($x=-4$mm) and the upper boundary; no-slip and isothermal ($T_{w}=293$K) conditions are set for $x\geq 0$mm on the wall; the extrapolation condition is set at the outflow boundary.

The 3D instantaneous flow fields of M1 and M2 are shown in Fig. 4. In the central $x-y$ planes, it clearly shows that the SW configurations include the leading edge SW, the separation SW, and the reattachment SW. In the first $x-z$ planes (on the wall), it can be noted that the wavelengths of the spanwise streaks on the ramp, colored by the Stanton number $St$ defined in Eq. 5, are about $4\sim 6$mm, which are consistent with the previous experimental observations (Roghelia et al., 2017a, b; Chuvakhov et al., 2017) and numerical results (Cao et al., 2021).

The aerothermodynamic characteristics are used to verify and validate the DNSs’ results, including the skin friction coefficient $C_{f}$, the surface pressure coefficient $C_{p}$, and the Stanton number $St$, which are defined as

$$C_{f}=\frac{\tau_{w}}{\frac{1}{2}\rho_{\infty}u_{\infty}^{2}},\quad C_{p}=\frac{p_{w}-p_{\infty}}{\frac{1}{2}\gamma Ma_{\infty}^{2}p_{\infty}},\quad St=\frac{q_{w}}{\rho_{\infty}u_{\infty}c_{p}(T_{aw}-T_{w})}$$

(5)

where $\tau_{w}$, $p_{w}$ and $q_{w}$ are the friction, pressure and heat flux on the wall, respectively. $c_{p}$ and $T_{aw}$ are the specific heat capacity and the adiabatic wall temperature, respectively.

The verifications include the time- and grid-convergences. For the time-convergences, as shown in Fig. 5(a), the non-dimensional simulation time $tu_{\infty}/L$ of M1 and M2 are both longer than $145$. It shows that the spanwise-averaged locations of the separation and reattachment points only oscillate slightly near $x_{S}/L\approx 0.5$ and $x_{R}/L\approx 1.3$ after $tu_{\infty}/L>70$, implying the separated DCR flows have been stabilized in both M1 and M2. For the grid-convergence, as shown in Fig. 5(b), there are only small discrepancy of the $C_{f}$ distributions between M1 and M2, implying the mesh resolution of M1 is suitable for the present simulations.

The validations include the comparisons of the present $C_{f}$, $C_{p}$, $S_{t}$distributions, and $u^{*}$, $T^{*}$ profiles with the accepted theoretical, numerical and experimental results. The laminar boundary layer before separation is self-similar satisfying the Blasius theoretical solutions (White & Majdalani, 2006). Fig. 5(b) shows the theoretical $C_{f}$ distribution, and Fig. 5(c) shows the theoretical normal profiles of $u^{*}$ and $T^{*}$ at $x/L=0.36$. It is noted that both the present $C_{f}$ distribution and normal profiles agree well with the theoretical solutions, validating the simulations’ accuracy. Furthermore, as shown in Fig. 5(d), the present $C_{p}$ and $St$ distributions are also compared with the previous experiment (Roghelia et al., 2017a; Chuvakhov et al., 2017) and DNS (Cao et al., 2021) results, and the excellent agreements validate the present DNSs.

As mentioned in Sec. 2, “adding the viscosity” and “filling the corner” are the key operations to realize the thought experiment. In terms of “adding the viscosity”, the implement method is changing the governing equations from Euler to Navier-Stokes. In terms of “filling the corner”, since it requires to fill the corner macroscopical continuously (fill the corner atom by atom, microscopically), the strict operation is to manipulate the flow fields with an enormous number of steps, the number of which is the order of magnitude of the Avogadro number. Obviously, it is impossible in numerical simulation, and the equivalent implement is to replace the filling mode from ‘macroscopically continuous’ to ‘macroscopically discrete’. The specific strategy is as follows.

For Route I (“adding the viscosity” $\rightarrow$ “filling the corner”), the simulation process can be denoted as D${}^{\text{I}}_{0}$ $\rightarrow$ C${}_{1}^{\text{I}}$ $\rightarrow$ C${}_{2}^{\text{I}}$ $\rightarrow$ … $\rightarrow$ C${}_{\text{N}}^{\text{I}}$, where the superscript ‘I’ denotes Route I, and the subscript ‘N’ represents the discrete number. Thus, D${}^{\text{I}}_{0}$ is the stable separated DCR flow in step 2, and C${}_{\text{N}}^{\text{I}}$ is the stable CCR flow after ‘N’ times of corner filling in step 3.

For Route II (“filling the corner” $\rightarrow$ “adding the viscosity”), the simulation process is denoted as C${}_{\text{N}}^{\text{II,inv}}$ $\rightarrow$ C${}_{\text{N}}^{\text{II}}$ ,where C${}_{\text{N}}^{\text{II,inv}}$ is the inviscid CCR flow in step 2, and C${}_{\text{N}}^{\text{II}}$ is the stale viscous flow in step 3, with the expectation of maintaining attached.

For D${}^{\text{I}}_{0}$ $\rightarrow$ C${}_{1}^{\text{I}}$ $\rightarrow$ C${}_{2}^{\text{I}}$ $\rightarrow$ … $\rightarrow$ C${}_{\text{N}}^{\text{I}}$ (Route I), we set N $=$ 4, and C₁, C₂, C₃, C₄ are four CCR configurations with curvature radiuses $R\approx 174.7,189.9,220.3,227.9$mm ($MO=23,25,29,30$mm), respectively. The four configurations’ meshes are set with the same resolution of M1 ($1020\times 250\times 200$). Three characteristic locations, $x_{S}$, $x_{R}$, and $x_{C_{f,min}}$, are considered to quantify the flow states variations, where $x_{C_{f,min}}$ denotes the locations of minimal $C_{f}$ values in attached states. During process D${}^{\text{I}}_{0}$ $\rightarrow$ C${}_{1}^{\text{I}}$ $\rightarrow$ C${}_{2}^{\text{I}}$ $\rightarrow$ C${}_{\text{3}}^{\text{I}}$ $\rightarrow$ C${}_{\text{4}}^{\text{I}}$, the evolutionary histories of spanwise-averaged $x_{S}$ and $x_{R}$ (red dashed lines) in C₁, C₁, C₃, and C₄ are shown in Fig.7(a), 7(b), 7(c), and 7(d), respectively, and the non-dimensional time $tu_{\infty}/L>100$ (see the upper red $x$-axes) for all cases. Obviously, $x_{S}$ and $x_{R}$ in C₁, C₂, and C₃ are all maintaining two branches for $tu_{\infty}/L>200$, implying C${}_{1}^{\text{I}}$, C${}_{2}^{\text{I}}$, and C${}_{\text{3}}^{\text{I}}$ can all be stabilized at separated states. However, for C₄ (the largest curvature radius), $x_{S}$ and $x_{R}$ ultimately merge into one branch (the $x_{C_{f,min}}$ branch— the blue solid line) at $tu_{\infty}/L\approx 100$, implying C${}_{\text{4}}^{\text{I}}$ ends up in the attached state. Flow fields of C${}_{1}^{\text{I}}$, C${}_{2}^{\text{I}}$, C${}_{\text{3}}^{\text{I}}$, and C${}_{\text{4}}^{\text{I}}$ are also shown in Fig. 9(a), 9(c), 9(e), and 9(b), respectively. It is noted that, on the ramp, the streaks in CCR are weaker than those in DCR.

For C${}_{\text{N}}^{\text{II,inv}}$ $\rightarrow$ C${}_{\text{N}}^{\text{II}}$ (Route II), the flows in Step 2 (inviscid flows in C${}_{\text{N}}$, N=4) are simulated with Euler equations, which are shown in Fig. 8(a) and denoted as C${}_{\text{1}}^{\text{II,inv}}$, C${}_{\text{2}}^{\text{II,inv}}$, C${}_{\text{3}}^{\text{II,inv}}$, and C${}_{\text{4}}^{\text{II,inv}}$, respectively. A series of compression waves distribute on the curved wall $\overline{MN}$ regions, and the $C_{p}$ distributions on $\overline{MN}$ satisfy the Prandtl-Meyer relations, as shown in Fig. 8(b), validating the present simulations.

To accomplish Step 3 in Route II (adding viscosity), Euler equations are replaced with Navier-Stokes form, with the expectation of obtaining the stable attached states. The processes are denoted as C${}_{\text{1}}^{\text{II,inv}}$ $\rightarrow$ C${}_{\text{1}}^{\text{II}}$, C${}_{\text{2}}^{\text{II,inv}}$ $\rightarrow$ C${}_{\text{2}}^{\text{II}}$, C${}_{\text{3}}^{\text{II,inv}}$ $\rightarrow$ C${}_{\text{3}}^{\text{II}}$, and C${}_{\text{4}}^{\text{II,inv}}$ $\rightarrow$ C${}_{\text{4}}^{\text{II}}$, whose evolutionary histories of $x_{C_{f,min}}$ are respectively shown in Fig.7(a), 7(b), 7(c), and 7(d) (non-dimensional time $tu_{\infty}/L>35$ for all cases, see the lower blue $x$-axes). $x_{C_{f,min}}$-variations in C₂, C₃, and C₄ are similar: being of only one branch, first moving downstream and then stabilizing after about $tu_{\infty}/L\approx 10$. This behavior implies C${}_{\text{2}}^{\text{II}}$, C${}_{\text{3}}^{\text{II}}$, and C${}_{\text{4}}^{\text{II}}$ all ultimately stabilize in attached states. However, in C${}_{\text{1}}$ (the smallest curvature radius), the $x_{C_{f,min}}$ bifurcates into two branches at $tu_{\infty}/L\approx 3$. The lower and upper branches repectively converge to the $x_{S}$ and $x_{R}$ branches obtained from D${}^{\text{I}}_{0}$ $\rightarrow$ C${}_{1}^{\text{I}}$, implying C${}_{\text{1}}^{\text{II}}$ finally stabilizes in the separated state. Flow fields of C${}_{1}^{\text{II}}$, C${}_{2}^{\text{II}}$, C${}_{\text{3}}^{\text{II}}$, and C${}_{\text{4}}^{\text{II}}$ are also shown in Fig. 9(a), 9(d), 9(f), and 9(b), respectively.