1. Introduction
Turbulent boundary layers in high-speed vehicles can be significantly affected by the presence of surface roughness. Although most flight systems are designed to have relatively smooth surfaces, roughness can still occur or develop for various reasons. For instance, surface topology may be altered by localised effects such as pitting, corrosion, spallation or contamination deposits (Kocher et al. Reference Kocher, Kreth, Schmisseur, LaLonde and Combs2022). Roughness can also extend across larger areas of the surface, especially where thermal protection systems (TPSs) are used. Tiled TPSs often exhibit seams at the interfaces between tiles, which interact with the incoming turbulent boundary layer, strongly impacting both drag and heat transfer (Ekoto et al. Reference Ekoto, Bowersox, Beutner and Goss2009). In contrast, ablative TPSs undergo material removal as a sacrificial layer, potentially creating regular or irregular roughness patterns (Peltier, Humble & Bowersox Reference Peltier, Humble and Bowersox2016; Wilder & Prabhu Reference Wilder and Prabhu2019). Due to the complexities of supersonic and hypersonic flows compared with incompressible flows (Candler Reference Candler2019), the additional impact of surface roughness has rarely been studied. While experimental campaigns have provided valuable insights into this topic (Bowersox Reference Bowersox2007), high-fidelity simulations remain limited, yet they could offer useful insights into critical aspects of the flow dynamics.
In the incompressible regime, the effects of roughness in wall-bounded turbulent flows have been extensively investigated through both experimental and numerical studies, as recently reviewed by Chung et al. (Reference Chung, Hutchins, Schultz and Flack2021)and Jiménez (Reference Jiménez2004). These studies have advanced our understanding of turbulent flows over rough surfaces and highlighted that predictive tools for total drag
$\tau _w$
and heat flux
$q_w$
per plane are available, primarily in the form of empirical correlations (Macdonald, Griffiths & Hall Reference Macdonald, Griffiths and Hall1998; Flack & Schultz Reference Flack and Schultz2010; Yang et al. Reference Yang, Stroh, Lee, Bagheri, Frohnapfel and Forooghi2023) and, to some extent, as more robust physics-based models (Yang et al. Reference Yang, Sadique, Mittal and Meneveau2016; Zhong, Hutchins & Chung Reference Zhong, Hutchins and Chung2023; Meneveau, Hutchins & Chung Reference Meneveau, Hutchins and Chung2024).
Characterising roughness relies on the relative importance of the different length scales involved, namely the boundary-layer thickness
$\delta _{99}$
, the roughness height
$k$
and the viscous length scale
$\delta _{\nu }=\nu _w/ u_{\tau }$
, where
$u_{\tau }=\sqrt {\tau _w/\rho _w}$
is the friction velocity, and
$\nu _w$
and
$\rho _w$
are the kinematic viscosity and density at the wall, respectively. At sufficiently high friction Reynolds numbers,
$Re_{\tau }=\delta _{99}/\delta _{\nu }$
, and roughness Reynolds numbers,
$k^+=k/\delta _{\nu }$
, the flow can be considered fully rough, where the drag per unit area becomes independent of the Reynolds number as pressure drag becomes the dominant contributor (Nikuradse Reference Nikuradse1933).
However, a quantity of engineering interest is often the relative increase in drag compared with a smooth wall, which remains Reynolds-number-dependent. For this reason, the Hama roughness function is preferred for quantifying the added drag, as it remains relatively independent of the Reynolds number. The Hama roughness function,
$\Delta U^+$
, describes the downward shift in the viscous-scaled mean velocity profile induced by roughness compared with a smooth wall. Its effectiveness as a measure of added drag is supported by the validity of outer-layer similarity.
The similarity of the outer flow implies that roughness affects only a small region near the wall, while turbulence in the outer layer remains independent of the wall condition (Jiménez Reference Jiménez2004), perceiving the differing surface only as a change in the mean drag (Flack, Schultz & Connelly Reference Flack, Schultz and Connelly2007). Most surface patterns exhibit outer-layer similarity, meaning that estimating the added drag translates into establishing a relationship between the geometrical characteristics of the roughness (e.g.
$k$
) and
$\Delta U^+$
, which generally varies for each surface pattern. In an attempt to have a semblance of universality, the equivalent sand-grain roughness
$k_s$
is defined as the roughness height that a sand-grain surface would need to produce the same
$\Delta U^+$
as the surface in question (Nikuradse Reference Nikuradse1933; Flack & Schultz Reference Flack and Schultz2014).
In the context of high-speed flows, it remains unclear whether theories developed for incompressible turbulence are applicable, and the literature on this topic is limited. This fact was emphasised in the literature survey by Bowersox (Reference Bowersox2007), who described the high-speed database available at the time as consisting solely of experimental studies. Addressing this gap through high-fidelity numerical simulations is the primary aim of this work, with the following open questions as key points of investigation.
First, compressibility transformations provide a framework for scaling the compressible rough-wall velocity profile to its incompressible counterpart, thereby recovering the outer-layer similarity described above. Several experimental studies have demonstrated that the theory of van Driest (Reference van Driest1951) effectively achieves Mach-number invariance when the wall is adiabatic (Ekoto et al. Reference Ekoto, Bowersox, Beutner and Goss2008; Williams et al. Reference Williams, Sahoo, Papageorge and Smits2021; Kocher et al. Reference Kocher, Kreth, Schmisseur, LaLonde and Combs2022). Recently, Modesti et al. (Reference Modesti, Sathyanarayana, Salvadore and Bernardini2022) performed direct numerical simulation (DNS) of supersonic diabatic turbulent channels, showing that more recent compressibility transformations (Volpiani et al. Reference Volpiani, Iyer, Pirozzoli and Larsson2020) can account for wall temperature effects when an equivalent roughness Reynolds number is used.
Another aspect highlighted in the survey by Bowersox (Reference Bowersox2007) is that the interaction between compressibility and wall roughness can produce shock and expansion waves generated by each roughness element. These waves traverse the boundary layer and extend into the free stream. This effect was observed in several experimental studies at Mach numbers from 2 to 2.9 (Latin & Bowersox Reference Latin and Bowersox2000; Ekoto et al. Reference Ekoto, Bowersox, Beutner and Goss2008; Kocher et al. Reference Kocher, Kreth, Schmisseur, LaLonde and Combs2022), where the wave structure was found to significantly impact both first- and second-order statistics, as well as to locally alter the wall shear stress. The geometries tested included three-dimensional (3-D) cubes, 2-D bars sand-grained surfaces and a diamond-like pattern, with the latter producing the most intense local distortions in the flow (Ekoto et al. Reference Ekoto, Bowersox, Beutner and Goss2008).
In this respect, it remains unclear whether, in certain configurations, compressibility effects may be strong enough to break outer-layer similarity. Peltier et al. (Reference Peltier, Humble and Bowersox2016) reached this conclusion when analysing the wave structure produced by flow over cross-hatch roughness at Mach 4.9, which consisted of a pattern of shocks and expansions that disturbed the entire boundary layer. They proposed a conceptual model suggesting that the lower part of the boundary layer is dominated by compression waves that enhance turbulence intensity, while expansion waves suppress it in the upper portion. This observation contrasts with the experimental study by Williams et al. (Reference Williams, Sahoo, Papageorge and Smits2021), who conducted experiments at hypersonic speeds (Mach 7.3) and observed a general agreement with similar studies in incompressible flows. It is notable that the friction Reynolds number,
$Re_{\tau }$
, in this study was much lower than in the experiments by Peltier et al. (Reference Peltier, Humble and Bowersox2016).
All previous experimental studies include a smooth-to-rough transition; however, almost none have assessed the adjustment length required to adapt to the new surface condition. In the incompressible literature, smooth-to-rough transitions (or vice versa) have been extensively studied, featuring the formation of an internal boundary layer,
$\delta _i$
, which develops within the incoming boundary layer and eventually merges with it (Elliott Reference Elliott1958; Antonia & Luxton Reference Antonia and Luxton1971; Rouhi, Chung & Hutchins Reference Rouhi, Chung and Hutchins2019; Li et al. Reference Li, de Silva, Chung, Pullin, Marusic and Hutchins2022). The experimental study by Kocher et al. (Reference Kocher, Kreth, Schmisseur, LaLonde and Combs2022) is, to our knowledge, the only one that attempts to characterise smooth-to-rough transitions in a supersonic boundary layer, although it considers only a limited number of streamwise locations. They observed a substantial change in turbulent statistics at least 19.5
$\delta _{ref}$
downstream of the transition, where
$\delta _{ref}$
is the reference boundary-layer thickness of the incoming smooth-wall boundary layer. Providing a deeper understanding of smooth-to-rough transitions in compressible turbulence is one of the key objectives of the present study.
In this work, we perform DNS of zero-pressure-gradient turbulent boundary layers at free-stream Mach numbers of
$M_{\infty }=[0.3,2]$
and friction Reynolds numbers up to
$Re_\tau \approx 1700$
over both smooth and rough walls, aiming to gain valuable insights into the flow physics. We focus on a single roughness geometry consisting of cubical roughness elements, preceded by an incoming smooth-wall boundary layer. This configuration allows us to study, for the first time, the streamwise development of a supersonic boundary layer adjusting to a rough surface, with the goal of assessing similarities and differences compared with the subsonic case. Compared with previous studies, the novelty of this analysis lies in the precise description of the flow fields, enabling an in-depth examination of all turbulent scales and a direct comparison with a subsonic case, which was impractical in earlier experimental studies. Turbulence and thermal statistics are then analysed downstream of the surface transition to determine whether any outer-layer similarity can be observed when compared with reference smooth-wall cases.
2. Methodology
We solve the compressible Navier–Stokes equations for a viscous, heat-conducting gas:

where
$\rho$
is the density,
$u_i$
denotes the velocity component in the ith Cartesian direction (
$i=1,2,3$
),
$p$
is the thermodynamic pressure,
$E=c_v T+u_i u_i/2$
the total energy per unit mass and

are the viscous stress tensor and the heat flux vector, respectively. The molecular viscosity
$\mu$
is assumed to follow Sutherland’s law, with a reference free-stream temperature
$T_{\infty }=220.0$
K. The thermal conductivity
$k$
is related to the viscosity through the Prandtl number
$Pr=0.72$
,
$k=c_p \mu /Pr$
, where
$c_p$
is the specific heat at constant pressure. This model is complemented by the equation of state for a calorically perfect gas. The system of equations is solved on a Cartesian grid using the in-house code STREAmS (Bernardini et al. Reference Bernardini, Modesti, Salvadore and Pirozzoli2021, Reference Bernardini, Modesti, Salvadore, Sathyanarayana, Della Posta and Pirozzoli2023), which has been extensively validated in numerous canonical flow configurations (Bernardini, Pirozzoli & Grasso Reference Bernardini, Pirozzoli and Grasso2011; Bernardini & Pirozzoli Reference Bernardini and Pirozzoli2011; Modesti & Pirozzoli Reference Modesti and Pirozzoli2016; Cogo et al. Reference Cogo, Salvadore, Picano and Bernardini2022, Reference Cogo, Baù, Chinappi, Bernardini and Picano2023). Convective terms are discretised using sixth-order, energy-preserving schemes applied in shock-free regions, while a high-order shock capturing scheme (weighted essentially non-oscillatory) is applied when shock waves are identified by the Ducros sensor (Ducros et al. Reference Ducros, Ferrand, Nicoud, Weber, Darracq, Gacherieu and Poinsot1999). Viscous terms are discretised using a locally conservative formulation (De Vanna et al. Reference De Vanna, Benato, Picano and Benini2021) with second-order accuracy.
The present database is composed of two simulations featuring compressible turbulent boundary layers over rough surfaces at
$M_{\infty }=0.3$
(RH_M03) and
$M_{\infty }=2$
(RH_M2), alongside a supersonic smooth-wall case for comparison (SM_M2), (see table 1). For all cases, the friction Reynolds number at the inflow is
$Re_{\tau ,in}=\delta _{in}/\delta _{\nu ,in}=400$
, where
$\delta _{in}$
and
$\delta _{\nu ,in}$
are the boundary-layer thickness and viscous length scale at the inflow boundary.
Table 1. Summary of parameters for DNS study. Domain lengths
$L_x,L_y,L_z$
are given in terms of the inflow BL thickness
$\delta _{in}$
.


Figure 1. Schematic of the computational set-up for a turbulent boundary-layer flow over cubical roughness. (a) Overview of the computational domain (not to scale). (b) Wall-parallel arrangement of 3-D cubes in the rough portion. (c) Cross-stream arrangement of 3-D cubes. The roughness elements have height
$k$
and spacing
$2k$
.
The rough-wall simulations consist of three regions: an initial smooth-wall part dedicated to the development of a turbulent boundary layer, obtained through a recycling/rescaling procedure, where the recycling plane is placed at
$x=40 \delta _{in}$
. The second part starts at
$x=55\, \delta _{in}$
, ends at
$x=147\delta _{in}$
and consists of cubic elements of side
$k$
, which are representative of the structured roughness patterns forming over ablative surfaces (Modesti et al. Reference Modesti, Sathyanarayana, Salvadore and Bernardini2022). Finally, a small smooth-wall part is placed at the end of the domain until
$x=150 \delta _{in}$
as a buffer region before the outflow. This specific design has been chosen after a series of preliminary tests in smaller computational domains, and is a compromise between computational cost and the need to adequately resolve each region. In particular, in the smooth section the location of the recycling/rescaling procedure is placed at
$x=40 \delta _{in}$
to achieve a fully developed turbulent boundary layer before the onset of roughness, followed by an additional
$15 \delta _{in}$
units as a buffer region to avoid any upstream influence of roughness on the recycling plane. The rough portion is extended as far as possible in order to adequately study the boundary-layer adjustment to the new surface, but it transitions back to a smooth one
$3 \delta _{in}$
units before the outflow (hence up to
$147 \delta _{in}$
) in order not to interfere with the boundary conditions. We highlight that the extension of the rough portion is unprecedented compared with previous studies (Kocher et al. Reference Kocher, Kreth, Schmisseur, LaLonde and Combs2022), considerably reducing the effects of the smooth-to-rough transition in the analysis of the wall-normal turbulence statistics near the end of the domain. The roughness elements are spaced by a distance
$2k$
in the wall-parallel directions, as shown in figure 1. The complexity of the geometry is handled using a ghost-point-forcing immersed boundary method (Piquet, Roussel & Hadjadj Reference Piquet, Roussel and Hadjadj2016; De Vanna, Picano & Benini Reference De Vanna, Picano and Benini2020), which is used to enforce no-slip adiabatic boundary conditions on the solid wall. Following the mesh convergence study of Modesti et al. (Reference Modesti, Sathyanarayana, Salvadore and Bernardini2022), which simulated the same roughness geometry in turbulent channel flows, the number of computational nodes per element is
$20,40,20$
in the three directions
$x,y,z$
, respectively. Throughout this study, mean flow statistics are collected at selected streamwise locations, reported in table 2. For rough-wall cases, we selected a station
$x/\delta _{in}$
far from the surface transition to avoid non-equilibrium effects, whereas smooth-wall statistics are collected at a similar
$Re_{\tau }$
. Table 2 also reports the computed roughness Reynolds number
$k^+$
and the ratio
$\delta _{99}/k$
at the selected stations for rough cases, where we note that at these stations we expect a fully rough regime as
$k^+_s\gt 80$
. Here, the ratio
$k_s/k=1.84$
is computed as discussed in § 4.3.
Table 2. Boundary-layer properties at the selected stations. Grid spacings are given in wall units according to the selected station. The values of
$\Delta y^+_{min}$
and
$\Delta y^+_{max}$
refer to the wall-normal spacing at the wall and at the boundary-layer edge, respectively.

Throughout this study, we use the symbols
$u$
,
$v$
and
$w$
to denote the streamwise, wall-normal and spanwise velocity components and the decomposition of any variable is conducted using either the standard Reynolds decomposition (
$f=\bar {f}+f'$
) or the density-weighted (Favre) representation (
$f=\tilde {f}+f^{\prime \prime }$
), with
$\tilde {f}=\overline {\rho f}/\bar {\rho }$
. Here, mean quantities are computed using different samples in time, exceeding 500
$\delta _{in}/ u_{\infty }$
for all cases, and space, considering the periodic spanwise direction
$z$
and a moving average over a small window
$\lambda _x$
in the streamwise direction
$x$
. The latter corresponds to four roughness periods, or equivalently,
$\lambda _x=1.44\delta _{in}$
, which is approximately half of the mean boundary-layer thickness after the onset of roughness. The presence of computational nodes inside the solid domain is taken into account by considering the intrinsic average of mean flow statistics only in the fluid plane area
$S_f$
, such that
$\bar {f} = 1/S_f \int _{fluid} f {\rm d}S$
.
3. Instantaneous flow
We start by providing an overview of the instantaneous flow organisation to build a qualitative understanding of the main flow features, particularly at the smooth-to-rough surface transition. First, we inspect the instantaneous density field of the supersonic case in figure 2 in a longitudinal and wall-parallel plane. Before the roughness, we find a canonical smooth-wall organisation which is abruptly disrupted by the occurrence of a shock wave caused by the first row of roughness elements. In the longitudinal plane, behind the shock wave, we also note a series of expansion and compression waves possibly caused by the roughness pattern below.

Figure 2. Instantaneous density field
$\rho /\rho _{\infty}$
for flow case case RH_M2. The flow domain is excluded in a portion around the transition region to highlight the roughness elements.
Figure 3 shows the instantaneous contours of streamwise velocity for both subsonic, RH_M03, and supersonic, RH_M2, cases, indicating the development of the turbulent boundary layer and its influence into the free stream. In contrast to the streamwise evolution of case RH_M03, which seems to be weakly affected by the onset of surface roughness, supersonic case RH_M2 shows emerging features that are clearly related to the interaction between compressibility and the roughness. A strong compression is observed around
$x/\delta _{in}=55$
, where the roughness starts. The first compression wave is not directly emanated by the first row of cubical elements, but by a sudden adjustment of the boundary-layer thickness due to an upward shift forced by the new surface condition. After the first compression, we observe a series of smaller compression/expansion waves that protrude into the free stream. The level of acoustic disturbances emanated from the boundary layer tends to be lower downstream, but they are still present until the end of the domain (not shown).

Figure 3. Instantaneous contours of the streamwise velocity
$u/u_{\infty }$
on a wall-normal plane. Panel (a) shows the RH_M03 case; panel (b) shows the RH_M2 case.
Similar patterns can be observed in the instantaneous contours of the density gradient, visualised as a numerical schlieren in figure 4(a). Here, we also report the time-averaged version of the same quantity over a short period of
$75\delta _{in}/u_\infty$
in figure 4(b). This averaged flow better highlights the formation of the first strong compression, primarily caused by the vertical shift of the boundary layer, although a subsequent distortion is emanated from the first cubical element. We classify this first compression as a shock wave given that the propagation angle is estimated as
$31.8^{\circ }$
, while Mach waves have a theoretical angle of
$30^{\circ }$
at
$M_{\infty }=2$
. Figure 4(b) also shows the sequence of compression waves emanated from each element and propagating into the free stream, which is visible throughout the domain and is in agreement with experimental studies that reported similar wave patterns (Ekoto et al. Reference Ekoto, Bowersox, Beutner and Goss2008; Kocher et al. Reference Kocher, Kreth, Schmisseur, LaLonde and Combs2022). Finally, from the time-averaged numerical schlieren we have a first glimpse of the formation of an internal boundary layer, visible as a dark blue region progressively growing as we move downstream.

Figure 4. Instantaneous (a) and time-averaged (b) contours of the numerical schlieren
$exp(-| \nabla \bar {\rho }|)$
on a wall-normal plane for the case RH_M2. In (b) the flow has been time averaged for a short time interval
$75\delta _{in}/u_\infty$
to filter out incoherent fluctuations.

Figure 5. Instantaneous contours of the streamwise velocity fluctuations normalised by the free-stream velocity
$u'/u_{\infty }$
on a wall-parallel plane at
$y/k = 0.97$
(close to the roughness crest). Panel (a) shows the RH_M03 case; panel (b) shows the RH_M2 case.
Additional insights can be gained from the streamwise velocity fluctuations in a wall-parallel plane close to the roughness crest, at
$y/k=0.97$
, figure 5. In both flow cases, the smooth-wall region features the presence of small-scale velocity streaks typical of near-wall turbulence. After the onset of roughness, both the subsonic and supersonic flow cases exhibit a breakdown of the near-wall streaks, although a roughness-induced coherence is apparent as high-speed velocity streaks are preferentially aligned in the gaps between cubes. An important difference between the subsonic and supersonic cases is that the breakdown of the near-wall cycle occurs earlier in the latter, as evident from the small separated region preceding the first row of cubes. This is attributed to the adverse pressure gradient imposed by the shock wave at this location.
4. Results
4.1. Added drag and boundary-layer development
We report a first assessment of the mean flow statistics, analysing the streamwise evolution of the boundary layer in figure 6. Figure 6(a) shows the mean streamwise evolution of the friction coefficient
$C_f=\tau _w/(1/2\rho _{\infty } u_{\infty }^2)$
, where the smooth-to-rough transition is clearly visible for both RH_M03 and RH_M2 cases downstream of
$x/ \delta _{in}\approx 55$
, resulting in an increase of the local drag induced by the roughness. After an initial overshoot, the friction coefficient
$C_f$
decreases following a similar trend in subsonic and supersonic conditions, although with different intensities, and in the last part of the domain it drops as a result of the rough-to-smooth transition right before the outflow.
Figure 6(b) shows the streamwise growth of the boundary-layer thickness
$\delta _{99}$
. The supersonic and subsonic cases exhibit similar growth rates, but the relative intensity is affected by the surface change. The supersonic case RH_M2 exhibits a sharp upward bump in the boundary-layer thickness after the onset of roughness, from which the boundary layer continues to grow. This effect is milder in subsonic case RH_M03, which in turn is only marginally affected by the surface change.
The increased boundary-layer growth for supersonic case RH_M2 is visible when compared with the smooth-wall counterpart SM_M2, and amounts to approximately 1.4 times the reference boundary-layer thickness
$\delta _{99,ref}$
computed before the transition at
$x/ \delta _{in}=45$
. This effect was already visible in the instantaneous flow visualisations in § 3, and we argue that it is directly related to the formation of an initial shock wave, representing a distinct feature of the smooth-to-rough supersonic transition, playing an important role in the development of the internal boundary layer (discussed in § 4.2).
The differences in the streamwise developments of typical lengths of the boundary layer are clearly noted in the profiles of the friction Reynolds number
$Re_{\tau }$
, figure 6(c). On the initial smooth wall, subsonic and supersonic rough-wall cases follow similar trends in
$Re_\tau$
. However, on its rough region, the supersonic case RH_M2 attains slightly higher values of
$Re_\tau$
than the subsonic counterpart RH_M03, as a result of the sharp growth of
$\delta _{99}$
. We highlight that, in order to compare rough and smooth walls at approximately matched
$Re_{\tau }$
values, the domain of the latter is twice as long, see figure 6(e).
To complete the picture, profiles of the wall density are shown in figure 6(d). Subsonic case RH_M03 shows a
$2\,\%$
variation compared with its free-stream value, while supersonic case RH_M2 shows variation of more than
$40\,\%$
, which is expected given the higher Mach number. In addition, we note that the wall density is slightly influenced by the smooth-to-rough transition, increasing by approximately
$3\,\%$
at the onset of the roughness, and then slowly decreasing as the post-shock effects attenuate.

Figure 6. Mean streamwise profiles of (a) skin friction coefficient
$C_f=\tau _w/(1/2\rho _{\infty } u_{\infty }^2)$
, (b) boundary-layer thickness
$\delta _{99}$
, (c) friction Reynolds number
$Re_{\tau }$
and (d) wall density
$\rho _w$
as a function of the streamwise coordinate
$x/ \delta _{in}$
. Panel (e) shows the friction Reynolds number
$Re_{\tau }$
over the whole domain of case SM_M2.
A topic of practical interest is the validity of classical compressibility transformations for predicting the friction coefficient of high-speed flows over rough walls. One of the most employed transformations is van Driest II (van Driest Reference van Driest1956), which allows one to map the compressible friction coefficient into the incompressible one, and it is known to be accurate for adiabatic walls (Pirozzoli & Bernardini Reference Pirozzoli and Bernardini2011),

where
$Re_{\theta }=\rho _{\infty } u_{\infty } \theta / \mu _w$
is the Reynolds number based on the momentum thickness
$\theta$
and
$F_c$
takes the expression

The transformed distribution is compared with the friction formula
$C_{f,i} = 0.024 Re_{\theta ,i}^{-1/4}$
by Smits, Matheson & Joubert (Reference Smits, Matheson and Joubert1983). Figure 7 shows the performance of the aforementioned scaling for all cases in the present database. Here, only specific portions of the domain are considered in order to be far enough from the inlet/outlet and the smooth-to-rough transition, which for the smooth cases results in the range
$20\lt x/\delta _{in}\lt 50$
and for the rough cases
$80\lt x/\delta _{in}\lt 130$
. In general, a very good collapse is observed when using the van Driest transformation for the smooth supersonic case SM_M2 and both in the smooth and rough portions of cases RH_M2 and RH_M03.

Figure 7. Transformed skin friction coefficient
$C_{f,i}$
(4.1) as a function of the Reynolds number based on the incompressible momentum thickness
$Re_{\theta ,i}$
. The dashed grey line represents the friction formula
$C_{f,i}=0.024Re_{\theta ,i}^{1/4}$
.
4.2. Internal boundary layer
This section aims to discuss the growth of the internal boundary layer (IBL), forming at the smooth-to-rough transition because the boundary layer progressively adjusts to the new surface condition. We focus on comparing the theory developed for incompressible boundary layers (Rouhi et al. Reference Rouhi, Chung and Hutchins2019) with supersonic flows and proposing improvements for predicting the growth of the IBL thickness
$\delta _{i}$
when compressibility effects are present.
Several definitions of IBL thickness are available for incompressible flows. One of the most popular is the one by Cheng & Castro (Reference Cheng and Castro2002) who determines
$\delta _{i}$
as the point where the velocity
$\tilde {u}/u_{edge}$
downstream of the transition point is
$99 \,\%$
of the upstream velocity. To compare velocity profiles at different streamwise locations, the wall distance is scaled by the local boundary-layer thickness
$\delta _{99}$
at each station.
In our database, we find that this method has some degrees of arbitrariness, which makes it very sensitive, especially in supersonic cases. On the other hand, we look for a definition for the IBL height that should not directly rely on the local boundary thickness value (both
$\delta _{99}$
or integral values) to scale profiles at different streamwise locations, and not be based on the difference between different streamwise locations (e.g. smooth reference).

Figure 8. Visualisation of the procedure for determining the IBL thickness
$\delta _i$
using the method by Elliott (Reference Elliott1958) for flow cases RH_M03 (a) and RH_M2 (b). The IBL thickness is the intersection (green cross) between two logarithmic regions indicated by dotted lines, intercepting the velocity profiles at the blue (upper logarithmic region) and red (lower logarithmic region) dots.
Among the different definitions developed for incompressible flows, the one from Elliott (Reference Elliott1958) is in line with requirements, and has been regarded by Rouhi et al. (Reference Rouhi, Chung and Hutchins2019) as one of the most consistent with turbulence statistics, while being not dependent on specific thresholds. Elliott (Reference Elliott1958) argues that, if we consider the velocity profiles in inner units right after the surface transition, they will show two logarithmic layers. The upper one would be reminiscent of the upstream surface condition while the lower one is adjusted with the new one. The intersection between these two logarithmic regions is identified as the location of
$\delta _{i}$
. Looking at figure 8, which shows inner-scaled velocity profiles at approximately
$x/\delta _{in}=60$
and the corresponding logarithmic fits, we qualitatively observe good evidence of this theory in both subsonic and supersonic regimes. This observation is confirmed by looking at figure 9, which shows the contours of the shear of the van Driest transformed velocity
$\partial u^+_{VD} / \partial \ln {y^+}$
(see § 4.3 for velocity transformations). Here, we used the van Driest (Reference van Driest1951) transformation to scale velocity profiles in order to have comparable ranges in the contours between the subsonic and the supersonic cases.We also remark that for this analysis we do not adjust for the virtual origin, given that at the smooth-to-rough transition the flow is not in equilibrium with the rough surface and defining a virtual origin would introduce some level of arbitrariness. A drawback of the definition by Elliott (Reference Elliott1958) is that the IBL growth (black dashed line in figure 9) is detectable only for a small streamwise extension after the roughness onset, especially for the supersonic case (
$x/\delta _{in}\lt 70$
). This is because the sharp interface demarcating the IBL only exists for a limited streamwise extent, before blending with the rest of the boundary layer.

Figure 9. Contour plots of
$\partial u^+_{VD} / \partial \ln {y^+}$
for flow cases RH_M03 (a) and RH_M2 (b). The dotted black lines indicate the IBL height
$\delta _i$
predicted using the method by Elliott (Reference Elliott1958).

Figure 10. Contours of
$(\bar {\rho } \widetilde {v^{\prime \prime 2}})/ (\rho _{\infty } u_{\infty })^2$
in a wall-normal plane for cases RH_M03 (a) and RH_M2 (b). The dotted black lines indicate the predicted IBL height
$\delta _i$
using the present method. The prediction using the algorithm from Elliott (Reference Elliott1958) is reported with dashed grey lines for reference.
Most definitions of IBL thickness proposed for incompressible flows (Rouhi et al. Reference Rouhi, Chung and Hutchins2019) that we tested were revealed to be very sensitive, inaccurate or inconsistent at supersonic Mach numbers. We remark that it is desirable to avoid algorithms based on specific thresholds, or dependent of reference wall-normal profiles upstream of the surface transition, as they may be highly sensitive to variations in the flow properties, such as the Mach number. On the other hand, we look for flow variables which can be directly related to the influence of roughness, without being perturbed by other factors (e.g. compressibility effects). In particular, we find the effects of acoustic disturbances and the shock wave emanated at the surface transition to be important, which substantially hampers the reliability of canonical definitions used at lower speeds. For these reasons, consider the Favre-averaged wall-normal momentum equation, which in the boundary-layer approximation reads (Pope Reference Pope2000)

In its integral form, (4.3) indicates that the free-stream pressure is equal to the sum of the static pressure and the wall-normal turbulent stress. From a preliminary analysis of the flow structure of the supersonic case, we argue that, since the local distortion emanated from each element tends to be almost normal to the flow right above the roughness crest (see figure 4
b), the wall-normal component of velocity fluctuations should be the least affected by compression waves in this region, while still being tightly related to the vertical displacement of turbulence caused by the roughness. Additionally, the wall-normal velocity is less affected by compressibility effects, since it is much smaller than the speed of sound, even at supersonic speed. Figure 10 shows the contours of wall-normal fluctuations
$(\bar {\rho } \widetilde {v^{\prime \prime 2}})$
in a longitudinal plane for both RH_M03 and RH_M2 cases. Comparing figures 10(a) and 10(b), the behaviour of wall-normal fluctuations seems to consistently indicate the presence of a region of high intensity (in yellow) starting from the smooth-to-rough transition and developing downstream, without spurious effects due to local pressure gradients or shock waves. In order to connect this observation to the detection of the IBL, we focus on the wall-normal profiles of
$\bar {\rho } \widetilde {v^{\prime \prime 2}}$
, reported in figure 11 at a representative streamwise location
$x=60 \delta _{in}$
.
Here, we observe an inflection point clearly separating two (approximately) linear regions, whose intersection is considered as a definition of the IBL height
$\delta _{i}$
. Repeating this analysis for different streamwise locations gives a measure of the IBL growth, reported as black dotted lines in figures 10(a) and 10(b). Note that the colour contours in figure 10 only serve as a qualitative assessment of the IBL presence, and no specific value should be matched by the dotted lines. To be more precise, we determine the slope and intersection point of these lines by selecting the relative height and location of minima and maxima of the function
$\partial (\bar {\rho } \widetilde {v^{\prime \prime 2}}) /\partial y$
. We argue that this procedure can be considered similar to finding the point of maximum concavity, although it is more reliable from a numerical standpoint. Figure 10 reveals a very good match between the present detection method (black dotted lines) and the one from Elliott (Reference Elliott1958) (grey dotted lines), along with a qualitative agreement with the underlying contours. We note that in the subsonic case (figure 10
a) the two methods nearly coincide both in terms of location and streamwise extent. On the contrary, in the supersonic case (figure 10
b) the prediction by Elliott (Reference Elliott1958) breaks down relatively early, whereas the proposed method is more robust.

Figure 11. Visualisation of the procedure for determining the IBL thickness
$\delta _i$
using the present method for flow cases RH_M03 (a) and RH_M2 (b). The IBL thickness is the intersection (red bullet) between two linear regions (blue and red dotted lines), representative of the slope change in the wall-normal Reynolds stress component
$\bar {\rho } \widetilde {v^{\prime \prime 2}}$
.
Figure 12 shows the predicted height of the IBL scaled by a reference boundary-layer thickness at
$x/\delta _{in}=45$
as a function of the streamwise distance from the surface transition using the algorithm of Elliott (Reference Elliott1958), in figure12(a), and the present method, figure 12(b). We observe that the sizes of the IBL evaluated with the two methods are comparable, but the proposed definition leads to smoother results compared with the one by Elliott (Reference Elliott1958), particularly around the smooth-to-rough transition region. Moreover, both methods exhibit two distinct regions with different power-law exponents
$\delta _i \propto x^{\alpha }$
, approximately before and after
$2.5\delta _{99,ref}$
, the first clearly having a lower exponent
$\alpha$
. We argue that the first region might be affected by an adverse pressure gradient imposed by the roughness right at the transition point, as well as a shock-wave/boundary-layer interaction for the supersonic case, which influences the boundary-layer statistics, resulting in a lower power-law exponent, see figure 12. Focusing on the second region, we find that the aforementioned methods are consistent with each other and estimate a similar exponent
$\alpha$
for supersonic, RH_M2, and subsonic, RH_M03, cases. In particular, the method of Elliott (Reference Elliott1958) yields
$\alpha =0.62$
for subsonic case RH_M03 and
$\alpha =0.58$
for supersonic case RH_M2. Similarly, the present method yields
$\alpha =0.58$
for subsonic case RH_M03 and
$\alpha =0.57$
for supersonic case RH_M2.

Figure 12. Predicted IBL height
$\delta _i/\delta _{99,ref}$
as a function of the streamwise distance from the surface transition, calculated using the method by Elliott (Reference Elliott1958) (a), and the present approach (b). Symbols represent DNS data at
$M_\infty =0.3$
(black) and
$M_\infty =2$
(blue). Dashed and dotted lines represent power-law extrapolations for before and after
$2.5\delta _{99,ref}$
, respectively. Both axes are normalised with a reference boundary-layer thickness
$\delta _{99,ref}$
at
$x/\delta _{in}=45$
.

Figure 13. Growth of the boundary layer,
$\delta _{99}/\delta _{99,ref}$
, (solid lines) and IBL,
$\delta _i/\delta _{99,ref}$
, (dotted lines) as a function of the streamwise distance from the surface transition: (a) using the method by Elliott (Reference Elliott1958) and (b) using the present method. Dotted lines are power-law extrapolations disregarding the region before
$2.5\delta _{99,ref}$
, as shown in figure 12.
It is important to relate these results to the growth of the external boundary layer, which is significant, especially for case RH_M2 (as discussed in § 4.1). Figure 13 shows the IBL growth estimated with the power-law fit and extended throughout the domain. This is compared with the growth of the external boundary layer. Here, the IBL growth of the subsonic and supersonic cases is very similarly when scaled with their respective reference boundary-layer thickness
$\delta _{99,ref}$
, especially the one predicted by the present method, figure 13(b). Figure 13 also shows the streamwise development of the external boundary layer, which is substantially thicker in the supersonic case, as a result of the sharp boundary-layer growth at the roughness onset. On the contrary, figure 14 shows the relative growth of IBL compared with the local boundary-layer thickness
$\delta _{99}$
. The figure reveals that the supersonic flow case exhibits a slower growth towards the equilibrium state. This is attributed to the sharp growth of the external boundary layer right after the roughness onset, as the absolute growth of the IBL seems independent from the Mach number. Looking at figure 14, we can see that, while the IBL for the subsonic case is at almost
$0.9 \delta _{99}$
of the boundary-layer height, the supersonic case still lies below the
$0.8 \delta _{99}$
. These estimates provide a measure of the degree of equilibrium reached by the boundary layer in its adjustment to the different surface condition. We consequently expect the supersonic case RH_M2 to retain mild effects of the smooth-to-rough transition in streamwise stations located downstream, as discussed in § 4.3.
We conclude that the formation and development of the IBL are similar across subsonic and supersonic cases, but when compared with the growth of the external boundary layer we observe a slower recovery for the supersonic case RH_M2 due to the sharp growth of the external boundary layer at the roughness onset. We believe that this is a compressibility effect induced by the shock wave located upstream of the roughness onset.

Figure 14. Growth of the IBL with respect to the local boundary-layer thickness,
$\delta _i / \delta _{99}$
, as a function of the streamwise distance from the surface transition: (a) using the method by Elliott (Reference Elliott1958) and (b) using the present method. Here,
$\delta _i$
is estimated with power-law extrapolations, as in figure 12.
4.3. Velocity statistics
In this section, we consider average velocity profiles at stations listed in table 2. In particular, rough-wall cases are all located at the station
$x/ \delta _{in}=140$
in order to minimise the out-of-equilibrium effects induced by the surface transition and the development of the IBL (§ 4.2). We use compressibility transformations to incorporate the effect of Mach number, comparing the roughness-induced velocity deficit with the incompressible counterpart. It is worth mentioning that rough-wall velocity profiles have been shifted by an effective virtual origin of the flow (Chung et al. Reference Chung, Hutchins, Schultz and Flack2021), which has been chosen as
$d=0.9 k$
. This parameter was found to reduce the uncertainty of the velocity shift and it is in agreement with the value used by Modesti et al. (Reference Modesti, Sathyanarayana, Salvadore and Bernardini2022) for the same roughness geometry. Among the plethora of compressibility transformations, we first report the classical velocity scaling proposed by van Driest (Reference van Driest1951), subscript ‘VD’, together with the two most recent and successful transformations proposed, namely the ones by Griffin, Fu & Moin (Reference Griffin, Fu and Moin2021), subscript ‘GFM’, and Hasan et al. (Reference Hasan, Larsson, Pirozzoli and Pecnik2023), subscript ‘HLPP’, which show improved accuracy in the case of strong density variations. All the above transformations (except the one by Griffin et al. (Reference Griffin, Fu and Moin2021)) can be expressed in terms of convolution integrals mapping the compressible velocity and wall distance into the incompressible ones (denoted with the subscript I),

The kernels functions
$f_I$
and
$g_I$
are reported in table 3. The subscript
$I$
indicates a general compressibility transformation, and it is replaced by the specific subscripts reported in table 3. For details on the approach of Griffin et al. (Reference Griffin, Fu and Moin2021) the reader can refer to the original paper.
Table 3. Compressibility transformations for wall distance and mean velocity according to (4.4), where
$R=\bar {\rho }/\bar {\rho }_w$
and
$M=\bar {\mu }/\bar {\mu }_w$
. In the transformation by Hasan et al. (Reference Hasan, Larsson, Pirozzoli and Pecnik2023),
$D^i=[1-\exp (-y_{TL} / A^{+})]^2$
and
$D^c=[1-\exp (-y_{TL} /(A^{+}+f(M_{\tau })))]^2$
are damping functions,
$A^+$
and
$\kappa$
are constants and
$M_{\tau }=u_{\tau }/\sqrt {\gamma R T_w}$
is the friction Mach number.


Figure 15. Mean velocity profiles for smooth and rough-wall cases obtained at stations listed in table 2, for different compressibility transformations: (a) untransformed, (b) van Driest (Reference van Driest1951, VD), (c) Griffin et al. (Reference Griffin, Fu and Moin2021, GFM), (d) Hasan et al. (Reference Hasan, Larsson, Pirozzoli and Pecnik2023, HLPP). The smooth-wall incompressible case of Sillero, Jiménez & Moser (Reference Sillero, Jiménez and Moser2013) at
$Re_{\tau }=1571$
is used as a reference.
Figure 15 shows the transformed mean velocity profiles for smooth- and rough-wall cases. For the smooth wall, the untransformed mean velocity profile of supersonic case SM_M2 in figure 15(a) shows some compressibility effects, which are accounted for by the compressibility transformations. We find only minor differences between the various transformations, which is expected on an adiabatic wall at this Mach number. Rough-wall cases show the typical downward shift compared with the smooth wall, indicative of higher drag.

Figure 16. Mean velocity deficit between smooth and rough cases
$\Delta \tilde {u}^+$
.
In order to assess the performance of each transformation in collapsing both smooth and rough velocity profiles, we consider the velocity deficit function

with
$\tilde {u}^+_S (y^+)$
and
$\tilde {u}^+_R (y^+)$
the mean velocity of the smooth- and rough-wall cases, respectively. If outer-layer similarity holds,
$\Delta u^+$
is nearly constant in the log layer (i.e.
$100\lt y^+\lt 0.3 Re_\tau$
), although this is only a necessary condition as second-order quantities may still be out of equilibrium (Li et al. Reference Li, De Silva, Rouhi, Baidya, Chung, Marusic and Hutchins2019). Figure 16 shows
$\Delta u^+$
of the untransformed and transformed velocity profiles. In general, we observe a nearly constant trend across the different compressibility transformations, but the most recent transformations of Griffin et al. (Reference Griffin, Fu and Moin2021) and Hasan et al. (Reference Hasan, Larsson, Pirozzoli and Pecnik2023) show a superior agreement with the subsonic reference. The velocity deficit is a convenient measure of the added drag, as unlike the relative friction coefficient, it does not depend on the Reynolds number (Chung et al. Reference Chung, Hutchins, Schultz and Flack2021). Throughout this work, we evaluate the velocity shifts at the nominal edge of the logarithmic region
$\Delta U^+ = \Delta u^+ (y^+\approx 0.3 Re_{\tau } )$
. Measuring at a different location closer to the wall does not affect the results as
$\tilde {u}^+$
is approximately constant in this region. A important difference with respect to the incompressible flow regime where
$\Delta U^+$
directly relates to the added drag is that, for compressible flows, the density ratio
$R=\rho _{\infty }/\rho _w$
also plays a role and the ratio between smooth and rough friction coefficients can be written as
$C_{fs}/C_{f}=R/R_s(1-\Delta U^+/\tilde {u}^+_s)^2$
(Modesti et al. Reference Modesti, Sathyanarayana, Salvadore and Bernardini2022).
As customary for rough walls, introducing the equivalent sand-grain roughness
$k_s$
helps the standardisation of the present roughness geometry, therefore, we evaluate its relation to the geometric roughness height
$k$
. The relation between
$k$
and
$k_s$
is at the basis of the Moody diagram, where the drag of a given rough geometry can be estimated by relating it to the drag produced by an equivalent sand-grain level. As discussed in Chung et al. (Reference Chung, Hutchins, Schultz and Flack2021), if the roughness Reynolds number is based on the equivalent sand-grain roughness height,
$k_s^+=k_s/\delta _\nu$
, and the separation of scales is large enough to yield the fully rough regime, the velocity deficit function
$\Delta U^+$
takes the form

Here,
$\kappa \approx 0.41$
is the von Kármán constant,
$A=5.2$
and
$B_s=8.5$
.
We first find the ratio between
$k_s$
and
$k$
by matching the incompressible case RH_M03 with the fully rough asymptote, (4.6), which yields
$k_s/k=1.84$
. This is in agreement with Modesti et al. (Reference Modesti, Sathyanarayana, Salvadore and Bernardini2022), who reported a value of 1.9 for the same roughness geometry in a turbulent channel flow. For the supersonic flow case we follow the approach of Modesti et al. (Reference Modesti, Sathyanarayana, Salvadore and Bernardini2022), namely we introduce the incompressible roughness height
$k_{I} = y_{I}(k)$
and then calculate
$k_{s,I}$
using the same ratio of the subsonic flow case
$k_{s,I}/k_I=1.84$
. Figure 17 shows
$\Delta U^+$
as a function of
$k_s^+$
for the RH_M2 case using both the untransformed and transformed values according to Hasan et al. (Reference Hasan, Larsson, Pirozzoli and Pecnik2023).
We observe that using a suitable compressibility transformation yields a slightly better agreement with the incompressible asymptote, although at the present Mach number the effect of the compressibility transformation is minor. We believe that this result is powerful as it supports a framework in which compressibility transformations can act as a tool to achieve Mach-number independency, thus maintaining the applicability of the incompressible theory of roughness. Knowing the Hama roughness function for an incompressible flow, one can, in principle, calculate
$\Delta U^+$
in the compressible regime by applying the inverse compressibility transformation (4.4).

Figure 17. Shift of the mean streamwise velocity
$\Delta \tilde {u}^+$
as a function of the inner-scaled equivalent sand-grain roughness height
$k_s^+$
. The case RH_M2 is also reported using velocity profiles transformed with the relation of Hasan et al. (Reference Hasan, Larsson, Pirozzoli and Pecnik2023). The dotted grey line indicates the Hama roughness function
$\Delta U^+ = 1/k \ln {k_s^+}+8.5-5.2$
, while open circles are data from Nikuradse (Reference Nikuradse1933).
Figure 18 reports the Reynolds stress components as a function of
$y^+$
, figure 18(a), and
$y/\delta _{99}$
, figure 18(b). Here, the wall-normal coordinates for rough-wall cases are again shifted by the virtual origin in order to assess the outer-layer similarity. In general, we note that both cases follow quite well the behaviour of the smooth-wall counterpart in the outer layer, with few discrepancies in the streamwise component
$\tau _{11}$
. Among the two rough cases, the subsonic case RH_M03 seems to better agree with the reference case of Sillero et al. (Reference Sillero, Jiménez and Moser2013) up to
$y^+\gt 300$
, while both profiles exhibit a plateau approaching the roughness crest, located at approximately
$y^+\approx k^+ \approx 60$
. For the same component, supersonic case RH_M2 shows a more intense level of fluctuations when compared with the corresponding smooth case SM_M2 at a similar
$Re_{\tau }$
. For this case, the influence of the smooth-to-rough transition extends further downstream, as discussed in § 4.2, and this is especially true for the second-order statistics. We attribute this effect to the modulation of the turbulent structures on the outer layer induced by the initial shock wave, that is expected to mainly affect the streamwise velocity fluctuations.
The Reynolds shear stress
$\tau _{12}$
of the rough-wall cases shows an excellent collapse even with the respective smooth-wall reference case. A clear plateau close to a value of
$\tau _{12}/\tau _w \approx 1$
in the log layer is established, meaning that a constant stress layer exists in all cases.

Figure 18. Turbulent velocity fluctuations
$\tau _{ij}=\widetilde {u_i^{\prime }u_j^{\prime }}$
scaled with the wall shear stress
$\tau _w$
as a function of the wall-normal distance in wall units (a)
$y^+, y^+-d^+$
and outer units (b)
$y/\delta _{99},(y-d)/\delta _{99}$
. Rough-wall cases are adjusted using a virtual origin shift
$d=0.9k$
.
4.4. Thermal statistics
We consider thermodynamic statistics at stations listed in table 2 for smooth and rough supersonic cases, SM_M2 and RH_M2, respectively. This aspect is particularly interesting in order to assess if the outer-layer similarity, extensively studied for velocity statistics, can be also observed in the temperature field. Figure 19 shows the mean,
$\tilde {T}$
, and root-mean-square,
$\tilde {T}_{rms}$
, temperature profiles as a function of
$y^+$
. Profiles are not shifted with a virtual origin in order to show their behaviour even below the roughness crest. Both
$\tilde {T}$
and
$\tilde {T}_{rms}$
of flow cases RH_M2 and SM_M2 are fundamentally different. In particular, the mean temperature profile
$\tilde {T}$
, figure 19(a), is consistently higher for the rough case, excluding the region very close to the wall where it settles to slightly lower values. Looking at the temperature fluctuation root mean square,
$\tilde {T}_{rms}$
, figure 19(b), we note even more significant differences between the smooth and rough cases. In particular, temperature fluctuations are mostly damped in the roughness layer, while an intense peak appears in the outer layer, right below the edge of the boundary layer. These findings are consistent with the results of Modesti et al. (Reference Modesti, Sathyanarayana, Salvadore and Bernardini2022) for turbulent channel flows over rough walls with an isothermal wall condition. They concluded that there is no evidence of outer-layer similarity for the temperature field, and roughness is able to influence the temperature fluctuations throughout, up to the edge of the boundary layer. Similar conclusions are evident for the present cases despite the fact that adiabaticity is enforced at the wall, which for smooth cases is a condition that yields a strong similarity between velocity and temperature fields even at higher Mach numbers (Cogo et al. Reference Cogo, Baù, Chinappi, Bernardini and Picano2023). This aspect suggests that the interaction between compressibility and roughness fundamentally alters the thermal field and this process dominates over the influence of the wall temperature condition.
From a physical point of view, we can explain the influence of roughness in altering the temperature fluctuation profile by considering the definition of the thermal production term
$\mathcal {P}_T= - \bar {\rho } \widetilde {v^{\prime \prime } T^{\prime \prime }} \partial \tilde {T} / \partial y$
, see for example Cogo et al. (Reference Cogo, Baù, Chinappi, Bernardini and Picano2023). In particular, of note is the behaviour of the wall-normal component of the Reynolds stresses
$\bar {\rho } \widetilde {v^{\prime \prime 2}}$
, figure 18, which influences the velocity–temperature fluctuation correlation
$\bar {\rho } \widetilde {v^{\prime \prime } T^{\prime \prime }}$
, and the mean temperature profile
$\tilde {T}$
, figure 19(a), whose gradient appears in the expression of
$\mathcal {P}_T$
. The temperature profile appears to stagnate around a nearly constant value, roughly
$\tilde {T}_w$
, up to the roughness crest, while the smooth profile is already at much lower and decreasing temperatures. Hence, the vanishing mean temperature gradient in the roughness sublayer dampens the near-wall thermal fluctuation production
$\mathcal {P}_T$
. Far from the wall, the rough-wall profile is forced to have a higher wall-normal temperature gradient
$\partial \tilde {T} / \partial y$
in order to reach the edge temperature, while wall-normal velocity fluctuations
$\bar {\rho } \widetilde {v^{\prime \prime 2}}$
are similar to those of the smooth-wall profile in the same region (see figure 18). Hence, we attribute the outer-layer peak of temperature fluctuations
$\tilde {T}_{rms}$
, figure 19(b), to the steeper mean temperature gradient
$\partial \tilde {T} / \partial y$
, which increases the thermal fluctuation production
$\mathcal {P}_T$
.

Figure 19. (a) Normalised temperature profiles
$\tilde {T}/T_{\infty }$
for cases RH_M2 and SM_M2 as a function of the wall-normal distance
$y^+$
. (b) Temperature fluctuations scaled with the wall temperature
$\tilde {T}_{rms}/T_{w}$
as a function of
$y^+$
. Vertical dashed grey lines represent the location of the roughness crest
$y^+=k^+$
.
We finally assess the coupling between the average thermal and kinetic fields through the well-known temperature–velocity quadratic relation (Busemann Reference Busemann1931; Crocco Reference Crocco1932), which has been improved over the years to account for finite heat transfers and high Mach numbers (Zhang et al. Reference Zhang, Bi, Hussain and She2014),

where
$T_{rg}=T_{\infty }+r_g U_{\infty }^2/(2 c_p)$
and
$r_g=2 c_p (T_w-T_{\infty })/U_{\infty }^2-2 \, Pr \, q_w/(U_{\infty } \tau _w)$
. Figure 20 compares smooth- and rough-wall profiles of the mean temperature as a function of the mean velocity with their respective estimates given by the Zhang et al. (Reference Zhang, Bi, Hussain and She2014) relation. As expected, the smooth DNS profile is well approximated by the temperature–velocity relation, (4.7). The supersonic rough-wall case RH_M2 closely follows its smooth counterpart in the region denoted by
$\tilde {u}/u_{\infty }\gt 0.5$
, showing a slight mismatch below this threshold. Despite this small difference, we conclude that the temperature–velocity relation of Zhang et al. (Reference Zhang, Bi, Hussain and She2014) appears to be a robust tool even in rough-wall supersonic boundary layers. As pointed out by Modesti et al. (Reference Modesti, Sathyanarayana, Salvadore and Bernardini2022), the fact that a quadratic relation exists between temperature and velocity invalidates the outer-layer similarity for the mean temperature by construction.

Figure 20. Temperature–velocity relation for cases RH_M2 (solid) and SM_M2 (dashed). Direct numerical simulation data (blue) are compared with (4.7), in grey.
5. Conclusions
We have studied compressible turbulent boundary layers over smooth and rough surfaces using DNS. In particular, we compared the flow dynamics of subsonic (
$M_{\infty }$
= 0.3) and supersonic (
$M_{\infty }$
= 2) turbulent boundary layers with a smooth-to-rough surface transition, the rough part being composed of 3-D cubical elements equally spaced in the streamwise and spanwise directions.
First, we analysed the instantaneous flow features of each case, noting clear differences in the supersonic case due to compressibility effects. In particular, the supersonic case RH_M2 features an oblique shock wave at the onset of roughness, followed by a pattern of subsequent local waves emanating from each element and extending into the free stream, consistent with the observations of previous experimental studies (Latin & Bowersox Reference Latin and Bowersox2000; Ekoto et al. Reference Ekoto, Bowersox, Beutner and Goss2008; Kocher et al. Reference Kocher, Kreth, Schmisseur, LaLonde and Combs2022).
For the first time, we extend the study of smooth-to-rough transition and the formation of the IBL to compressible flows. In order to detect the edge of the IBL, we compare the classical method by Elliott (Reference Elliott1958) with a newly proposed approach based on the wall-normal velocity fluctuations. The two methods are similar, but the latter is revealed to be more robust for supersonic flows, leading to a sharper determination of the IBL.
Interestingly, we find that the absolute growth of the IBL is not affected by the Mach number, probably because for our flow case the Mach number in the IBL does not exceed 1.5. However, relative to the local boundary-layer thickness, the supersonic case shows a slower recovery, due to the increased thickness of the external boundary layer over the roughness.
Regarding the velocity statistics, we observe outer-layer similarity for the mean velocity and the Reynolds stresses in both the subsonic and supersonic regime, although minor deviations in the streamwise Reynolds stress component are visible for the supersonic case, in agreement with the slower recovery rate suggested by the IBL analysis.
Regarding the temperature statistics, we find that outer-layer similarity does not hold, neither for the mean temperature nor for the temperature fluctuations, corroborating the analysis of Modesti et al. (Reference Modesti, Sathyanarayana, Salvadore and Bernardini2022). We point out that this cannot be attributed solely to the thermal boundary condition, as the present dataset features adiabatic walls, as opposed to the cold walls of Modesti et al. (Reference Modesti, Sathyanarayana, Salvadore and Bernardini2022). Nonetheless, we note that the lack of universality in the mean temperature does not complicate its modelling. The present data indicate that the classical temperature–velocity relation also holds over rough surfaces, allowing the temperature field to be predicted from the mean velocity.
Regarding the prediction of added drag, we confirm the finding of Modesti et al. (Reference Modesti, Sathyanarayana, Salvadore and Bernardini2022), namely, that it is possible to use compressibility transformations to map the compressible velocity shift and equivalent sand-grain roughness height onto the incompressible counterparts. This has the important implication that added drag in supersonic flows can, in principle, be predicted from knowledge of
$\Delta U^+$
and
$k_s^+$
in the incompressible regime, for the same roughness type. While this framework is very powerful, several questions on its applicability for supersonic flows remain to be answered. For instance, Mach-number effects have been barely explored, similarly the thermal boundary conditions, and the roughness shape. Therefore, additional studies are required to full assess the applicability of this framework and its accuracy in estimating the added drag.
Another point that deserves further scrutiny is the development of the IBL and its dependency on the Mach number, the roughness shape and the flow blockage induced at the transition onset. The present data would suggest that the boundary-layer growth induced by the roughness depends on the intensity of the initial shock and the pressure jump it imposes. This could potentially alter the universality of the flow downstream.
Acknowledgements
We acknowledge that the results reported in this paper have been achieved using the EuroHPC JU Extreme Scale Access Infrastructure resource Leonardo Booster based at CINECA, Casalecchio di Reno, Italy, under project EUHPC-E02-044. We also acknowledge the CINECA award under the ISCRA initiative, for the availability of high-performance computing resources and support.
Funding
This research received financial support from ICSC – Centro Nazionale di Ricerca in ‘High Performance Computing, Big Data and Quantum Computing’, funded by European Union – NextGenerationEU. We acknowledge financial support under the National Recovery and Resilience Plan (NRRP), Mission 4, Component 2, Investment 1.1, Call for tender No. 104 published on 2.2.2022 by the Italian Ministry of University and Research (MUR), funded by the European Union – NextGenerationEU – Project Title ‘ADvanced Modeling of hIgh-speed aerodynamics for MaRs Entry (ADMIRE)’, PRIN project number ‘20224NY2YK’ – CUP C53C24000830006 – Grant Assignment Decree No. 20435 adopted on 06-11-2024 by the Italian Ministry of University and Research (MUR).
Declaration of interests
The authors report no conflict of interest.
Data availability statement
The data that support the findings of this study are available upon reasonable request.