Hypersonic shock wave and turbulent boundary layer interaction in a sharp cone/flare model

Fulin TONG, Junyi DUAN, Jiang LAI, Dong SUN, Xianxu YUAN,*

a State Key Laboratory of Aerodynamics, China Aerodynamics Research and Development Center, Mianyang 621000, China

b Computational Aerodynamics Institute, China Aerodynamics Research and Development Center, Mianyang 621000, China

c LHD, Institute of Mechanics, Chinese Academy of Sciences, Beijing 100190, China

d School of Engineering Science, University of Chinese Academy of Sciences, Beijing 100049, China

KEYWORDS Shock wave;Hypersonic turbulent boundary layer;Heat transfer;Skin friction decomposition;Direct numerical simulation

Abstract A direct numerical simulation of hypersonic Shock wave and Turbulent Boundary Layer Interaction(STBLI) at Mach 6.0 on a sharp 7° half-angle circular cone/flare configuration at zero angle of attack is performed.The flare angle is 34°and the momentum thickness Reynolds number based on the incoming turbulent boundary layer on the sharp circular cone is Reθ = 2506.It is found that the mean flow is separated and the separation bubble occurring near the corner exhibits unsteadiness.The Reynolds analogy factor changes dramatically across the interaction, and varies between 1.06 and 1.27 in the downstream region, while the QP85 scaling factor has a nearly constant value of 0.5 across the interaction.The evolution of the reattached boundary layer is characterized in terms of the mean profiles,the Reynolds stress components,the anisotropy tensor and the turbulence kinetic energy.It is argued that the recovery is incomplete and the near-wall asymptotic behavior does not occur for the hypersonic interaction.In addition,mean skin friction decomposition in an axisymmetric turbulent boundary layer is carried out for the first time.Downstream of the interaction, the contributions of transverse curvature and body divergence are negligible,whereas the positive contribution associated with the turbulence kinetic energy production and the negative spatial-growth contribution are dominant.Based on scale decomposition, the positive contribution is further divided into terms with different spanwise length scales.The negative contribution is analyzed by comparing the convective term, the streamwise-heterogeneity term and the pressure gradient term.

1.Introduction

Understanding the physical phenomena associated with Shock wave and Turbulent Boundary Layer Interaction(STBLI)is of practical importance in the aeronautical and aerospace industries, because the shock-induced fluctuations of pressure and heat transfer are particularly severe, and are often responsible for thermal and structural failure of high-speed aircraft.Comprehensive review articles on the interaction have been published by Green,1Dolling2and Gaitonde,3who extensively discuss the physics as well as the computational and experimental methods used to address the problem.Depending on the complexity of the interaction, Gaitonde3has divided the canonical configurations used in STBLI into seven categories:impinging shock, compression ramp, cylinder flare or double cone, swept ramp (or corner), single sharp fin, double fin and internal flow.

There exists a wide body of experiments and numerical investigations by means of Direct Numerical Simulation(DNS)and Large-Eddy Simulation (LES) in the published literature for the nominally two-dimensional configurations,impinging shocks and compression ramps.Significant advances in the fundamental physics of STBLIs have been made over the past few decades, in areas such as the source of low-frequency unsteadiness,4–7turbulence amplification8,9and Go¨rtler-like vortices.10,11For example,as reviewed by Clemens and Narayanaswamy,4low-frequency unsteadiness in the interaction region is driven either by the fluctuations in the upstream Turbulent Boundary Layer(TBL)or by the intrinsic instability of the downstream separated flow, which is dependent on the size of the separated flow.They concluded that the downstream mechanism is mainly dominant for large flow separation, while the upstream and downstream mechanisms are both present in weakly separated flows.Fang et al.9reviewed the mechanism of turbulence amplification in STBLIs and suggested that turbulence amplification in the upstream part of the interaction zone is caused by the interaction between the deceleration of the mean flow and streamwise velocity fluctuations, rather than the free shear layer.Using Dynamic Mode Decomposition (DMD), Priebe et al.10found that the reconstructed flow using the five low-frequency DMD modes exhibited similarities with the counter-rotating Go¨rtlerlike vortices.Recently, the presence of Go¨rtler-like vortices in impinging shock wave and turbulent boundary layer interactions was experimentally demonstrated for the first time by Zhuang et al.11They discussed the absence of these vortices in previous studies, and proposed a possible mechanism for their generation in the reflected interactions.Despite remarkable progress, the majority of previous research has focused on supersonic STBLIs, and a comprehensive understanding of hypersonic interactions requires further study.

Until now,due to the difficulty of reliable and accurate turbulence measurements at high Mach number, there have been very limited experimental studies of hypersonic STBLIs including both mean and statistical turbulence quantities.An example is the particularly important measurement of turbulence fluctuations across the interaction region taken in the early experiment by Mikulla and Horstman12for two Mach 7 axisymmetric reflected interactions.Comparison of the separated and attached flows showed that a coupling between turbulent energy and separation bubble unsteadiness was responsible for substantial differences in turbulence lifetimes.Later, similar experiments were performed at Mach 8 by Bookey et al.13on an 8° compression corner, at Mach 7 by Schrijer et al.14on a double ramp flow and, more recently, at Mach 10 by Brooks et al.15on a hollow cylinder flare.These experiments demonstrate a strong amplification throughout the interaction.The recent particle image velocimetry data acquired by Schreyer et al.16at Mach 7.2 for 8° and 33° compression corners show a different amplification process for the fluctuation components, where the peak locations are further moved away from the wall in the strong interaction.

Recently,a preliminary attempt to perform a DNS analysis of hypersonic STBLI was carried out by Priebe and Martin,17who analyzed the evolution of turbulence on an 8° compression ramp at Ma∞= 7.2 and Reθ= 3500.In their attached STBLI, the Reynolds stress components were found to be amplified by factors of 1.8–2.5, and the turbulence became more isotropic in the interaction region.A quantitative analysis of the heat transfer showed that the commonly used Reynolds analogy did not apply in the interaction, and the heat transfer scaled better with pressure, rather than skin friction.

The objective of this paper is to systematically investigate hypersonic STBLI on a sharp cone/flare by means of DNS.To the best of the authors’knowledge, no other DNS studies of hypersonic interactions on the axisymmetric configuration have been reported so far.The geometry and the inflow conditions model used in the present study closely approximate recent experiments on a 7° half-angle circular cone/flare in the Air Force Research Laboratory(AFRL)Mach-6 Ludwieg tube by Running et al.18The flare angle is set to 34° and the mean flow in the interaction is separated.In this study, we comprehensively analyze the recovery of the reattached boundary layer and the generation of mean skin friction and heat transfer in the downstream region.

2.Numerical methodology

2.1.Governing equations

We use a high-order finite-difference DNS code, OpenCFDSC, to solve the three-dimensional compressible conservative Navier-Stokes equations in conical coordinates without any modeling.This code has been extensively validated for supersonic and hypersonic flows including STBLIs over a compression ramp19and incident shock interaction over a wavy-wall.20In this study,the governing equations are non-dimensionalized by the inflow parameters and the reference length scale(in millimeters), and are written as

with T∞being the inflow static temperature and Ts=110.4 K.Unless otherwise stated, the subscript‘‘∞”stands for a variable in the inflow freestream.

2.2.Numerical scheme

In the DNS code, the inviscid terms are calculated using the WENO-SYMBO scheme21with limiters22and the Steger-Warming vector flux splitting method.WENO-SYMBO is a bandwidth-optimized, 4th-order Weighted Essentially Non-Oscillatory (WENO) scheme that uses a symmetric collection of candidate stencils in the linear part of the original WENO method.Employing an absolute limiter in the smoothness measurement and a relative limiter based on the total variation, together, to further reduce numerical dissipation and to stabilize the simulation,WENO-SYMBO can properly resolve turbulence in the smooth region, while shock waves are still well captured.The viscous terms are discretized by an 8thorder central difference scheme, and the 3rd-order Runge-Kutta method23is used to perform the time integration.

2.3.Computation overview

Fig.1 shows the sharp cone/flare model used in the present simulation,which closely approximates the recent experiments conducted in the AFRL Ludwieg tube by Running et al.18An instantaneous iso-surface of density gradient colored by wallnormal distance is included.In Fig.1, the variable zois the location of the domain inlet, zaand zbdenote the start and end points of the blowing and suction region, zrefis the reference location in the upstream fully developed TBL.The sharp circular cone has a length of 609 mm and a half vertex angle of 7°, followed by the 34° flare with a length of 76 mm.The approach flow has a Mach number of Ma∞= 6.0, a unit inflow Reynolds number of Re∞=10.8×103/mm and a static temperature of T∞= 65 K.Following the DNS by Li et al.20and Sivasubramanian and Fasel,24the computational domain(see Fig.1) does not include the tip of the cone, which ranges from z = 156 mm to z = 675 mm.The domain height in the wall-normal direction is set to be 40 mm, and the domain size in the azimuthal direction corresponds to a span of φ=11.5°.Similar to the laminar-to-turbulent transition method of Pirozzoli et al.,25a region of strong random blowing and suction disturbances is enforced along the wall ranging from z = 180 mm to z = 185 mm.As a result, the inflow laminar boundary layer is triggered, and fully turbulent conditions are produced at the reference location zref= 580 mm, just upstream of the 34°flare.At zref,the boundary layer thickness is estimated to be δ=5.88 mm and the friction Reynolds number and the momentum thickness Reynolds number are Reτ= 343 and Reθ= 2506, respectively.

The computation grid and boundary conditions are illustrated in Fig.2.The grid is plotted at intervals of every ten and five points in the x and y directions, respectively.Unless otherwise specified, x and y denote the streamwise direction along the wall and the wall-normal direction,respectively,with usand unthe corresponding velocity components.In the x direction, a total number of 4819 points is used to discretize the domain, which is bounded by a laminar boundary layer at the domain inlet and a supersonic outflow boundary at the domain outlet.Specifically, 2820 points are equally distributed in the interaction region 576 mm 660 mm,49 points are located with progressively coarsening grid spacing to inhibit the reflection of disturbances.In the following analysis, the superscript‘‘+”refers to normalization with local wall units taken at zref.The profile of the imposed laminar boundary layer is determined from an auxiliary DNS of the sharp cone under the same inflow conditions.In the y direction, 319 points are clustered toward the wall, ensuring that there are 190 points located inside the boundary layer.Thus, the grid resolutions at the first point from the wall and at the edge of the boundary layer are Δy+= 0.59 and Δy+= 4.24, respectively, comparable with those in recent DNS studies on a hypersonic compression ramp flow by Priebe and Martin.17At the bottom boundary, the wall is modeled with no-slip isothermal boundary conditions, and the wall temperature is set to be Tw= 296 K, corresponding to a cold wall with wall-to-recovery temperature ratio Tw/Tr= 0.62.The velocity disturbances in the wall blowing and suction region(see Fig.1, za< z < zb) are similar to those in previous simulations by Li et al.,20,26where random numbers between –β and β are spatially randomly distributed and kept temporally constant.In order to greatly accelerate the transition process over the sharp cone, we choose β = 0.4, much larger than the value used by Li et al.20In addition,as suggested by Pirozzoli et al.,27nonreflecting boundary conditions are used at the upper boundary to inhibit the reflection of spurious disturbances.In the φ direction, periodic boundary conditions are enforced due to the homogeneity in the azimuthal direction,with 360 points being equally distributed.The grid spacing increases from Δrφ+= 2.31 at zrefto Δrφ+= 3.55 at z = 660 mm, caused by the rapid increase of the radius of the 34° flare angle.

Fig.1 Sketch of sharp cone/flare model and computational domain.

Fig.2 A two-dimensional sketch of computational grid and boundary conditions.

We performed a sensitivity study to assess the chosen grid resolution and the domain extent in the azimuthal direction.The selected mesh (Nx× Ny× Nφ= 4819 × 319 × 360) is labeled as Grid-B, and another three grids are labeled as Grid-A, Grid-C and Grid-D.Grid-A is generated by coarsening 30%of the grid spacing in both the x and φ directions,corresponding to Nx×Ny×Nφ=3917×319×240.In contrast,from Grid-B to Grid-C,the grid is refined by 30%in both the x and φ directions,so Nx×Ny×Nφ=5719×319×480.For Grid-D, the azimuthal domain extent only is increased by 50%, so Nx× Ny× Nφ= 4819 × 319 × 540.Note that the variation of the grid resolution in the x direction is applied in the interaction region, whereas the grid resolutions in the transition region and the buffer region are not changed.Fig.3 shows a comparison of mean wall quantities for Grid-A to Grid-D.In the results that follow, the mean refers to the average in time and in the azimuthal direction.It is clear that the curves for mean wall pressure pw(see Fig.3(a)) agree very well,and the streamwise locations of mean separation and reattachment point are insensitive to the variations of grid resolution and domain size (see Fig.3(c)).Slight variations are found only in the mean skin fiction coefficient Cfand heat transfer coefficient Chon the flare z > 620 mm.The results with Grid-A exhibit a deviation of approximately 10% from those with Grid-B, while the relative differences between the other three grids are confined within 5%, indicating that Grid-B is sufficient to obtain converged results and the selected azimuthal domain size is reasonable.

The adequacy of the adopted azimuthal domain extent is further assessed by analyzing the two-point autocorrelation function in the φ direction, which is defined as25,27

Fig.3 Grid and domain sensitivity study with respect to mean wall quantities evolution.

in which the bar denoting averages in time and in the φ direction,and sφand s representing the spacing and the coordinate in the azimuthal direction, respectively.Fig.4 shows the twopoint autocorrelations for the fluctuations of density, temperature and velocity components at two wall-normal locations in the interaction region, where y+= 6.7 and y+= 85 at z=617 mm.Both figures show that all the correlations decay rapidly and the fluctuations are decorrelated as the spacing is increased over half of the azimuthal width Lφ/2, confirming that the azimuthal domain used in the present simulation is wide enough and the turbulence fluctuations are not inhibited.

Fig.4 Distributions of two-point autocorrelation coefficient at different wall-normal locations in the interaction region z = 617 mm.

3.Results and discussion

3.1.Upstream turbulent boundary layer

Fig.5 Instantaneous numerical schlieren visualization of upstream TBL.

An instantaneous numerical schlieren visualization of the upstream TBL is shown in Fig.5, with the domain rotated 7° clockwise.The variable, defined as ds= 0.8exp[-10(|?ρ|-|?ρ|min)/ (|?ρ|max-|?ρ|min)], is a nonlinear function of the magnitude of the density gradient |?ρ|, often used in previous DNS studies of STBLIs.17,28With this transformation, the dark regions in the figure denote large density gradients in the flow field; the visualization exhibits a pattern consistent with observations in previous experiments.As observed by Prebie and Martin17in a DNS of a Mach 7.2 turbulent boundary layer and by Bookey et al.13in experiments on Mach 8 flows, sharp interfaces are clearly observed at the boundary layer edge, which separates the turbulent structures from the outer irrotational fluid, and the inclined turbulent bulges with high intermittency are highlighted inside the boundary layer.

Fig.6 Turbulence statistics obtained at the reference location zref.

Fig.6 shows a quantitative comparison of turbulence statistics taken at the reference location zrefwith previous experimental and numerical results.In Fig.6(a), the van Direst transformed mean streamwise velocity U+vdprofile agrees very well with the DNS data obtained by Priebe and Martin17at Reθ= 3350, except for notable differences in the wake region due to the different Reynolds number used.The classic law of the wall is only satisfied in the region y+< 5, and a narrow logarithmic region is observed for 40 < y+< 80, where the von Karman and log law constants are set to 0.41 and 5.8,respectively.Fig.6(b) compares the profiles of the densityscaled Reynolds stress components Rijwith previous compressible DNS data at similar Reynolds numbers.In Fig.6(b),uτis the friction velocity; δ is the reference boundary layer thickness.A satisfactory agreement with the numerical results reported by Pirozzoli et al.29in a spatially developing boundary layer at Mach 2 and Reτ= 278–358, and by Subbareddy and Candler30in a Mach 6 cold-wall boundary layer flow at Reθ≈ 6000, is obtained, confirming that the generation method of the incoming turbulent boundary layer used in the present study is reliable.

3.2.Instantaneous and mean flow fields

To give a qualitative overview of the interaction region,Fig.7 and Fig.8 show contours of instantaneous density gradient and mean pressure, respectively.In Fig.7, the foot of the strong shock penetrates into the boundary layer,and the thickness of the disturbed boundary layer decreases significantly,with a slow increase in the recovery region.Downstream of the interaction, the pressure is strongly amplified and the density gradient becomes much steeper (as reflected by the darker region at z>610 mm),due to the strong adverse pressure gradient caused by the large deflection angle and the high inflow Mach number.As is observed in Fig.8, a small region of separation bubble occurs at 605 mm < z < 610 mm.

Computed instantaneous streamwise velocity fields in the zr mid-plane are shown in Fig.9,which shows two uncorrelated samples of the DNS flow field.The pink dash-dot line shows where the streamwise velocity u=0,and the mean separation and reattachment points are denoted by symbols S and R,respectively.It is seen from Fig.9(a)that the separation bubble becomes much larger than the mean separation region, with significant reversed flows observed downstream of the reattachment point, whereas the bubble in Fig.9(b) is very small,with only a few locations with reversed flow concentrated below the mean separation region.This demonstrates the unsteady motion of the separated flow.Moreover, the variation in instantaneous bubble size (not shown here) reveals a breathing motion of the separation bubble, as previously observed by Priebe and Martin31in their DNS studies of a 24° compression ramp in Mach 2.9 flow.

Fig.7 Contours of instantaneous density gradient with five selected streamwise locations labeld as R1-R5.

Fig.8 Contours of mean pressure superimposed with in-plane streamlines(Inset is the enlargement of the separation region).

Fig.9 Contours of two instantaneous streamwise velocity fields in the z-r mid-plane with mean streamwise velocity u=0 isoline in pink.

To characterize the relation between mean wall quantities across the interaction zone, we analyze the distributions of the Reynolds analogy32and the QP85 scaling33factors,respectively.Here, the Reynolds analogy factor is defined as RAF = 2Ch/Cf≈1, and the QP85 scaling factor is defined as QP85 = (Ch/Chu)×(pwu/pw)0.85≈1, where pwuand Chuare the mean wall pressure and heat transfer coefficient,respectively, at zref.The Reynolds analogy distribution is shown in Fig.10(a) as a function of z.Upstream of the interaction,RAF attains a nearly constant value of 1.07,agreeing well with what was previously shown by Roy and Blottner,32who found RAF = 0.9–1.3 in hypersonic zero-pressure-gradient boundary layers.In the separation region, RAF varies significantly and the Reynolds analogy is not satisfied, as a consequence of the small negative Cf.Past the reattachment point(z > 610 mm), RAF rapidly decreases, attaining a minimum value of RAF=1.06 at z=628 mm.Subsequently,it follows a slight increase in the recovery region, with RAF = 1.27 at z = 660 mm.In Fig.10(b), it is seen that QP85 is nearly constant in the upstream TBL,where QP85 ≈1.0,consistent with the experimental results of Murray et al.,34whereas the QP85 scaling is not valid in the separation region,with QP85 varying between 0.6 and 1.8.However, it is worth noting that QP85 relaxes towards a nearly constant value (QP85 ≈0.5) in the recovery region at z > 625 mm, unlike in the DNS results of Priebe and Martin,17who found that QP85 monotonously decreased downstream of the interaction in a hypersonic compression ramp flow.This discrepancy is mainly ascribed to the different distribution of mean wall pressure.In our simulation,pwis nearly constant in the downstream region, while the results of Priebe and Martin17exhibited a continuous increase.

3.3.Evolution of boundary layer

To obtain a better understanding of the recovery of the reattached boundary layer, we further compare turbulence statistics in terms of mean velocity and temperature profiles,Reynolds stress tensor components,anisotropy invariant maps and Turbulent Kinetic Energy (TKE), at five selected streamwise locations on the flare, denoted by R1–R5as shown in Fig.7.The first location R1is at the mean reattachment point z=610 mm,while R2–R5are located at z=617 mm,626 mm,635 mm and 650 mm, respectively.In the following analysis,the reference boundary layer thickness δ is taken at zreffor direct comparison.

Fig.11(a) compares the mean streamwise velocity profiles at the five selected locations in the reattached boundary layer with that at zrefin the upstream boundary layer, providing a general view of the recovery of the velocity profiles downstream of the interaction.It is seen that the upstream velocity profile is in good agreement with the DNS data of Priebe and Martin17at Mach 7.2 flow, characterized by typical TBL behavior.At R1, very close to the reattachment point, the velocity significantly decreases for y/δ < 0.5, and the profile becomes much less full, compared to that of zref.For R2–R5,we observe that only the velocities in the small inner region are fully recovered, whereas most parts of the downstream profiles strongly deviate from the upstream profile for y/δ > 0.1.For instance, the velocities at y/δ = 0.1 and 0.3 are about 0.60U∞and 0.63U∞, respectively, for R5, which are much lower than 0.72U∞and 0.84U∞for zref.It is suggested that the recovery process of the reattached boundary layer on the flare is not fully completed.Note that as the wallnormal location moves further away from the wall (above the shock wave), the velocities at R3–R5sharply increase,finally approaching the freestream value.The corresponding temperature profiles are shown in Fig.11(b), where the DNS data at zrefand R5from the present study are compared with predicted values obtained from the modified Crocco-Busemann relation given by Walz.35According to the modified Croco-Busemann relation, symbols in Fig.11(b) represent the predicted temperatures at locations zrefand R5, respectively.As previously observed in supersonic turbulent boundary layers by Volpiani et al.36and Pirozzoli et al.,27the computed data at zrefmatches well the theoretical relation, indicating its reliability in characterizing the temperature–velocity relationship in hypersonic flows.We observe that the theoretical relation is not equally valid downstream of the interaction at R5, where significant deviations are clearly observed and the predicted values are much smaller than the DNS data over most of the boundary layer.

Fig.10 Distribution of mean wall quantities.

Fig.11 Profiles of mean streamwise velocity and mean temperature at various streamwise locations.

The profiles of the Reynolds stress components at R1-R5are plotted in Fig.12.The results at zrefare also included for comparison.As previously observed by Fang et al.9in a Mach 2.25 reflected interaction and by Loginov et al.37in compression ramp interactions, the figure highlights the large amplification of the Reynolds stress components downstream of the interaction.In particular, the components ～u′′u′′, ～v′′v′′, ～w′′w′′and-～u′′v′′at R1are amplified by factors of about 3.0, 16.0, 9.0 and 17.0, respectively, consistent with the PIV measurements of two hypersonic STBLIs by Schreyer et al.,16who found that the wall-normal fluctuations were more strongly amplified than the streamwise component.At R2–R5,all components follow a continuous decrease as the reattached boundary layer develops downstream.Note that the strong spikes in Fig.12(b) and (d) are the direct result of the shock wave, which is located far away from the boundary layer edge, see Fig.7.Importantly, we observe that the recovery of the streamwise component is much slower than those of the other three components.In Fig.12(a),from zrefto R1,the peak of ～u′′u′′moves outwards and appears at y/δ ≈0.1, as often observed previously in separated STBLIs.This behavior is presumably related to the reattachment of the detached shear layer developing above the separation bubble in the interaction zone.From R1to R2, the peak is shifted toward the wall instead,appearing in the very-near-wall region y/δ < 0.01; this can be ascribed to the decay of the reattached shear layer and the regeneration of the near-wall structures.From R3to R5,the magnitude gradually decreases, with the peak location changing slightly.It is seen that the profile at R5still significantly deviates from that of the upstream TBL,indicating that the recovery of the streamwise component is far from complete.Clearly, a different behavior is observed in Fig.12(b)-Fig.12(d), where the profiles of the other three components quickly recover,and the differences between R5and zrefare less significant.

Fig.12 Distribution of Reynolds stress tensor components at various locations.

Following the work of Lumley,38we show the Reynolds stress anisotropy invariant maps at various streamwise locations in Fig.13 to further analyze the amplification of the Reynolds stress components; the second IIband third IIIbinvariants of the anisotropy tensor bijare defined as

At zref,the anisotropy invariant map exhibits a typical characteristic of zero-pressure-gradient TBLs, as also observed in DNS studies by Sun et al.39In the near-wall region, twocomponent turbulence is attained and the anisotropy attains a maximum value in the buffer layer at y/δ ≈0.029 (y+≈10),which is in agreement with the numerical findings of Grilli et al.8and Pirozzoli et al.40As the wall-normal location moves outward, it is seen that the most turbulent state follows an axisymmetric expansion state, and the anisotropy of the flow in the outer region is further reduced,falling close to an isotropic state.At R1, a substantial reduction of the anisotropy of the velocity fluctuations in the near-wall region occurs, as highlighted in Fig.13(a).Turbulence in the immediate vicinity of the wall is characterized by a two-component axisymmetric state, and it closely resembles axisymmetric compression state in the inner part of the boundary layer y/δ < 0.025.As y+is further increased,the anisotropy of the flow begins to increase and turbulence develops along the limit line of the axisymmetric expression state,attaining its maximum anisotropy at y/δ ≈0.175.Fig.13(b)-Fig.13(d) shows a reversal tendency in the near-wall region at R2–R5, where the invariant trace is pulled toward the top right of the map, corresponding to a onecomponent turbulence state, suggesting that the anisotropy is significantly greater.This behavior is attributed to the regeneration of the elongated high-and low-speed streaks in the nearwall region.

Fig.14 shows the distributions of the density-scaled turbulence kinetic energy k+as a function of wall-normal distance y+in local wall units at various streamwise locations.Upstream of the interaction, the TKE at zrefattains a maximum value of k+= 5.2 at y+= 15, which is in close agreement with the computed values given by Patel et al.41for an incompressible TBL (k+= 4.5 at y+=15) and by Pirozzoli et al.27for a spatially developing supersonic TBL (k+=5.25 at y+=14.3).As expected, the TKE at R1is greatly amplified, indicating turbulence amplification after passing through the interaction.Moving downstream, the TKE decreases rapidly.It is also interesting to note that the TKE profile at zrefobeys the near-wall asymptotic behavior, k+≈0.2035y+2, but this asymptotic behavior is not satisfied in the near-wall region at R1–R5.Such a scenario shows essential differences from the findings of Pirozzoli et al.27and Tong et al.20in their DNS of supersonic STBLI flow.They argued that the asymptotic consistency was held at different constant values downstream of the interaction.

Fig.13 Reynolds stress anisotropy invariant maps at various streamwise locations.ISO:isotropic;AC:axisymmetric compression;2CA:two-component axisymmetric; AE: axisymmetric expansion; 2C: two components.

Fig.14 Profiles of turbulent kinetic energy at various streamwise locations with the red solid lines being the correlation k+ =0.2035y+2.

3.4.Mean skin friction decomposition

Useful information for understanding and predicting mean skin friction generation can be obtained using the Renard-Deck(RD)identity42to decompose the mean skin friction into physics-based contributions.We extend the method of Li et al.,43for the first time, to investigate the generation of the mean skin friction in hypersonic STBLI on a cone/flare.The mean skin friction decomposition in conical coordinates is expressed as

where Cf,V,Cf,T,Cf,Gand Cf,Rrepresent,respectively,the contributions of direct viscous dissipation, TKE production, thespatial growth of the flow, and the combined effect of transverse curvature and body divergence.Recalling the compressible RD identity proposed by Li et al.,43we highlight that the only difference in our decomposition is the transverse curvature and body divergence contribution Cf,R, due to the body radius r.The explicit expressions for the four terms are given as

Table 1 Fractal contributions to mean skin friction at various streamwise locations.

Fig.15 Pre-multiplied integrands of Cf,V/Cf as a function of y/δ at various streamwise locations.

Here, σxx, σyxand ?p/?x denote the normal viscous stress,the shear viscous stress and the pressure gradient in the x direction, respectively, and σθθ, σrrand σrzare the components of the viscous stress in conical coordinates.

Considering that the mean skin friction is significantly increased downstream of the interaction, as previously shown in Fig.3(b),we quantitatively compare the decomposed results at the five locations R1–R5in the reattached boundary layer with those at the reference location zrefin the upstream TBL to investigate the shock interaction on the skin friction generation.Table 1 lists the fractional contributions of the above four terms at different streamwise locations.The relative error,quantifying the deviation from the raw Cf, is defined as(Cf,V+ Cf,T+ Cf,G+ Cf,R– Cf) / Cf.As is apparent from Table 1, the relative errors are confined to within ±6.82%,demonstrating the high accuracy of this decomposition.In good agreement with the decomposition analysis by Fan et al.44of a spatially developing zero-pressure-gradient supersonic flat plate turbulent boundary layer,the mean skin friction at zrefis dominated by Cf,Vand Cf,T, which contribute approximately 49.09%and 42.57%,respectively,to Cf,whereas the contributions of Cf,Gand Cf,Rare rather small and negligible.However, the decomposed results at R1–R4show a completely different picture.Here, the positive contribution of Cf,Tplays the leading role, and rapidly decreases from 2994.16% at R1to 54.65% at R4.Consistent with recent studies of reflected interactions by Zhang et al.45and Tong et al.46, Cf,Gbecomes negative and contributes significantly to Cf, in particular –2945.62% at R1and –240.81% at R2.It is suggested that the small positive Cfdownstream of the interaction is mainly produced by the very large positive Cf,Tand the very large negative Cf,G.At R5, the negative Cf,Gbecomes positive and the contribution of Cf,Vbecomes dominant again, showing close similarities with the decomposition in the upstream TBL.It is also emphasized that the influence of the direct shock interaction on the contribution of Cf,Rcan be negligible, where Cf,R/Cfis very small in the downstream region.

To obtain a more intuitive understanding, we also show profiles of the pre-multiplied integrands of the significant components in the wall-normal direction at various streamwise locations, where the area below the curve is in proportion to the ratio Cf,V/Cf,Cf,T/Cfor Cf,G/Cf.The results at zrefare also included for reference.

In Fig.15(a), the pre-multiplied integrand of Cf,V/Cfat zrefpeaks at y/δ = 0.02 (y+= 6.8) and the contribution of Cf,Vmainly comes from the inner region y/δ<0.1.At R1,it is seen in Fig.15(b) that the profile shape is different, with the peak located at y+= 61.5 and most of Cf,Vis generated in the region 0.04 < y/δ < 0.4, indicative of the great importance of the outer region just after passing through the interaction.At R2–R5(see Fig.15(c)-Fig.15(f)),the profiles have recovered and shown the similarities with that of the upstream TBL,but the peak reappears in the very-near-wall region, being approximately fixed at a wall-normal location of y/δ = 0.003(y+= 1.1), which is independent of the relaxation process.

Regarding the pre-multiplied integrand of Cf,T/Cf, as displayed in Fig.16, a similar evolution process in the reattachment region is observed.As shown in Fig.16(a), the peak at zrefis located in the near-wall region at y/δ = 0.043 (y+=14.8), while a secondary peak is identified in the logarithmic region at y/δ = 0.31 (y+= 104).Such double-peak behavior is very similar to that seen in the decomposition analysis of Fan et al.44in zero-pressure-gradientflat plate turbulent boundary layers;they suggested that the secondary peak could be ascribed to the generation of large-scale energy-containing motions in the logarithmic and outer regions.Additionally,the computed secondary peak location y2+ndin the present study is consistent with the linear relationship y+2nd=0.43Reτ– 10.11 proposed by Fan et al.44In Fig.16(b)-Fig.16(c), it is seen that the profiles at R1and R2are characterized by a single outer peak,appearing at y/δ=0.16 and 0.20,respectively,and Cf,Tis mainly generated in the region y/δ=0.043–0.31.This is mainly related to the outer large-scale structures downstream of the interaction and will be quantitatively analyzed in the following.Going from R3to R5,as the reattached boundary layer has developed further downstream, the profiles shown in Fig.16(d)-Fig.16(f) relax toward the upstream doublepeaked distribution at zref.Note that the strong spikes in the figures are caused by the shock wave located above the boundary layer.Clearly, the outer peak decreases rapidly, with its location moving slightly outward,whereas the inner peak reoccurs at y/δ < 0.01 and significantly increases.It is seen from Fig.16(f) that the inner peak is close to that at zref, with its location being lower,while the outer peak is much smaller than the upstream value, suggesting that the relaxation process is incomplete inside the entire boundary layer.

Fig.16 Pre-multiplied integrands of Cf,T/Cf as a function of y/δ at various streamwise locations.

Similar to the analysis of Duan et al.,47we further break down the Cf,Tcontribution into various circumferential length scales.At a given y location, the cospectra of the Reynolds shear stress in Eq.(9) are defined as.

Correspondingly, the integrand function ksyφRD/Cfin the λ-y plane on logarithmic scales represents the contribution of a given circumferential wavelength λ at a given location y to Cf.

Fig.17 shows the integrand function ksyφRD/Cfat three streamwise locations,showing variation with both inner scales(λ+and y+).The spectra are normalized by the maximum value.At zref(see Fig.17(a)),the most significant contribution appears at y+=15.4 and λ+=87.1.In Fig.17(b),we observe that the peak location moves outward and the peak wavelength becomes much larger, so that y+= 54.9 and λ+=177.8 at R1.Such behavior is consistent with the above analysis of Fig.16(a)-Fig.16(b).At R5, Fig.17(c) shows the strongest contribution at y+= 3.3 and λ+= 25.9, suggesting that the Cf,Tcontribution is mainly associated with the small-scale fluctuations in the very-near-wall region.It is seen that the peak wavelength at R5is still much smaller than that at zref, also indicating that the recovery process is incomplete.Recalling the analysis of the pre-multiplied integrands in Fig.16, we see that the generation of Cf,Tin the reattached boundary layer is characterized by the transformation from large-scale turbulence in the outer region to small-scale turbulence in the nearwall region.

Fig.17 Local Cf,T contribution at three streamwise locations.

Fig.18 Pre-multiplied integrands of Cf,G/Cf as a function of y/δ at various streamwise locations.

According to Eq.(8), the contribution Cf,Gcan be decomposed as.

Here,Cf,G1,Cf,G2and Cf,G3denote the convective term,the streamwise-heterogeneity term and the pressure-gradient term,respectively.Therefore,the pre-multiplied integrand of Cf,G/Cfconsists of three corresponding parts: Cf,G1/Cf, Cf,G2/Cfand Cf,G3/Cf, as show in Fig.18.At R1, the positive Cf,G(see Fig.18(a)), confined to the outer region, is mainly related to the positive Cf,G1, whereas the magnitudes of the other two terms are rather small and can be neglected.However, the pre-multiplied integrands at R1and R2reveal a different scenario in Fig.18(b)-Fig.18(c).The curves of Cf,Gexhibit a sign switching, being negative in the near-wall region and positive in the outer region.Such behavior is the combined result of the very large positive Cf,G1and the very large negative Cf,G3, consistent with the recent decomposition analysis of Fan et al.48on adverse pressure gradient turbulent boundary layers developing on flat plates and airfoils.As seen from Fig.18(d)-Fig.18(f), Cf,G3rapidly decreases due to the reduction of the pressure gradient in the downstream region,and the contribution Cf,G1is again dominant over most of the boundary layer.

4.Conclusions

A DNS of a hypersonic sharp cone/flare configuration has been performed.The numerical study has focused on the fundamental mechanisms dictating the hypersonic interaction,including mean and instantaneous properties, turbulence evolution and the decomposition of mean skin friction.

(1) The simulation shows that the mean flow is separated and a substantial unsteadiness of the moving separation bubble is observed near the corner.The Reynolds analogy and QP85 scaling are not satisfied in the interaction zone.Downstream of the interaction, RAF is about 1.06–1.27, very close to the upstream value, whereas QP85 is maintained at a nearly constant value of 0.5,half of the upstream value.

(2) The study shows that the turbulence is strongly amplified, especially for the wall-normal component of the Reynolds stresses.The reduced turbulence anisotropy slowly increases as the reattached boundary layer relaxes downstream.The downstream TKE does not obey follow near-wall asymptotic behavior, which is different from the findings in previously studied supersonic interactions.

(3) The mean skin friction decomposition exhibits similarities with the reflected interaction, where the generation in the downstream region is dominated by the positive contribution related to TKE production and the negative contribution associated with spatial growth of the flow.The positive contribution is mainly linked to the counterbalance between the outer large-scale and inner small-scale turbulence.The convective term and the pressure gradient term contribute most to the negative spatial growth contribution.The combined effect of transverse curvature and body divergence on the generation is negligible.

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgements

This study was co-supported by the National Natural Science Foundation of China (Nos.11972356 and 91852203) and the National Key Research and Development Program of China(No.2019YFA0405300).