1 Introduction
Superhydrophobicity is a property attributed to surface roughness (ridges, grooves, posts or random textures) and surface chemistry which maintains large contact angles for sessile drops, thus producing low wettability, known as the Cassie–Baxter state (Cassie & Baxter Reference Cassie and Baxter1944). The interface meniscus creates an air mattress that acts like a lubricant for the outer flow (Rothstein Reference Rothstein2010). When the interface fails, the liquid fills the surface cavities and the superhydrophobic effect is lost. This is referred to as the Wenzel state (Wenzel Reference Wenzel1936).
Nature provides numerous examples of superhydrophobic surfaces (SHS), which can be exploited for practical applications. For example, the lotus leaf is believed to take advantage of superhydrophobicity for a self-cleaning mechanism (Barthlott & Neinhuis Reference Barthlott and Neinhuis1997). Frictional drag reduction is central to the performance of marine vessels, and anti-biofouling, anti-icing and microfluidic devices (Furstner et al. Reference Furstner, Barthlott, Neinhuis and Walzel2005; Genzer & Efimenko Reference Genzer and Efimenko2006; Fang et al. Reference Fang, Li, Wang and Qiao2008; Jung et al. Reference Jung, Dorrestijn, Raps, Das, Megaridis and Poulikakos2011). Any impact on skin friction drag reduction substantially improves the overall performance and yields savings in fuel cost (Choi & Kim Reference Choi and Kim2006). In the present work we focus on drag reduction using SHS in two canonical flow configurations: laminar Couette flow and turbulent channel flow.
With recent developments in three-dimensional printing and microfabrication processes, it is possible to create surfaces exhibiting superhydrophobic characteristics when coupled with chemical treatments. Laminar flows over SHS have been studied both numerically and experimentally. SHS have been shown to achieve drag reduction (Ou, Perot & Rothstein Reference Ou, Perot and Rothstein2004; Ou & Rothstein Reference Ou and Rothstein2005; Choi & Kim Reference Choi and Kim2006; Joseph et al. Reference Joseph, Cottin-Bizonne, Benoit, Ybert, Journet, Tabeling and Bocquet2006; Maynes et al. Reference Maynes, Jeffs, Woolford and Webb2007; Woolford et al. Reference Woolford, Prince, Maynes and Webb2009; Emami et al. Reference Emami, Tafreshi, Gad-el Hak and Tepper2011). Analytical models relate the slip lengths to various surface parameters such as groove width, pitch and height (Lauga & Stone Reference Lauga and Stone2003; Ybert et al. Reference Ybert, Barentin, Cottin-Bizonne, Joseph and Bocquet2007) or the slip velocities to geometry (Seo & Mani Reference Seo and Mani2016). In general the SHS are considered to be simple grooved geometries, and numerically the interface is typically assumed to be flat and represented using zero-shear boundary conditions. Others have included the effect of viscosity on the interface (Vinogradova Reference Vinogradova1995; Belyaev & Vinogradova Reference Belyaev and Vinogradova2010; Nizkaya, Asmolov & Vinogradova Reference Nizkaya, Asmolov and Vinogradova2014). Several authors have investigated the effect of the curvature due to the meniscus and modified the analytical solutions to take curvature into account (Cottin-Bizonne et al. Reference Cottin-Bizonne, Barrat, Bocquet and Charlaix2003; Sbragaglia & Prosperetti Reference Sbragaglia and Prosperetti2007b ; Wang, Teo & Khoo Reference Wang, Teo and Khoo2014; Li, Alame & Mahesh Reference Li, Alame and Mahesh2017).
Turbulent flows over textured surfaces have been studied extensively in the past. Experimentally, it becomes difficult to conduct measurements near the wall and to maintain a stable interface, but drag reduction and slip lengths have been investigated. Some past work reported that SHS had no effect on turbulent statistics (Zhao, Du & Shi Reference Zhao, Du and Shi2007; Peguero & Breuer Reference Peguero, Breuer and Eckhardt2009), while others reported otherwise (Gogte et al.
Reference Gogte, Vorobieff, Truesdell, Mammoli, van Swol, Shah and Brinker2005; Henoch et al.
Reference Henoch, Krupenkin, Kolodner, Taylor, Hodes, Lyons, Peguero and Breuer2006; Daniello, Waterhouse & Rothstein Reference Daniello, Waterhouse and Rothstein2009; Jung & Bhushnan Reference Jung and Bhushnan2009; Aljallis et al.
Reference Aljallis, Sarshar, Datla, Sikka, Jones and Choi2013; Bidkar et al.
Reference Bidkar, Leblanc, Kulkarni, Bahadur, Ceccio and Perlin2014; Park, Sun & Kim Reference Park, Sun and Kim2014; Srinivasan et al.
Reference Srinivasan, Kleingartner, Gilbert, Cohen, Milne and McKinley2015). Investigation of interface stability was studied using post-processed pressure fluctuations (Seo, Garcia-Mayoral & Mani Reference Seo, Garcia-Mayoral and Mani2015). Rosenberg et al. (Reference Rosenberg, Van Buren, Fu and Smits2016) showed that the turbulent skin friction is reduced over air- and liquid-impregnated surfaces (SLIPS) for Taylor–Couette flows. Numerically, the interface is assumed to be flat and modelled using zero-shear boundary conditions (Martell, Perot & Rothstein Reference Martell, Perot and Rothstein2009; Frohnapfel, Hasegawa & Kasagi Reference Frohnapfel, Hasegawa and Kasagi2010; Martell, Rothstein & Perot Reference Martell, Rothstein and Perot2010; Park, Park & Kim Reference Park, Park and Kim2013; Jelly, Jung & Zaki Reference Jelly, Jung and Zaki2014; Türk et al.
Reference Türk, Daschiel, Stroh, Hasegawa and Frohnapfel2014) and homogenised slip length models instead of zero-shear boundary conditions (Min & Kim Reference Min and Kim2004; Fukagata, Kasagi & Koumoutsakos Reference Fukagata, Kasagi and Koumoutsakos2006; Busse & Sandham Reference Busse and Sandham2012). The profiles of mean velocity, turbulence intensities and Reynolds shear stresses were characterised for the inner part of turbulent boundary layers over several SHS with varying textures and a range of
$Re_{\unicode[STIX]{x1D70F}}$
(Ling et al.
Reference Ling, Srinivasan, Golovin, McKinley, Tuteja and Katz2016). Jung, Choi & Kim (Reference Jung, Choi and Kim2016) studied the effect of anisotropy in the slip-length models for different interface heights in idealised SHS. Recently, Fu et al. (Reference Fu, Arenas, Leonardi and Hultmark2017) used the level-set method to study SHS with varying viscosity ratios over spanwise and streamwise grooves. Rastegari & Akhavan (Reference Rastegari and Akhavan2018) applied the Boltzmann method in their DNS to study the effect of the liquid–gas interface in longitudinal grooves by modelling it as a stationary, curved and free shear boundary; the meniscus shape was determined using the Young–Laplace equation. The sustainable pressure bounds of SHS were further investigated by Rastegari & Akhavan (Reference Rastegari and Akhavan2019). Fairhall, Abderrahaman-Elena & García-Mayoral (Reference Fairhall, Abderrahaman-Elena and García-Mayoral2019) showed that drag reduction is proportional to the difference between the virtual origin of the mean flow and the virtual origin of the overlying turbulence.
Most past work on flow over SHS has considered idealised geometries such as grooves or posts. Seo & Mani (Reference Seo and Mani2018) recently studied turbulent flow over SHS idealised as random slip/no-slip patches. To the best of our knowledge, none of the past numerical work has simulated a multiphase flow over realistically rough surfaces as presented in this paper. The main goal of the present work is to perform DNS of (i) laminar Couette flow and (ii) turbulent channel flow, where the bottom wall is a realistically rough surface. We aim to explore the effect of interface height on slip, drag reduction, near-wall flow field and turbulence statistics. The rest of the paper is organised as follows: § 2 describes simulation details including the numerical method, parameters and problem formulation. Results are described in § 3 which include flow visualisations, steady and mean flow field properties and drag reduction in laminar Couette flow. The mean flow statistics, scaling laws and flow structure are also presented for turbulent flow. Finally, the work is summarised in § 4.
2 Simulation details
2.1 Numerical method
Direct numerical simulation (DNS) is performed using a mass-conserving volume-of-fluid (VOF) methodology on structured grids to study the effect of an air–water interface over a realistically rough surface. The governing equations are solved using the finite-volume algorithm developed by Mahesh, Constantinescu & Moin (Reference Mahesh, Constantinescu and Moin2004) for the incompressible Navier–Stokes equations. The governing equations for the momentum and continuity are given by the Navier–Stokes equations:
where
$u_{i}$
and
$x_{i}$
are the
$i$
th component of the velocity and position vectors respectively,
$p$
denotes pressure,
$\unicode[STIX]{x1D70C}$
is density and
$\unicode[STIX]{x1D707}$
is the viscosity of the fluid. The fluids are assumed to be immiscible. Additionally in (2.1),
$\unicode[STIX]{x1D6FF}_{i1}$
is the Kroenecker delta,
$K_{i}$
is the body force which is only active in the liquid phase and
$F_{st,i}$
the surface tension force. The algorithm is robust and emphasises discrete kinetic energy conservation in the inviscid limit which enables it to simulate high-Reynolds-number flows without adding numerical dissipation. The solution is advanced in time by an implicit scheme using successive over-relaxation (SOR). A predictor–corrector methodology is used where the velocities are first predicted using the momentum equation and then corrected using the pressure gradient obtained from the Poisson equation yielded by the continuity equation. The Poisson equation is solved using a multigrid pre-conditioned conjugate gradient method (CGM) using the Trilinos libraries (Sandia National Labs). The multigrid pre-conditioner uses a Chebyshev smoother with a third-order polynomial and a maximum number of levels set to
$4$
. The implicit time advancement uses the Crank–Nicholson discretisation with a linearisation of the convection terms.
The volume fraction is represented by a colour function
$c$
to keep track of two different fluids. The colour function
$c$
varies between the constant value of one in a filled cell to zero in an empty cell, with an intermediate value between zero and one to define an interface cell where
$0\leqslant c\leqslant 1$
. The reconstruction and advection steps are based on a set of analytic relations proposed by Scardovelli & Zaleski (Reference Scardovelli and Zaleski2000). The governing equations for the colour function material derivative are given by

where the advection term is neglected in the following simulations since we assume that the interface is stationary due to an infinite surface tension. The density and viscosity are evaluated as
where the subscript ‘
$l$
’ denotes the liquid phase and ‘
$g$
’ the gas phase. The surface is represented by obstacle cells which are masked out. At the beginning of a simulation run, the fluid and obstacle cells are flagged accordingly:
this step is performed once. The wetted masked cells (cells that share a face between a fluid and obstacle cell) enforce a zero face-normal velocity
$v_{N}|_{mask}=0$
. The cell-centred velocities satisfy a no-slip boundary condition, with the exception of corner cells that take a weighted average of the neighbouring cell-centred values. The algorithm has been validated with experimental results for a variety of flows involving superhydrophobicity (Li, Alame & Mahesh Reference Li, Alame and Mahesh2016; Li et al.
Reference Li, Alame and Mahesh2017) and fully wetted roughness (Ma, Alamé & Mahesh Reference Ma, Alamé and Mahesh2019). In this study, we enforce a zero face-normal velocity at the interface
$v_{N}|_{interface}=0$
. The condition models a high-surface-tension regime with a stable flat interface. This was done to focus on the effect of varying interface heights
$h$
for a finite-viscosity lubricant. The assumption made is valid for flow regimes where the interfacial surface tension dominates the interface dynamics. Further discussion of the validity of our assumptions is presented in appendix A. The statistics of the turbulent channel flow were averaged over a period of 300 flow-through times after the discharge had reached a steady-state value.
2.2 Surface generation
The roughness used in the present work is obtained from a real surface manufactured at UT Dallas (courtesy Professor Wonjae Choi), with a three-dimensional (3-D) surface profile measurement using a 20X objective lens obtained from MIT (courtesy Professor Gareth McKinley). The sample is aluminium 6061 sandblasted using 150 grit, etched for 25 s, boehmetized for 30 min and hydrophobised using Ultra Ever Dry top coat in isopropanol. Figure 1 provides an illustration of the scanned surface data coloured with height. The surface statistics and power spectral density (PSD) of the surface height are provided in appendix B.
Illustration of the real rough surface. The contour legend describes the height of the surface profile.

We begin with a pre-processing step by reading the scanned surface data. The number of pixels in the scan width and height are stored as the number of nodes in the streamwise and spanwise directions respectively. The values of the roughness height and spatial location are then interpolated to cell centres given our domain of choice. The cell centre values are then written to a new file with a structured data format. Any obstacle cell which shares an edge with a fluid cell is tagged as a boundary cell. Boundary cells can either be an edge cell (if the boundary cell borders exactly one fluid cell) or a corner cell (if the boundary cell shares a corner with two or more fluid cells). The discretised surface is checked with the original data and the errors in the surface statistics are presented in table 5, appendix B. The momentum equations are solved inside the fluid domain while the pressure is solved everywhere. The weighted average applied at the corner cells does not affect the pressure equation since we use collocated grids where the face-normal velocities are set to zero at the boundaries independent of the cell centre value. This ensures a proper pressure jump recovery at the obstacle walls where the values inside the obstacle domain do not affect the pressure values in the fluid domain.
2.3 Problem description
Simulations are performed for two canonical problems: (i) laminar Couette flow and (ii) turbulent channel flow, where the surface described in § 2.2 is used as the bottom wall. In the experiments performed by Ling et al. (Reference Ling, Srinivasan, Golovin, McKinley, Tuteja and Katz2016), the tunnel pressure (which controls the interface location) is increased, which compresses the air layer into the SHS and in turn exposes more asperities, thereby reducing the extent of drag reduction. The aim of this paper is to model this effect over an idealised flat interface numerically by progressively increasing the height
$h$
and measuring flow properties for each interface location in different flow regimes. The maximum interface height is non-physical in a realistic scenario, but it serves the purpose of providing the largest amount of slip that is theoretically achievable. It also helps describe the trend between limiting cases. The problem description is given in the following sections.
Illustration of the computational domain for the laminar Couette simulation, roughness and the instantaneous velocity field.

Case names, interface location, domain extents and grid resolution for the laminar Couette flow problem. L denotes the laminar cases. S and R denote a smooth and rough wall respectively. Fully wetted roughness is denoted by FW. I1, I2 and I3 represent the interface height at three locations:
$0$
,
$S_{q}$
and
$S_{p}$
respectively.

2.3.1 Laminar Couette
The height
$H$
of the top wall was chosen such that the root-mean-square (RMS) roughness height
$S_{q}$
is around 2 % of
$H$
. The original surface is scaled to achieve the roughness height ratios described above. The reference system is chosen such that the origin coincides with the arithmetic mean elevation of the roughness. The schematic diagram shown in figure 2 illustrates the flow domain. No-slip boundary conditions are prescribed on the bottom surface and a constant velocity
$U_{\infty }$
in the streamwise
$x$
-direction is prescribed at the top wall. The streamwise (
$x$
) and spanwise (
$z$
) directions are periodic; a non-uniform grid is used in the wall-normal (
$y$
) direction with clustering in the rough wall region. The interface location was varied from the maximum valley depth
$S_{v}$
all the way up to the maximum peak height of the roughness
$S_{p}$
over
$18$
increments. Table 1 gives the grid details. The Reynolds number
$Re=U_{\infty }H/\unicode[STIX]{x1D708}=740$
, where
$U_{\infty }$
and
$H$
are taken to be unity and the liquid phase being the reference material property. A smooth planar Couette flow (Case L-S) is used as a baseline such that the reference shear stress
$\unicode[STIX]{x1D70F}_{o}=\unicode[STIX]{x1D707}U_{\infty }/H=\unicode[STIX]{x1D707}_{w}$
, where
$\unicode[STIX]{x1D707}_{w}$
is the reference viscosity in the water phase. First, a fully wetted case (L-RFW) is simulated to baseline the effect of roughness on drag when compared to the smooth wall. Case L-RI1 denotes
$h=0$
, Case L-RI2 denotes
$h=S_{q}$
and Case L-RI3 denotes
$h=S_{p}$
. The viscosity ratio
$\unicode[STIX]{x1D707}_{r}=\unicode[STIX]{x1D707}_{a}/\unicode[STIX]{x1D707}_{w}=1/50$
is used to represent an air–water interface. The change in shear stress due to the roughness and
$h$
is used to compute the drag reduction defined using the following relation:
Illustration of the computational domain for the turbulent channel flow, interface location embedded within the roughness and the instantaneous velocity field.

Case names, interface location, domain extents and the grid resolution in wall units for the turbulent channel flow problem. T denotes the turbulent cases. S and R denote a smooth and rough wall respectively. Fully wetted roughness is denoted by FW. I1, I2 and I3 represent the interface height at three locations:
$0$
,
$S_{q}$
and
$S_{p}$
respectively.

2.3.2 Turbulent channel
A schematic diagram describing the turbulent channel domain is given in figure 3. No-slip boundary conditions are applied on both the top smooth wall and the bottom rough wall with periodicity in the streamwise (
$x$
) and spanwise (
$z$
) directions; non-uniform grids are used in the wall-normal (
$y$
) direction where the grid is clustered near the rough wall region. The grid details are given in table 2. A constant body force in the liquid phase is applied such that the friction Reynolds number is
$Re_{\unicode[STIX]{x1D70F}}=u_{\unicode[STIX]{x1D70F}}\unicode[STIX]{x1D6FF}/\unicode[STIX]{x1D708}=180$
where
$u_{\unicode[STIX]{x1D70F}}$
is the wall friction velocity,
$\unicode[STIX]{x1D6FF}=(L_{y}-y_{o})/2$
the channel half-height and
$y_{o}$
the reference bottom plane. Four cases were considered: (i) fully wetted rough channel for Case T-RFW, (ii, iii) two-phase rough channel with
$h=0$
for Case T-RI1,
$h=S_{q}$
for Case T-RI2, and (iv)
$h=S_{p}$
for Case T-RI3. The reference plane
$y_{o}$
is taken to be the arithmetic mean elevation of the roughness for Case T-RFW, and the location of the interface for Cases T-RI1, T-RI2 and T-RI3. The viscosity ratio is that of an air–water interface given by
$\unicode[STIX]{x1D707}_{r}=\unicode[STIX]{x1D707}_{a}/\unicode[STIX]{x1D707}_{w}=1/50$
. The solver was validated (not shown here) for the flat smooth channel (Kim, Moin & Moser Reference Kim, Moin and Moser1987). The original surface is scaled such that
$S_{q}^{+}\approx 1.6$
. The original surface scan was not large enough to cover the bottom wall after scaling, and the roughness patch had to be tiled in random orientations to minimise any directional bias and create a larger area. The required domain extents were then extracted from the tiled surface as shown in figure 4(a). The computational domain required after scaling was twice as long in the streamwise (
$x$
) direction and
$35\,\%$
longer in the spanwise (
$z$
) direction. The ratio of
$L_{x}$
to
$L_{z}$
is
$4:3$
in the laminar Couette case compared to the
$2:1$
ratio of the turbulent channel case. In order to ensure this does not affect the roughness height distribution, a probability density function (p.d.f.) distribution is plotted for both cases in figure 4(b) and is compared to a Gaussian distribution. No appreciable difference is observed between the two cases. The p.d.f. distribution is negatively skewed when compared to a Gaussian which is also calculated in the surface statistics presented in table 4 of appendix B.
(a) Turbulent channel domain extent (dashed line) and the tiled surface at different orientations with their physical boundaries (solid red), not drawn to scale. (b) The probability density function (p.d.f.) distribution of the real surface height for the laminar Couette (red dash dot) and turbulent channel (blue dashed) compared to a Gaussian (black solid line) of the same root-mean-square height.

For the sake of brevity, the streamwise mean velocity is denoted by
$U$
where the overline (denoting temporal averaging) and angle brackets (denoting spatial averaging) are dropped, e.g.
The bulk velocity is defined as follows:
Similarly, the Reynolds stresses are denoted by
$u_{i}^{\prime }u_{j}^{\prime }$
dropping the angle brackets and overline. Given that the channel is under a constant pressure gradient, at a statistically stationary state, the average friction wall velocity is given by
$u_{\unicode[STIX]{x1D70F}}=(\unicode[STIX]{x1D6FF}K_{1})^{1/2}$
and the average shear stress by
$\unicode[STIX]{x1D70F}_{w}=\unicode[STIX]{x1D6FF}K_{1}$
. It also holds that
$\unicode[STIX]{x1D70F}_{w}=(\unicode[STIX]{x1D70F}_{w}^{T}+\unicode[STIX]{x1D70F}_{w}^{B})/2$
where
$\unicode[STIX]{x1D70F}_{w}^{T}$
is the top wall shear stress and
$\unicode[STIX]{x1D70F}_{w}^{B}$
the bottom wall shear stress. The top wall is flat therefore
$\unicode[STIX]{x1D70F}_{w}^{T}$
is calculated directly by averaging
$\unicode[STIX]{x1D707}(\unicode[STIX]{x2202}U/\unicode[STIX]{x2202}y)_{y=L_{y}/2}$
and
$\unicode[STIX]{x1D70F}_{w}^{B}$
is calculated indirectly to avoid averaging over the masks using
$\unicode[STIX]{x1D70F}_{w}^{B}=2\unicode[STIX]{x1D70F}_{w}-\unicode[STIX]{x1D70F}_{w}^{T}$
. The bottom wall friction velocity is then calculated using
$u_{\unicode[STIX]{x1D70F}}^{B}=(\unicode[STIX]{x1D70F}_{w}^{B})^{1/2}$
. Results are plotted against the channel height in wall units
$y^{+}=u_{\unicode[STIX]{x1D70F}}y/\unicode[STIX]{x1D708}$
. If the bottom wall friction velocity is used as a reference, then a distinction is made explicitly. For example
$y^{+}(u_{\unicode[STIX]{x1D70F}}^{B})$
denotes the channel height in wall units based on the bottom wall friction velocity.
Two simulations were performed at different resolutions to quantify the effect of grid size. The refined grid is
${\sim}3.5$
as fine as the previous grid. Case T-RI2 was used as a baseline for the grid refinement comparison. No appreciable difference (less than
$1\,\%$
) in the mean velocity profiles, bulk velocity and Reynolds stresses is observed in figure 5. The slip velocity increased by
$3.76\,\%$
and the bottom wall shear stress
$\unicode[STIX]{x1D70F}_{w}^{B}$
decreased by
$2.3\,\%$
. We report results from the finer grid in this paper.
Grid refinement comparison between (a) mean velocity profiles and (b) the Reynolds stresses. Black solid lines represent the fine grid and symbols the coarse grid.

3 Results
3.1 Laminar Couette flow
3.1.1 Steady-state flow field
Initially, the fully wetted Case L-RFW is considered. Figure 6 shows the flow field after it is fully developed. The wall-normal velocity (figure 6
a) and the vorticity magnitude along with surface pressure (figure 6
b) are shown. A wall-normal velocity component into the flow is induced due to the surface asperities. Additional vorticity is generated due to the surface roughness, and large variations of pressure on the surface are evident due to the presence of peaks and valleys. The penetration effect of the surface roughness is illustrated in figure 7, where the percent change in instantaneous streamwise velocity
$(u(y)-u_{o}(y))/u_{o}(y)$
is shown for the fully wetted rough case (Case L-RFW) compared to the smooth channel case (Case L-S) at four wall-parallel planes varying from
$y=0.02H$
to
$y=0.4H$
. The baseline streamwise velocity
$u_{o}(y)$
represents Case L-S and
$u(y)$
represents Case L-RFW. Notice that it is not until the location
$y=0.4H$
that the change in velocity is less than
$1\,\%$
, suggesting that the surface roughness effects can penetrate up to that distance.
Laminar Couette flow (Case L-RFW) with (a) wall-normal velocity contours normalised by the maximum wall-normal velocity
$v_{max}$
and (b) vorticity magnitude line contours normalised by the maximum vorticity
$\unicode[STIX]{x1D714}_{max}$
with surface pressure (normalised by
$p_{max}$
) on the roughness for the range shown in the colour bar.

Percent change in the streamwise velocity field for the fully wetted laminar Couette flow (Case L-RFW) compared to the baseline smooth wall (Case L-S) at wall-parallel planes: (a)
$0.02H$
, (b)
$0.08H$
, (c)
$0.16H$
and (d)
$0.40H$
from the bottom surface.

Laminar Couette flow: drag reduction as a function of interface height normalised by the RMS roughness height
$S_{q}$
.

3.1.2 Mean flow field properties
Simulations are performed for each of the interface heights
$h$
varying from
$S_{v}$
to
$S_{p}$
.
$DR$
is shown for all the interface heights in figure 8. The increase in
$DR$
is not linear when
$h$
is varied from
$S_{v}$
to
$S_{p}$
. Note that the fully wetted case has negative
$DR$
indicating that the absence of the interface has increased drag due the exposed asperities. The presence of the interface produces nearly the same drag reduction for
$h/S_{q}\leqslant -0.32$
which we will refer to as region I. This suggests that the value of
$DR$
is insensitive to
$h$
in region
$I$
. This holds in general for any surface with negative skewness (
$S_{sk}=-0.32$
) since that region
$I$
holds most of the valleys. Figure 9 shows a 2-D slice of the surface roughness; note that the valleys dominate over the peaks about the reference line. As the interface fills up more of the cavities, the slip area becomes significant enough to cause drag reduction. In region
$II$
, for
$-0.32<h/S_{q}<2.15$
, the increase in drag reduction is rapid since the increase in interface height exposes fewer asperities to the outer flow. The slip is enhanced due to a much larger area of air–water interface. In region
$III$
beyond
$h/S_{q}\geqslant 2.15$
,
$DR$
hits a plateau and becomes insensitive to the interface height since it covers most of the asperities.
$DR$
is therefore sensitive to
$h$
in the vicinity of mean roughness
$S_{q}$
. It is therefore evident from figure 8 that we can classify the interface cases into three distinct regions:
$$\begin{eqnarray}\text{Region}=\left\{\begin{array}{@{}ll@{}}I,\quad & \text{if }h/S_{q}\leqslant -0.32;\\ II,\quad & \text{if }-0.32<h/S_{q}<2.15;\\ III,\quad & \text{if }h/S_{q}\geqslant 2.15.\end{array}\right.\end{eqnarray}$$
One can extract a representative case from each of the three regions shown above. The baseline smooth wall is denoted as Case L-S. For the rough wall, one can extract a representative case from each of the three regions shown above. The fully wetted case in region
$I$
is denoted as Case L-RFW, the interface at the mean elevation height of the roughness
$h=0$
denoted by Case L-RI1, the interface at
$h=S_{q}$
in region
$II$
is represented by Case L-RI2 and the interface at
$h=S_{p}$
in region
$III$
is represented by Case L-RI3. The velocity profiles are extracted and compared in figure 10. Case L-RFW exhibits an increase in velocity gradient when compared to Case L-S indicating an increase in drag. Once the interface is introduced, the effect is reversed and the velocity gradient decreases for Cases L-RI1, L-RI2 and L-RI3. The effect is more pronounced in Case L-RI3 since it corresponds to the interface being at the highest peak where most of the asperities are covered and the largest slip effect is achieved.
A 2-D slice of the surface roughness (solid blue line) to highlight negative skewness about the reference line (dashed red line).

Mean streamwise velocity
$u$
as a function of the wall-normal distance
$y$
, where
$u$
and
$h$
are normalised with
$U_{\infty }$
and
$H$
respectively.

Gas fraction
$\unicode[STIX]{x1D719}_{g}$
as a function of interface height
$h$
normalised with the RMS roughness height
$S_{q}$
. The red symbols represent the data and the solid blue line represents the nonlinear fit.

The increase in interface height reduces the amount of rough surface area exposed to the flow. As a result, the flow is subjected to an increase in slip area. The asperities exposed to the outer flow can be represented by a solid fraction
$\unicode[STIX]{x1D719}_{s}$
, which is found by calculating the area of the rough surface above the interface normalised by the projected area of the bottom wall. It is evident that there exists a relationship between the interface height
$h$
and gas fraction
$\unicode[STIX]{x1D719}_{g}$
defined by
$\unicode[STIX]{x1D719}_{g}=1-\unicode[STIX]{x1D719}_{s}$
. This is useful since
$\unicode[STIX]{x1D719}_{g}$
is not known a priori and
$h$
is prescribed as an initial condition. A simple nonlinear fit relates
$\unicode[STIX]{x1D719}_{g}$
to
$h/S_{q}$
of the form of
$c[1+\tanh (ax+b)]$
and is described by the following equation:
This equation can be applied to any general rough surface: the coefficients do change for different surfaces, but the overall fit is general since any surface roughness can be represented by a bearing area curve (BAC). Region
$I$
represents the index of the deepest valleys where the interface is retained, region
$II$
represents the core interface retention index where the maximum amount of air is trapped within the cavities, and region
$III$
the upper zone index related to the largest asperities that contribute to drag. Figure 11 shows a comparison between the actual data and (3.2) for
$\unicode[STIX]{x1D719}_{g}$
as a function of
$h/S_{q}$
. The negative skewness
$S_{sk}=-0.32$
coincides with the transition between regions
$I$
and
$II$
.
Laminar Couette flow: drag reduction
$DR$
as a function of gas fraction
$\unicode[STIX]{x1D719}_{g}$
.

Alternatively,
$DR$
can be represented as a function of slip area instead of interface height by using
$\unicode[STIX]{x1D719}_{g}$
as shown in figure 12. Based on the definition of
$DR$
in (2.7), it can be shown that
$DR$
is related to the slip length
$b_{eff}$
by the following equation:
$$\begin{eqnarray}DR=\frac{1}{\displaystyle 1+\frac{H}{b_{eff}}},\end{eqnarray}$$
therefore,
Philip (Reference Philip1972a
,Reference Philip
b
) obtained an analytic solution for the normalised slip lengths
$b_{eff}/H$
for periodic grooves oriented parallel and perpendicular to the flow respectively:
Ybert et al. (Reference Ybert, Barentin, Cottin-Bizonne, Joseph and Bocquet2007) showed that for a post geometry,
$b_{eff}$
scales with the solid fraction as
where
$\unicode[STIX]{x1D6FC}$
is a prefactor that depends on the geometry. Davis & Lauga (Reference Davis and Lauga2010) were able to use superposition of point sources, where the infinite series is interpreted as a Riemann sum to obtain an analytical solution that agrees with Ybert et al. (Reference Ybert, Barentin, Cottin-Bizonne, Joseph and Bocquet2007) in the asymptotic limit of small surface coverage. A linear regression is performed on the numerically obtained data using the scaling given by (3.7) to obtain the following expression:
Effective slip length
$b_{eff}$
normalised with the channel height
$H$
as a function of (a) the gas fraction
$\unicode[STIX]{x1D719}_{g}$
and (b) the square root of gas fraction area
$\unicode[STIX]{x1D719}_{g}$
. In (a), the solid black line is for longitudinal grooves, the solid dashed line for transverse grooves, the thick solid blue line for post geometry and the symbols for the random roughness. In (b), the solid blue line is the power-law fit and the symbols represent the random roughness.

Figure 13(a) shows a comparison between the different solutions obtained for the longitudinal and transverse grooves and posts to that of random roughness. The solution for post geometry best approximates the data for a random rough geometry for
$\unicode[STIX]{x1D719}_{g}<0.85$
but starts to diverge as gas fraction increases. In the limit of large
$\unicode[STIX]{x1D719}_{g}$
, the transverse groove solution captures the slip effect in random roughness more accurately. We present a power-law fit to the current data to obtain an expression for
$b_{eff}$
as a function of
$\unicode[STIX]{x1D719}_{g}$
:
This formula provides a simple expression for slip length over the rough surface given the amount of gas fraction present when
$\unicode[STIX]{x1D719}_{g}>0$
and is shown in figure 13(b). The fully wetted roughness
$b_{eff}$
remains an outlier.
3.2 Turbulent channel flow
3.2.1 Mean velocity profiles
Mean velocity profiles and Reynolds stresses are computed for five cases (T-S, T-RFW, T-RI1, T-RI2 and T-RI3). Figure 14(a) shows the mean velocity profile
$U$
normalised by the average friction velocity
$u_{\unicode[STIX]{x1D70F}}$
of each corresponding case, as a function of the wall-normal distance
$y$
normalised by the channel half-height
$\unicode[STIX]{x1D6FF}=L_{y}/2$
. The smooth channel (Case T-S) is also shown for reference. The presence of roughness (Case T-RFW) causes a small slip effect: it shows an insignificant decrease of 0.6 % in peak value of
$U/u_{\unicode[STIX]{x1D70F}}$
and about a 3 % shift in its centreline location away from the rough wall. The presence of an interface (Case T-RI1) shows a further increase in slip, a slight increase in the centreline peak value of
$U/u_{\unicode[STIX]{x1D70F}}$
and a 1 % shift in its location towards the SHS wall. The presence of an interface for Case T-RI1 adds just enough slip to offset the effect of roughness. The slip effect is more pronounced for Case T-RI2 and is largest for Case T-RI3 when the interface location covers all the roughness. The mean peak velocity
$U/u_{\unicode[STIX]{x1D70F}}$
increases and the mean profile shifts towards the SHS wall by 5 % and 6 % for Cases T-RI2 and T-RI3 respectively when compared to Case T-S.
Mean profile of (a) velocity normalised with the average friction velocity
$u_{\unicode[STIX]{x1D70F}}$
as a function of wall-normal distance normalised by the channel half-height
$\unicode[STIX]{x1D6FF}$
and (b) close-up near the wall region of the velocity normalised by the bottom wall friction velocity
$u_{\unicode[STIX]{x1D70F}}^{B}$
as a function of wall-normal distance
$y^{+}(u_{\unicode[STIX]{x1D70F}}^{B})$
. Symbols for each case are: Case T-S (▫), Case T-RFW (▵), Case T-RI1 (○), Case T-RI2 (♢), Case T-RI3 (▿). The symbols are not representative of the grid resolution.

Figure 14(b) shows a close-up view near the SHS wall of the mean velocity profile
$U$
normalised by the bottom wall friction velocity
$u_{\unicode[STIX]{x1D70F}}^{B}$
as a function of the wall-normal distance
$y^{+}(u_{\unicode[STIX]{x1D70F}}^{B})$
. The scaling with
$u_{\unicode[STIX]{x1D70F}}^{B}$
describes a more accurate picture in terms of slip and drag reduction. Case T-RFW shows a slip velocity at the wall due to the presence of roughness. The mean velocity profile is 5.8 % lower than the baseline case in the viscous wall region (
$y^{+}<50$
), indicating an overall increase in drag. This effect is not apparent when the mean velocity profile is scaled with the average
$u_{\unicode[STIX]{x1D70F}}$
. Case T-RI1 sees a further increase in slip at the wall due to the presence of an interface: the mean velocity profile is around 14 % higher within the viscous sublayer and extends into the buffer layer where the two velocity profiles of Case T-S and Case T-RI1 intersect at
$y^{+}\sim 22$
. Cases T-RI2 and T-RI3 exhibit the largest slip at the wall: the two profiles intersect at
$y^{+}\sim 28$
. The intersection of the velocity profiles is due to asymmetry caused by the slip effect which shifts the profile towards the SHS. This trend has been observed in Martell et al. (Reference Martell, Perot and Rothstein2009, Reference Martell, Rothstein and Perot2010) for longitudinal grooves and post geometries. In terms of scaling laws, the law of the wall
$u^{+}=y^{+}+C$
still holds (
$C$
is the constant shift that represents the normalised slip velocity), with the exception of Case T-RFW and Case T-RI3. Case T-RFW shows a deviation from the law of the wall at
$y^{+}<3$
, where the Reynolds stress is negligible compared to the viscous stress, but follows it for
$3<y^{+}<10$
. This implies that shifting the profile by the normalised slip velocity
$u_{s}^{+}$
will cause the profile to move further below the baseline case. This behaviour in the profile for the fully wetted case has been demonstrated in the literature (Yuan & Piomelli Reference Yuan and Piomelli2014). Case T-RI3 shows a complete deviation from the law of the wall due to the large amount of slip which modifies the slope of the velocity profile such that
$u^{+}=\unicode[STIX]{x1D6FC}y^{+}+C$
where
$\unicode[STIX]{x1D6FC}<1$
.
Semi-log plot of the mean profile for (a) velocity normalised with the bottom wall friction velocity
$u_{\unicode[STIX]{x1D70F}}^{B}$
of each respective case and (b) velocity shifted by the corresponding slip velocity all normalised by the bottom wall friction velocity
$u_{\unicode[STIX]{x1D70F}}^{B}$
as a function of the wall-normal distance
$y^{+}(u_{\unicode[STIX]{x1D70F}}^{B})$
. Symbols for each case are: Case T-S (▫), Case T-RFW (▵), Case T-RI1 (○), Case T-RI2 (♢), Case T-RI3 (▿). The symbols are not representative of the grid resolution.

A semi-log plot is shown in figure 15(a) where, unlike figure 14(a), the mean velocity profile
$U$
is normalised with the bottom wall friction velocity
$u_{\unicode[STIX]{x1D70F}}^{B}$
. The slip effect is more pronounced in the near-wall region and the difference in peaks are more apparent. The mean velocity profile of Case T-RFW shows a positive slip effect as mentioned earlier: the mean velocity profile is lower than Case T-S over all regions ranging from the viscous sublayer through the log law. The roughness reduces the overall mass flux, indicating an increase in drag. The profile shifts down and away from the wall by roughly 6 % from the centreline location of the baseline. The peak velocity at the centreline is lower than Case T-S by 2 %. Case T-RI1 shows a further increase in slip and the velocity profile is shifted further up indicating a drag reduction. The presence of an interface at that specific height is not enough to overcome the effect of roughness. This is evident from the velocity profile which is still 1 % below the baseline case. The profile shifts back towards the SHS by around 2 % from the centreline location of Case T-RFW. The near-wall slip seems to affect the velocity profile only within the viscous wall region for Cases T-RFW and T-RI1 when compared to Case T-S; however the log-law region shows a collapse in the data where the difference is within 2 %–3 %. Cases T-RI2 and T-RI3 show a large increase in slip near the wall with Case T-RI3 having the largest slip given that the roughness is fully covered by an interface. The peaks for those cases are 8 % and 10 % higher than Case T-S respectively. The profiles for both Cases T-RI2 and T-RI3 shift closer to the SHS wall such that their centrelines are 8.75 % and 11 % away from Case T-S respectively. Overall, the largest difference is clearly seen in the viscous wall region but the log-law region shows a collapse for Cases T-RI2 and T-RI3. This collapse holds until
$y^{+}=40$
and the deviation in the slope of the log-law region in Case T-RI3 becomes more apparent for
$y^{+}>50$
when compared to Case T-RI2. Overall, fully wetted roughness exhibits a decrease in mass flux whereas the presence of an interface increases mass flux. This is evident by the downward shift in the log-law region for Case T-RFW and an upward shift for Cases T-RI1, T-RI2 and T-RI3. We can conclude that the trend for an increase in mass flux (more fluid mass moving) directly correlates with drag reduction and vice versa. This is due to the fact that for all the cases,
$Re_{\unicode[STIX]{x1D70F}}$
and pressure gradient are held constant.
It has been shown that for structured geometries (grooves and posts) not only does the gas fraction
$\unicode[STIX]{x1D719}_{g}$
matter, but also the gap spacing (Ou & Rothstein Reference Ou and Rothstein2005; Daniello et al.
Reference Daniello, Waterhouse and Rothstein2009; Martell et al.
Reference Martell, Perot and Rothstein2009). For a random rough geometry with varying interface heights, the gas fraction is indeed increasing, but the gap spacing is also altered since more pockets are being filled with air. Notice that between Cases T-RI2 and T-RI3 there was not much increase in peak centreline velocity aside from the large increase in slip near the wall which does not substantially alter the total mass flux. However a significant increase in mass flux is observed between Cases T-RI1 and T-RI2. This can be attributed to the change in gap spacing. Although the gas fraction increases with increasing interface height, more gaps and surface valleys are covered up, since the surface roughness is dominated by valleys due to negative skewness. As the interface covers nearly all the valleys then a sharp increase in drag reduction is observed. Case T-RI1 is dominated by small gap size features which may be ineffective in reducing drag; this is in agreement with the literature on longitudinal grooves and posts.
The velocity profiles can be corrected by offsetting them with the slip velocity as shown in figure 15(b). Close to the wall we see a good collapse in the viscous sublayer with an early departure
$y^{+}>1$
from the
$u^{+}=y^{+}$
. With the roughness fully covered by the interface in Case T-RI3, a large deviation from the law of the wall is observed when compared to other cases. This implies that we should expect the structures of the wall-normal turbulence to remain intact for Cases T-RFW, T-RI1 and T-RI2, and Case T-RI3 to be fundamentally different. Away from the near wall, the relative velocity
$(U-u_{s})/u_{\unicode[STIX]{x1D70F}}^{B}$
decreases with roughness and increasing interface height (equivalently with increasing
$\unicode[STIX]{x1D719}_{g}$
). This trend has also been observed in Türk et al. (Reference Türk, Daschiel, Stroh, Hasegawa and Frohnapfel2014) for structured geometries. The profile in the log-law region is given by
where
$\unicode[STIX]{x1D705}$
is the von Kármán constant and
$B$
is the intercept. The value of
$\unicode[STIX]{x1D705}$
decreases from
$0.41$
to
$0.38$
for Case T-RI3 but remains the same for the other cases. The value of
$B$
decreases with increasing
$\unicode[STIX]{x1D719}_{g}$
, going from
$B\approx 5.5$
to
$B\approx 2$
for Case T-RI2 and
$B\approx -4$
for Case T-RI3. Similar trends have been observed in the literature (Busse & Sandham Reference Busse and Sandham2012; Yuan & Piomelli Reference Yuan and Piomelli2014; Busse, Thakkar & Sandham Reference Busse, Thakkar and Sandham2017). The decrease in
$B$
implies an increase in friction as discussed in the literature through surface manipulation (Luchini, Manzo & Pozzi Reference Luchini, Manzo and Pozzi1991; Jiménez Reference Jiménez1994; Garcia-Mayoral & Jiménez Reference Garcia-Mayoral and Jiménez2011). Therefore a decrease in
$B$
is associated with an increase in friction due to roughness and a decrease in friction due to SHS, and in order to differentiate the two, we consider the following argument. This increase in friction that is typically associated with surface roughness is offset by the drag reducing slip velocity
$u_{s}$
due to the presence of an interface in SHS. Roughness induces slip, but the velocity profile in the viscous wall region does not necessarily follow the law of the wall. If we were to compare it to a smooth channel, the profile of the rough channel would intersect the profile of the smooth channel somewhere in the viscous sublayer (
$y^{+}<5$
) to merge into the log region, which is shifted below the baseline indicating an increase in drag. If we offset the profile by the slip velocity, the near-wall region does not collapse but moves further below. If we take the SHS with structured geometries, the interface is flush with the top location of the roughness. The slip effect is also present but the velocity profile in the viscous wall region obeys the law of the wall with some offset
$u^{+}=y^{+}+C$
; therefore if the velocity profile is shifted by the slip velocity, then a collapse in the near-wall region is observed. What we see in our simulation is somewhere in between. This is simply due to the fact that even when we have an interface, some roughness protrudes. Our analysis shows a combination of both behaviours where the law of the wall holds to a certain extent in the viscous sublayer before any appreciable deviation is observed. Also, with increasing interface height, the profile in the viscous wall region tends to move closer to the baseline case and away from the fully wetted roughness, with the exception of Case T-RI3 in which the log-law region seems to extend to the vicinity of the wall, disrupting the near-wall cycle. In general, rough surfaces tend to shift both the law of the wall and the logarithmic region away from baseline while SHS tend to shift the law of the wall closer and simultaneously moving the logarithmic region away from baseline.
3.2.2 Scaling laws
Figure 16(a) shows the slip effect as a function of
$\unicode[STIX]{x1D719}_{g}$
. The slip velocity
$u_{s}$
is normalised by the mean bulk velocity
$U_{b}$
for each corresponding case. As
$\unicode[STIX]{x1D719}_{g}$
increases, the slip effect is more pronounced:
$u_{s}$
exhibits a steady increase. Slip velocities can reach as much as 68 % of the bulk velocities which can obscure some of the effects due to roughness. This is another reason why it would be important to plot the mean profile offset by the slip as shown in figure 15(b).
$DR$
(%) as a function of
$u_{s}$
is not an accurate predictor of drag reduction since fully wetted roughness increases drag while achieving positive slip velocities as shown in figure 16(a). The slip length
$b_{s}$
normalised by the average roughness gap width
$L$
is reported in figure 16(b) as a function of
$\unicode[STIX]{x1D719}_{g}$
. The average roughness gap
$L$
is not a straightforward property to obtain since there is no periodicity (typically associated with structured geometry such as grooves and posts). We use a two-step process to calculate
$L$
. First we find the profile peak count HSC, the number of profile peaks that exceed a pre-selected threshold (e.g. arithmetic mean elevation of the roughness), calculated over the entire streamwise length for all the spanwise slices. Second we calculate the mean peak spacing
$S_{m}$
, the mean spacing between profile peaks, averaged over all the spanwise slices to obtain
$L$
. We refer the reader to figure 28 in appendix B for an example of profile peak count. There is a good agreement between our data and a fit based on the scaling presented by Ybert et al. (Reference Ybert, Barentin, Cottin-Bizonne, Joseph and Bocquet2007) for posts given by the following equation:
Mean flow field properties as a function of
$\unicode[STIX]{x1D719}_{g}$
. (a) Percentage slip velocity
$u_{s}$
normalised by the mean bulk velocity
$U_{b}$
for each corresponding case, (b) slip length
$b_{s}$
normalised by the average gap width of the roughness
$L$
where the solid blue line represents the scaling law of Ybert et al. (Reference Ybert, Barentin, Cottin-Bizonne, Joseph and Bocquet2007) and the dashed black line that of Sbragaglia & Prosperetti (Reference Sbragaglia and Prosperetti2007a
), (c) percentage change in bulk velocity
$\unicode[STIX]{x0394}U_{b}$
normalised by the mean bulk velocity
$U_{b,0}$
of the baseline case, (d) percentage shear stress reduction
$\unicode[STIX]{x0394}\unicode[STIX]{x1D70F}_{w}$
based on the ratio of the bottom wall
$\unicode[STIX]{x1D70F}_{w}^{B}$
to the average channel
$\unicode[STIX]{x1D70F}_{w}$
and (e) correlation between the percentage shear stress reduction and percentage change in bulk velocity compared to a linear fit (solid blue line).

The solution by Sbragaglia & Prosperetti (Reference Sbragaglia and Prosperetti2007a
) which does not require a fit gives good agreement for
$\unicode[STIX]{x1D719}_{g}<0.7$
. Slip length is not a good indicator of
$DR$
since both fully wetted roughness and SHS produce a positive
$b_{s}$
. As mentioned earlier during the discussion of the mean velocity profiles, overall mass flux appears to be a good indicator for drag reduction given that our channels are run at a constant
$Re_{\unicode[STIX]{x1D70F}}$
and pressure gradient. Drag reduction would imply that more fluid mass is moving for the same conditions, i.e. a larger mass flux and a larger change in bulk velocity; the opposite is also true. Figure 16(c) shows the percent change in mean bulk velocity
$\unicode[STIX]{x0394}U_{b}$
normalised by
$U_{b,0}$
of the baseline Case T-S. A negative change in
$\unicode[STIX]{x0394}U_{b}$
indicates an increase in drag, i.e. a lower mass flux. Once an interface is present (Case T-RI1) drag is reduced but is not enough to offset the effect roughness since
$\unicode[STIX]{x0394}U_{b}$
is still negative. As the interface height increases,
$\unicode[STIX]{x0394}U_{b}$
increases to become a positive value indicating an increase in mass flux and a larger
$DR$
effect. This description is exactly what we see when compared to the analysis of the mean velocity profile
$U/u_{\unicode[STIX]{x1D70F}}^{B}$
done earlier.
$\unicode[STIX]{x0394}U_{b}$
does not exhibit a steady increase and plateaus for large values of
$\unicode[STIX]{x1D719}_{g}$
; this was also observed by Türk et al. (Reference Türk, Daschiel, Stroh, Hasegawa and Frohnapfel2014). The reduction in the shear stress on the bottom wall
$\unicode[STIX]{x1D70F}_{w}^{B}$
compared to the average shear stress of the channel
$\unicode[STIX]{x1D70F}_{w}$
provides a straightforward result for
$DR$
as shown in figure 16(d). Note the increase in drag by around 2.5 % for Case T-RFW. The roughness has
$S_{q}^{+}\approx 1.6$
which is in the hydrodynamically smooth regime. Busse et al. (Reference Busse, Thakkar and Sandham2017) reported that the surface property
$S_{z5\times 5}^{+}$
is a more suitable measure of the sand-grain roughness
$k^{+}$
where
$S_{z5\times 5}^{+}\approx 11$
for our surface. This explains why we see a drag increase for the current configuration. A reduction of 1 % in drag is obtained once an interface is introduced in Case T-RI1. A 7 %
$DR$
is achieved for Case T-RI2 and peaks at 8 % for Case T-RI3. The results are in agreement with the experimental results of Ling et al. (Reference Ling, Srinivasan, Golovin, McKinley, Tuteja and Katz2016) obtained for turbulent boundary layers over SHS at higher
$Re_{\unicode[STIX]{x1D70F}}$
, where
$k_{rms}^{+}\approx 0.45$
–
$0.75$
and
$DR$
ranged from 9 % to 12 %. One can therefore reasonably assume that the change in shear stress and the change in bulk velocity correlate well with each other. Figure 16(e) shows
$\unicode[STIX]{x0394}\unicode[STIX]{x1D70F}_{w}$
as a function of
$\unicode[STIX]{x0394}U_{b}$
compared to a linear fit given by
The roughness function
$\unicode[STIX]{x0394}U^{+}$
depends on the three-dimensional topography of a surface where Bradshaw (Reference Bradshaw2000) suggested a scaling given by
$\unicode[STIX]{x0394}U^{+}\sim (k_{rms}^{+})^{\unicode[STIX]{x1D6FC}}$
. Since the interface height essentially modifies
$k_{rms}^{+}$
, then a similar scaling argument can be made for
$\unicode[STIX]{x0394}U^{+}$
as a function of
$b_{s}^{+}$
. Figure 17(a) shows a good agreement with a power-law behaviour where the roughness function goes to zero as the slip length diminishes. The power-law formula obtained from fitting the data is given by
A semi-analytical formula by Ybert et al. (Reference Ybert, Barentin, Cottin-Bizonne, Joseph and Bocquet2007) relates the slip velocity to the cavity width, in terms of wall units
$u_{s}^{+}=C_{Y}L^{+}$
. We simply use the slip length
$b_{s}^{+}$
instead since slip is a direct consequence of
$L^{+}$
. Figure 17(b) shows the linear approximation near the wall for
$C_{Y}=0.65$
where
The linear scaling deviates for
$b_{s}^{+}>10$
which was also shown for
$L^{+}>10$
in Seo et al. (Reference Seo, Garcia-Mayoral and Mani2015) for a value of
$C_{Y}=0.535$
. This is expected since the linear relationship is based on Stokes’ flow which becomes less accurate as slip increases. It is worth noting that by definition
$u_{s}^{+}=b_{s}^{+}$
for the cases where the air–water interface is aligned with the top of the roughness and where the total stress is equal to the viscous stress. However, as discussed in Ling et al. (Reference Ling, Srinivasan, Golovin, McKinley, Tuteja and Katz2016),
$u_{s}^{+}<b_{s}^{+}$
for the cases when
$-(u^{\prime }v^{\prime })^{+}>0$
or when the total stress
$\unicode[STIX]{x1D70F}_{t}^{+}>\unicode[STIX]{x1D707}(\text{d}U/\text{d}y)^{+}$
due the random nature of the roughness and the variable interface heights.
Roughness function
$\unicode[STIX]{x0394}U^{+}$
as a function of slip length
$b_{s}^{+}$
is shown in (a) and is compared to a power-law fit (solid blue line) and (b) shows the slip velocity scaling
$u_{s}^{+}$
as a function of
$b_{s}^{+}$
with a comparison to a linear approximation (solid black line).

3.2.3 Reynolds stresses
The Reynolds stress profiles are shown for the different cases as a function of the wall-normal distance. Figure 18(a) shows the streamwise component of the Reynolds stress normalised by
$u_{\unicode[STIX]{x1D70F}}$
. Case T-RFW seems to slightly increase the peak of
$u^{\prime }u^{\prime }$
while simultaneously shifting it away from the SHS wall. The presence of an interface seems to dampen that effect where Case T-RI1 shows a decrease in the peak while shifting it towards the SHS wall. This effect is further amplified in Case T-RI2 where the peak is clearly damped and the profile moves closer to the wall. Case T-RI3 exhibits the largest slip and hence is closest to the SHS wall where there exists a sharp rise in the streamwise component; however the peak is larger than the other cases, breaking the symmetry of the profile completely. The overall trend in damping the peak of
$u^{\prime }u^{\prime }$
while shifting it closer to the slip wall has been observed in the literature for longitudinal grooves and posts. For reference,
$u^{\prime }u^{\prime }$
is scaled by
$u_{\unicode[STIX]{x1D70F}}^{B}$
in figure 18(b). We get a reasonable collapse in the profiles near the SHS wall with the exception of Case T-RI3 indicating a change in the near-wall behaviour. Figure 18(c) shows a
$\log$
plot of
$u^{\prime }u^{\prime }$
as a function of the wall-normal distance
$y^{+}(u_{\unicode[STIX]{x1D70F}}^{B})$
in wall units scaled by
$u_{\unicode[STIX]{x1D70F}}^{B}$
. Large differences are observed in the near-wall region due to the slip effect. The shift away from the SHS wall for Case T-RFW and the shift towards the SHS wall for Cases T-RI1 and T-RI2 are more evident.
Streamwise component of the normal Reynolds stress normalised by (a) the average friction velocity
$u_{\unicode[STIX]{x1D70F}}^{2}$
as a function of
$y/\unicode[STIX]{x1D6FF}$
, (b) the bottom wall friction velocity
$(u_{\unicode[STIX]{x1D70F}}^{B})^{2}$
as a function of
$y/\unicode[STIX]{x1D6FF}$
and (c) the bottom wall friction velocity
$(u_{\unicode[STIX]{x1D70F}}^{B})^{2}$
as a function of
$y^{+}(u_{\unicode[STIX]{x1D70F}}^{B})$
on a
$\log$
scale. Symbols for each case are: Case T-S (▫), Case T-RFW (▵), Case T-RI1 (○), Case T-RI2 (♢), Case T-RI3 (▿). The symbols are not representative of grid resolution.

Wall-normal component of the normal Reynolds stress normalised by (a) the average friction velocity
$u_{\unicode[STIX]{x1D70F}}^{2}$
as a function of
$y/\unicode[STIX]{x1D6FF}$
, (b) the bottom wall friction velocity
$(u_{\unicode[STIX]{x1D70F}}^{B})^{2}$
as a function of
$y/\unicode[STIX]{x1D6FF}$
and (c) the bottom wall friction velocity
$(u_{\unicode[STIX]{x1D70F}}^{B})^{2}$
as a function of
$y^{+}(u_{\unicode[STIX]{x1D70F}}^{B})$
on a
$\log$
scale. Symbols for each case are: Case T-S (▫), Case T-RFW (▵), Case T-RI1 (○), Case T-RI2 (♢), Case T-RI3 (▿). The symbols are not representative of grid resolution.

Spanwise component of the normal Reynolds stress normalised by (a) the average friction velocity
$u_{\unicode[STIX]{x1D70F}}^{2}$
as a function of
$y/\unicode[STIX]{x1D6FF}$
, (b) the bottom wall friction velocity
$(u_{\unicode[STIX]{x1D70F}}^{B})^{2}$
as a function of
$y/\unicode[STIX]{x1D6FF}$
and (c) the bottom wall friction velocity
$(u_{\unicode[STIX]{x1D70F}}^{B})^{2}$
as a function of
$y^{+}(u_{\unicode[STIX]{x1D70F}}^{B})$
on a
$\log$
scale. Symbols for each case are: Case T-S (▫), Case T-RFW (▵), Case T-RI1 (○), Case T-RI2 (♢), Case T-RI3 (▿). The symbols are not representative of grid resolution.

The wall-normal Reynolds stress component
$v^{\prime }v^{\prime }$
is shown in figure 19. Initially
$v^{\prime }v^{\prime }$
is normalised by
$u_{\unicode[STIX]{x1D70F}}$
as shown in figure 19(a). Case T-RFW amplifies the peak stress and shifts it away from the slip wall. Comparing it to Case T-RI2, it is evident that the interface has a damping effect. It is interesting to note that once the interface was introduced initially, Case T-RI1 showed a further amplification in the peak from Case T-RFW and not the opposite. This is likely to be due to the fact that the interface height barely covers any of the roughness, and therefore the inhomogeneity between slip and no-slip due to the interface and random rough patches cause larger fluctuations. This is of course damped out once the interface covers more of the surface and larger slip areas are present. Case T-RI3 shows the largest shift towards the SHS wall as expected but the peak stress does not not seem to follow any further damping effect with increasing slip area. Figure 19(b) shows
$v^{\prime }v^{\prime }$
scaled by
$u_{\unicode[STIX]{x1D70F}}^{B}$
and no collapse in the data is observed. The
$\log$
plot in figure 19(c) shows the velocities going to zero since at the interface an infinite surface tension is assumed. Therefore the wall-normal velocity
$v_{N}=0$
at the interface and a no-slip boundary condition is applied elsewhere over the rough surface. Note that the profiles of Case T-RI1 and Case T-RI2 intersect around
$y^{+}\approx 28$
where Case T-RI2 exhibits larger wall-normal stresses for
$y^{+}<28$
and then Case T-RI1 tends to become larger for the region above that. Case T-RFW does not exhibit an observable shift away from the wall while it remains clear that when an interface is present, the shift towards the SHS wall remains considerable.
The streamwise component of Reynolds stress is shown in figure 20, which follows a similar trend in the behaviour to that described above. Case T-RFW shows an increase in
$w^{\prime }w^{\prime }$
when scaled by
$u_{\unicode[STIX]{x1D70F}}$
as shown in figure 20(a). The addition of an interface (Case T-RI1) further increases the peak in
$w^{\prime }w^{\prime }$
. The same reasoning applies to the spanwise component as described earlier for the wall-normal component since the roughness has no preferential direction. Therefore the interplay between slip and no-slip due to the interface and the protruding roughness holds here too. Similarly as soon as the interface covers a large portion of the roughness as in Case T-RI2, the peak is damped out. Case T-RI3 is highly skewed towards the rough wall with a peak velocity that is much larger than the rest of the cases. Similar to what was observed in the wall-normal stresses, Cases T-RI1 and T-RI2 intersect each other as shown in figure 20(c). The location however is smaller, where
$y^{+}\approx 13$
. Anywhere below that, Case T-RI2 exhibits higher stresses than Case T-RI1 and the opposite holds true when
$y^{+}>13$
. The spanwise slip is evident due to the presence of the interface where Case T-RI3 still exhibits the largest slip effect.
The Reynolds shear stress is an important quantity to examine since its behaviour is closely related to turbulence levels and the structure of the near-wall turbulence. Figure 21(a) shows how the presence of an interface (Case T-RI2) reduces shear while the presence of roughness enhances it. Note that since the wall-normal shear stress tends to increase when an interface is introduced (Case T-RI1), that effect translates here where we see a further enhancement in mixing instead of a reduction in the peak. Similar to the discussion above, it is not until the interface covers a large portion of the roughness that the damping effect takes place. Therefore it is evident that there are competing effects between the interface suppressing vertical velocity fluctuations and the asperities doing the opposite by enhancing them. Case T-RI3 exhibits a large gradient near the SHS wall since the wall-normal velocity is zero yet the streamwise component sees a large slip effect as was shown in figure 18(a). As we move away from the wall, the shear stress has to balance out with the top wall given that we are running a constant pressure gradient, which explains why the profiles are parallel in that region. Figure 21(b) shows the
$u^{\prime }v^{\prime }$
component normalised by
$u_{\unicode[STIX]{x1D70F}}^{B}$
, and a good collapse of the data is observed with the exception of Case T-RI3. This indicates that, overall, the near-wall turbulence is not fundamentally changed for Cases T-RFW, T-RI1 and T-RI2 whereas the turbulent structures of Case T-RI3 are different. Figure 21(c) shows
$-u^{\prime }v^{\prime }$
plotted on a
$\log$
scale in wall units (normalised by
$u_{\unicode[STIX]{x1D70F}}^{B}$
) along with the viscous stress
$\unicode[STIX]{x1D707}(\text{d}U/\text{d}y)$
and total stress
$\unicode[STIX]{x1D70F}_{t}$
which shows the total shear stress budget. We can see the cases where
$-u^{\prime }v^{\prime }>0$
and where
$\unicode[STIX]{x1D70F}_{t}>\unicode[STIX]{x1D707}(\text{d}U/\text{d}y)$
leading to
$u_{s}^{+}<b_{s}^{+}$
as discussed earlier. Case RFW seems to collapse onto the baseline while Cases T-RI1 and T-RI2 shift towards the slip wall. The shear stresses go to zero near the wall due to the infinite surface tension that keeps the interface flat.
Reynolds shear stress normalised by (a) the average friction velocity
$u_{\unicode[STIX]{x1D70F}}^{2}$
as a function of
$y/\unicode[STIX]{x1D6FF}$
, (b) the bottom wall friction velocity
$(u_{\unicode[STIX]{x1D70F}}^{B})^{2}$
as a function of
$y/\unicode[STIX]{x1D6FF}$
. (c) The viscous stress, Reynolds shear stress and total stress normalised by the bottom wall friction velocity
$(u_{\unicode[STIX]{x1D70F}}^{B})^{2}$
as a function of
$y^{+}(u_{\unicode[STIX]{x1D70F}}^{B})$
on a
$\log$
scale. Symbols for each case are: Case T-S (▫), Case T-RFW (▵), Case T-RI1 (○), Case T-RI2 (♢), Case T-RI3 (▿). In (c) the solid line represents total stress
$\unicode[STIX]{x1D70F}_{t}^{+}$
, dashed lines the viscous stress
$\unicode[STIX]{x1D707}(\text{d}U/\text{d}y)^{+}$
and dashed-dotted lines the Reynolds shear stress
$-(u^{\prime }v^{\prime })^{+}$
. The symbols are not representative of grid resolution.

Earlier we discussed how
$u_{s}$
is not an accurate predictor of
$DR$
since both fully wetted roughness and SHS result in a positive slip velocity.
$\unicode[STIX]{x0394}U_{b}$
on the other hand correlates with
$\unicode[STIX]{x0394}\unicode[STIX]{x1D70F}_{w}$
. One can show that the Reynolds shear stress is tied to the change in bulk velocity which is manifested in additional turbulent losses. This is done by applying a triple integration to the averaged transport equation for the streamwise momentum equation. This was demonstrated in Hasegawa, Frohnapfel & Kasagi (Reference Hasegawa, Frohnapfel and Kasagi2011) and Türk et al. (Reference Türk, Daschiel, Stroh, Hasegawa and Frohnapfel2014) where the following identities are obtained:
For the baseline case (T-S) where we have no-slip walls:
This leads to the final form given by the following:
where
$u_{0}^{\prime }v_{0}^{\prime }$
denotes the shear stress of the baseline case. Take Case T-RFW as an example. The fully wetted roughness enhances vertical velocity fluctuations and so does having an interface at a small height location (Case T-RI1). There
$-u^{\prime }v^{\prime }>-u_{0}^{\prime }v_{0}^{\prime }$
and therefore
$\int _{0}^{\unicode[STIX]{x1D6FF}}(1-y/\unicode[STIX]{x1D6FF})(-u^{\prime }v^{\prime }+u_{0}^{\prime }v_{0}^{\prime })>0$
which happens also to be larger than the
$u_{s}$
caused by the roughness or the presence of the interface. This gives
$\unicode[STIX]{x0394}U_{b}<0$
, indicating an increase in drag. For Cases T-RI2 and T-RI3,
$u^{\prime }v^{\prime }$
is damped when compared to baseline. Therefore
$-u^{\prime }v^{\prime }<-u_{0}^{\prime }v_{0}^{\prime }$
and the integral term ends up coming out to be negative which results in
$\unicode[STIX]{x0394}U_{b}>0$
and hence drag reduction. One implication of this result is that a superhydrophobic surface might fail in the sense of reducing drag at high pressure in spite of the interface itself being stable.
3.2.4 Flow structures
The near-wall shear stress
$u^{\prime }v^{\prime }$
collapsed when scaled with the bottom wall shear velocity indicating that the near-wall behaviour remains unchanged with the exception of Case T-RI3. In this section we look more closely at the near-wall region. Figure 22(a–d) shows instantaneous pressure fluctuations in wall units where
$p^{+}=p/(u_{\unicode[STIX]{x1D70F}}^{B})^{2}$
. Figure 22(a) shows the pressure contours along the rough wall while figures 22(b)–22(d) are taken at the interface location (the reference plane) for each case. The axes are shown in wall units scaled by
$u_{\unicode[STIX]{x1D70F}}^{B}/\unicode[STIX]{x1D708}$
. The contours change noticeably from figures 22(a) to 22(d). The instantaneous behaviour in pressure for figures 22(a) and 22(d) shows similar large length scales for fully wetted and fully covered roughness on the order of
$100\unicode[STIX]{x1D708}/u_{\unicode[STIX]{x1D70F}}$
. Larger pressure intensities are visible in figure 22(d) due to slip at the interface. Local contribution due to the presence of protruding roughness elements is observed in figures 22(b) and 22(c). Asperities cause stagnation in front of them as observed by the high pressure values. Low pressure values are seen in their wake. This effect is reduced as more asperities are covered by the interface as shown in figure 22(c). This behaviour was also observed in Seo et al. (Reference Seo, Garcia-Mayoral and Mani2015).
Instantaneous contours of
$p^{+}$
for (a) Case T-RFW on the roughness, (b) Case T-RI1, (c) Case T-RI2 and (d) Case T-RI3 on the interface location. Range of the contours is from
$-5$
to
$5$
in wall units.

Time-averaged contours of
$p^{+}$
at
$y^{+}=15$
for (a) Case T-RFW, (b) Case T-RI1, (c) Case T-RI2 and (d) Case T-RI3. Range of contours is from
$-2$
to
$2$
in wall units.

Instantaneous greyscale contours of
$u/U_{b}$
at
$y^{+}=15$
for (a) Case T-RFW, (b) Case T-RI1, (c) Case T-RI2 and (d) Case T-RI3. Range of contours is from
$0.2$
to
$1$
.

We extract slices in the wall-normal plane such that
$y^{+}=15$
and examine the time-averaged pressure distribution. Figure 23(a–d) shows striking visual difference between Case T-RFW, Cases T-RI1 and T-RI2, and Case T-RI3. The mean pressure is zero for a smooth channel; this is only observed in figure 23(d) when the interface fully covers the roughness. Figure 23(a) shows a variation in the mean pressure due to the signature of the rough wall. The contribution of stagnation pressure at the interface is overlaid with the contribution of asperities. Large variations in mean pressure interspersed across the domain are observed in figures 23(b) and 23(c) due to the interaction between the flow at the interface and around the asperities.
Since the pressure footprint can be related to near-wall velocity, we plot the instantaneous velocity
$u$
normalised by
$U_{b}$
in greyscale on the same
$y$
-slice. The presence of an interface affects the spanwise streak motion, and the distance between streaks can be visually seen to be around
$100$
wall units in figure 24(a). Figure 24(b) shows regions with larger streamwise velocities where the separation distance between streaks becomes larger than
$100$
wall units in the spanwise direction. The slightly larger distance between the streaks is maintained in Case T-RI2 but that coherent structure is completely destroyed in Case T-RI3 as shown in figure 24(d). This behaviour for Case T-RI3 is likely to be due to the lower near-wall shear yielding lower values of
$Sq^{2}/\unicode[STIX]{x1D716}$
where
$S$
,
$q^{2}$
and
$\unicode[STIX]{x1D716}$
are the mean shear, twice the turbulence kinetic energy and turbulent dissipation respectively. As shown by Rogers & Moin (Reference Rogers and Moin1987) and Lee, Kim & Moin (Reference Lee, Kim and Moin1990) for homogeneous shear flow, only high values of
$Sq^{2}/\unicode[STIX]{x1D716}$
as encountered in the near-wall region produce streaks.
Instantaneous contours of
$\unicode[STIX]{x1D714}_{y}^{+}$
at
$y^{+}=15$
for (a) Case T-RFW, (b) Case T-RI1, (c) Case T-RI2 and (d) Case T-RI3. Range of contours is from
$-0.5$
to
$0.5$
in wall units.

Illustration of the near-wall vortex structure. The asymmetry in slip due to the random roughness is illustrated in (a) whereas (b) shows the effect of having an interface cover most of the roughness and its effect on slip and induced secondary flow interaction.

The instantaneous wall-normal vorticity
$\unicode[STIX]{x1D714}_{y}^{+}$
is shown in figure 25(a–d) to illustrate pairs of counter-rotating vortices in regions of low- and high-momentum fluid streaks. The
$\unicode[STIX]{x1D714}_{y}^{+}$
contours show similar behaviour to the streamwise velocity with varying interface height. A striking feature in the
$\unicode[STIX]{x1D714}_{y}^{+}$
plots is the slight asymmetry of counter-rotating vortex pairs. Figure 26(a) schematically shows how, due to the random nature of the slip no-slip behaviour over rough SHS, each vortex in a counter-rotating vortex pair experiences different slip areas, producing an inhomogeneous spanwise slip
$w_{s}$
. This causes asymmetric velocity profiles to interact in the wall-normal direction. As the interface height increases, the slip area increases and less of the solid protrudes, and therefore
$w_{s}$
becomes more homogeneous. The slip effect is however more amplified in this case (figure 26
b) and the spanwise velocity is much larger, penetrating further up in the wall-normal direction, with tertiary vortices set up above the pair of counter-rotating vortices. This could possibly explain why in figures 24(d) and 25(d) we observe less coherence and more violent mixing of those structures. While such mixing is associated with drag increase, the slip velocity in the streamwise direction is more dominant which offsets this deleterious effect to give a net positive
$DR$
. This is explained by the change in bulk velocity as a function of slip velocity and the integral difference of shear stresses given by (3.17).
4 Summary
DNS of laminar Couette flow at
$Re=740$
and turbulent channel flow at
$Re_{\unicode[STIX]{x1D70F}}=180$
are performed, where the bottom wall is a realistically rough SHS. The surface scan is reproduced computationally, and the surface statistics are verified again with the experiments. Simulations are also performed for a smooth wall to serve as a baseline. Simulations of the fully wetted case and an air–water interface at various heights are compared to the smooth channel. The effects of roughness and interface heights are discussed in detail.
Simulations of laminar Couette flow show a penetration effect up to
$40\,\%$
in the wall-normal direction due to the roughness. Various interface heights
$h$
were considered and a nonlinear dependence of drag reduction
$DR$
on
$h$
is observed. The dependence can be categorised into three distinct regions. The drag is sensitive to the interface location in region
$II$
described in § 3.1.2 where
$h$
is in the range of
$-0.32<h/S_{q}<2.15$
. The negative skewness
$S_{sk}=-0.32$
of the roughness profile indicates that the surface contains more valleys than peaks and asperities. More than half of the surface roughness is filled with gas when
$h$
is within the vicinity of
$S_{q}$
since a large number of valleys become wetted. The solid fraction
$\unicode[STIX]{x1D719}_{s}$
decreases and the gas fraction
$\unicode[STIX]{x1D719}_{g}$
increases with increasing
$h$
where fewer asperities are exposed to the outer flow. Therefore it can be shown that
$DR$
is a function of
$\unicode[STIX]{x1D719}_{g}$
. To demonstrate the relation between
$h$
and
$\unicode[STIX]{x1D719}_{g}$
, we first calculate
$\unicode[STIX]{x1D719}_{s}$
by measuring the amount of non-wetted area above the interface. This is similar to the bearing area curve (BAC). Based on the definition
$\unicode[STIX]{x1D719}_{g}=1-\unicode[STIX]{x1D719}_{s}$
, the gas fraction is obtained and plotted as a function of
$h$
. A nonlinear fit is given by
$\unicode[STIX]{x1D719}_{g}=0.5[1+\tanh (0.95(h/S_{q})-0.875)]$
and shows good agreement with the data. This is useful since
$\unicode[STIX]{x1D719}_{g}$
is not known a priori and
$h$
is prescribed as an initial condition. Effective slip can be directly related to
$DR$
using the definition
$b_{eff}/H=DR/(1-DR)$
. A power law using a linear regression fit is obtained and corrected for the fully wetted case. The relation is given by
$b_{eff}/H=0.5(\unicode[STIX]{x1D719}_{g})^{5/2}+0.02$
which shows a good agreement with the data and provides a useful model for the slip length given
$\unicode[STIX]{x1D719}_{g}$
. The results are compared to previous work done on structured geometry. It is observed that random rough surfaces behave like post geometries for
$\unicode[STIX]{x1D719}_{g}<0.85$
and like transverse grooves in the upper limit.
Based on the observations of the three distinct regions made in the laminar Couette flow, four interface heights are chosen for turbulent DNS channel flow. Simulations of a fully wetted rough case
$h=S_{v}$
and the interface heights at
$h=0$
,
$S_{q}$
and
$S_{p}$
are performed for turbulent channel flow and the results are discussed. The mean velocity profile is normalised by two quantities: the average channel wall friction velocity
$u_{\unicode[STIX]{x1D70F}}$
and the bottom wall friction velocity
$u_{\unicode[STIX]{x1D70F}}^{B}$
. The mean velocity profiles show the effect of roughness where a reduction in mass flux is obtained. The presence of an interface increases mass flux. The velocity profiles are offset by the slip velocities and show a good agreement in the law of the wall where the data collapse. Case T-RI3 exhibits an early departure from the law of the wall indicating a change in the turbulent structures. Various mean flow properties are extracted and plotted, where we show that
$u_{s}$
can become a large fraction of the bulk velocity
$U_{b}$
(up to 68 %), and the slip length
$b_{s}$
maintains the scaling law proposed by Ybert et al. (Reference Ybert, Barentin, Cottin-Bizonne, Joseph and Bocquet2007). We show that
$u_{s}$
is not a good indicator of
$DR$
by itself since roughness induces a positive
$u_{s}$
while increasing drag. A more reasonable quantity to describe
$DR$
is the change in bulk velocity
$\unicode[STIX]{x0394}U_{b}$
since it implicitly contains the information from
$u_{s}$
and the additional turbulent losses in the form of a weighted Reynolds shear stress
$\unicode[STIX]{x0394}U_{b}=u_{s}-\int _{0}^{\unicode[STIX]{x1D6FF}}(1-y/\unicode[STIX]{x1D6FF})(-u^{\prime }v^{\prime }+u_{0}^{\prime }v_{0}^{\prime })\,\text{d}y$
. This has implications for
$DR$
where the surface might fail in reducing drag although the interface itself is stable for high pressure.
$\unicode[STIX]{x0394}U_{b}$
shows a good correlation with
$\unicode[STIX]{x0394}\unicode[STIX]{x1D70F}_{w}$
. Although
$u_{s}$
and
$b_{s}$
continually increase with increasing
$\unicode[STIX]{x1D719}_{g}$
, the change in wall shear stress
$\unicode[STIX]{x0394}\unicode[STIX]{x1D70F}_{w}$
plateaus at large gas fraction and so does the
$\unicode[STIX]{x0394}U_{b}$
. We discuss scaling laws for these quantities and correlate them with each other.
The Reynolds stresses are also examined showing an overall behaviour consistent with previous work on structured geometries. Slip tends to shift the profiles towards the SHS wall whereas roughness pushes it away from wall. Asperities enhance negative shear stress and therefore momentum mixing, while the interface suppresses them. The streamwise Reynolds stress
$u^{\prime }u^{\prime }$
and Reynolds shear stress
$u^{\prime }v^{\prime }$
show a good collapse in the data when normalised by
$u_{\unicode[STIX]{x1D70F}}^{B}$
indicating that the near-wall turbulence remains fundamentally unchanged with the exception of Case T-RI3. This prompted a further investigation into the nature of the near-wall turbulent structures. We looked at instantaneous pressure contours in the near-wall region and at
$y^{+}=15$
where time-averaged pressure
$p^{+}$
, instantaneous streamwise velocity
$u/U_{b}$
and wall-normal vorticity
$\unicode[STIX]{x1D714}_{y}^{+}$
are examined.
Pressure fluctuations in the near-wall region exhibit a competing effect between large-scale turbulent fluctuations and a contribution due to stagnation pressure in front of the asperities. The instantaneous behaviour in pressure shows similar large-length-scale fluctuations for Case T-RFW and Case T-RI3 on the order of
$100\unicode[STIX]{x1D708}/u_{\unicode[STIX]{x1D70F}}$
. Larger pressure intensities are visible in Case T-RI3 due to the slip effect at the interface. Local contribution due to the presence of protruding roughness elements is observed in Case T-RI1 and Case T-RI2. Asperities cause stagnation in front of them as observed by the high pressure values. Low pressure values are seen in their wake. This effect is reduced as more asperities are covered by the interface. Time-averaged pressure fluctuations show that Case T-RI3 resembles a smooth channel since the asperities are completely covered leading to a zero mean pressure variation. This is not observed for Cases T-RFW, T-RI1 and T-RI2 where the effect of asperities and the interface are clearly seen as large variations in mean pressure interspersed across the domain. At
$y^{+}=15$
, Case T-RI3 does indeed alter the near-wall turbulence where we see a complete loss of coherent streaks, as observed from the pressure fluctuations, streamwise velocity and wall-normal vorticity. A physical mechanism is proposed to explain the observed trends in flow structure.
Acknowledgements
This work was supported by the United States Office of Naval Research (ONR) MURI (Multidisciplinary University Research Initiatives) program under grant N00014-12-1-0874 managed by Dr K.-H. Kim. Computing resources were provided by the Minnesota Supercomputing Institute (MSI). We are grateful to Professor W. Choi at University of Texas Dallas and Professor G. H. McKinley at MIT for providing us with the scanned surface data used in the present work. The authors would like to thank Dr P. Kumar and Dr Y. Li for their helpful discussions and suggestions.
Appendix A. Validity of assumptions
In practice, for the superhydrophobic surface to sustain its drag reducing properties, the surface tension must be strong enough to maintain the presence of an air–water interface. This implies that the capillary pressure must be larger than the background turbulent pressure fluctuations. The balance between surface tension and the external pressure results in a meniscus shape and a contact angle at the wall contact boundaries. In the study, the interface is assumed to be flat, which it may not be, and that the interface is always sustained. We investigate the range of validity of our assumptions by using scaling arguments of the driving mechanisms in interfacial physics and comparing their orders of magnitude. The asterisk is used to denote dimensional quantities.
A.1 Small interface deflection approximation
Let
$s^{\ast }$
represent the interface deflection and
$w^{\ast }$
the average cavity width of the rough surface. The Young–Laplace equation gives
where
$\unicode[STIX]{x0394}p_{c}^{\ast }$
is the capillary pressure across the interface. Assuming the interface is pinned at the contact points of the cavity width, then
$R^{\ast }$
is the radius of the interface. We can then relate
$w^{\ast }$
to
$s^{\ast }$
given that
$w^{\ast }$
represents the chord of a circular segment such that
$w^{\ast }=2\sqrt{s^{\ast }(2R^{\ast }-s^{\ast })}$
. Substitute for
$R^{\ast }$
using (A 1) to obtain the following relation:
For a flat interface,
$s^{\ast }/w^{\ast }\ll 1$
where
$s^{\ast }/w^{\ast }$
represents the ratio of interface deflection to cavity width. Assume the maximum deflection to be no larger than 10 % such that the maximum deflection (contact) angle is less than
${\sim}3^{\circ }$
so we obtain
which gives
Therefore the maximum sustained pressure given a cavity width is
A.2 Interface stability approximation
In a realistic environment, the turbulent pressure fluctuations play an important role in determining whether the interface breaks or remains intact. In order for the surface to maintain its drag reducing properties, capillary pressure must be strong enough to maintain the air–water interface and overcome turbulent pressure fluctuations. Using similar scaling arguments as before, we know that the turbulent pressure fluctuations scale as follows:
and the capillary pressure as
In a stable configuration,
$\unicode[STIX]{x0394}p_{c}^{\ast }\gg p_{rms}^{\ast }$
must be satisfied. Therefore we obtain the following relation:
The above equation can be rearranged such that
which gives an upper bound on the friction velocity
Therefore in terms of
$Re_{\unicode[STIX]{x1D70F}}$
,
A.3 Range of validity
For the following analysis, we take water as a reference fluid at standard conditions:
$\unicode[STIX]{x1D70C}^{\ast }=997~\text{kg}~\text{m}^{-3}$
,
$\unicode[STIX]{x1D707}^{\ast }=8.94\times 10^{-4}~\text{Pa}~\text{s}$
and
$\unicode[STIX]{x1D70E}^{\ast }=7.2\times 10^{-2}~\text{N}~\text{m}^{-1}$
. In our numerical simulation, the surface
$S_{q}$
is approximately
$1/90$
th of the channel half-height
$\unicode[STIX]{x1D6FF}^{\ast }$
which gives
$\unicode[STIX]{x1D6FF}^{\ast }\sim O(10^{-4}\text{m})$
. From a design perspective, there exists a top down approach (the present study) where a surface is given and we estimate the range of validity of
$Re_{\unicode[STIX]{x1D70F}}$
. In a bottom up approach, we can find the upper limit of the maximum allowable
$w^{\ast }$
that sustains an interface given an
$Re_{\unicode[STIX]{x1D70F}}$
.
A.3.1 Top down approach
In our numerical experiment,
$w^{\ast }$
is of
$O(10~\unicode[STIX]{x03BC}\text{m})$
therefore (A 5) yields
$\unicode[STIX]{x0394}p_{c}^{\ast }<5.7~\text{kPa}$
suggesting that the interface can sustain pressures up to that value before the assumption of flat interfaces breaks down. Set the calculated pressure as the upper limit for
$p_{rms}^{\ast }$
and substitute (A 6) in (A 5) to obtain
In terms of
$Re_{\unicode[STIX]{x1D70F}}$
we have
Therefore, the assumption of a flat interface is valid for
$Re_{\unicode[STIX]{x1D70F}}<270$
. For the assumption of a stable interface, we use (A 11) to obtain
$Re_{\unicode[STIX]{x1D70F}}\ll 300$
. It is clear from these results that the assumption of a flat interface puts a more stringent requirement on the allowable
$Re_{\unicode[STIX]{x1D70F}}$
which can also be seen by comparing (A 11) to (A 13).
A.3.2 Bottom up approach
Given a range of
$Re_{\unicode[STIX]{x1D70F}}$
, we can estimate the largest allowable cavity width between roughness peaks. It is helpful to define terms in wall units such that (A 11) is rewritten as
Therefore
$w^{+}$
for interface stability is given as
The maximum sustained capillary pressure can also be written in wall units,
where
$Ca=\unicode[STIX]{x1D707}^{\ast }u_{\unicode[STIX]{x1D70F}}^{\ast }/\unicode[STIX]{x1D70E}^{\ast }$
is the ratio of viscous to capillary stresses known as the capillary number. For small interface deflections,
$w^{+}$
is therefore
An example of such bottom up calculation is given in table 3.
Power spectral density (PSD) contour plot where
$k_{x}$
and
$k_{z}$
are the wavenumbers in the streamwise and spanwise directions respectively.

Maximum allowable average cavity widths in wall units
$w^{+}$
for a range of
$Re_{\unicode[STIX]{x1D70F}}$
. Here
$w_{IS}^{+}$
and
$w_{SD}^{+}$
represent the average cavity width satisfying the interface stability and small deflection conditions respectively. The maximum allowable capillary pressure
$\unicode[STIX]{x0394}p_{c}^{+}$
in wall units is also shown.

It is important to note that for this analysis,
$w^{+}$
represents an average cavity width of the random rough surface. It does not say anything about the largest value that is prone to failure first. As the height of the interface increases,
$w^{+}$
increases and more asperities are covered up. Therefore the most realistic numerical simulations would be with an interface below
$S_{q}$
of the roughness where typically the gas fraction
$\unicode[STIX]{x1D719}_{g}<0.6$
. As mentioned earlier, the goal behind our numerical experiment was to investigate the effect of the interface height on the drag reducing properties of SHS. Also it is worth mentioning that although the above analysis gives the upper limit of allowable
$Re_{\unicode[STIX]{x1D70F}}$
for a given
$w^{+}$
, it is known from the literature that adding hierarchical structures to the same size posts can resist destabilisation. Hence, for the same geometry, the maximum allowable
$Re_{\unicode[STIX]{x1D70F}}$
can be larger due to the added multiscale roughness.
Appendix B. Surface statistics
The power spectral density (PSD) of the surface height obtained from the scan is shown in figure 27, where the visible cross-pattern is due to the aliasing effects at the non-periodic boundaries of the unfiltered surface. The original surface statistics are verified with the experimental values provided through private communication. The surface statistics are reported in table 4. Table 5 compares the values of the original surface statistics for the scaled turbulent channel roughness with the step-wise masked representation as used by the solver. Figure 28 illustrates the peaks as they are identified given a threshold.
Statistical parameters of the scanned surface used in the present work.
$N_{x}$
and
$N_{z}$
are the number of points in the streamwise and spanwise directions respectively.
$N_{p}$
denotes the total number of peaks,
$P_{k}$
the peak location and
$h$
the roughness height.

A 2-D slice of a random spanwise location of the original surface roughness (solid blue line) illustrating the identified peaks (solid blue triangles) using the arithmetic mean elevation (dashed red line) as a threshold. The mean distance between peaks is used to obtain the average roughness gap
$L$
.

Comparison of the statistical parameters of the original surface scaled for the turbulent channel flow with the step-wise distribution of the surface used in the present work.













































