Numerical Prediction of Turbulent Drag Reduction with Different Solid Fractions and Distribution Shapes over Superhydrophobic Surfaces

: The exploration of superhydrophobic drag reduction has been and continues to be of signiﬁcant interest to various industries. In the present work, direct numerical simulation (DNS) is utilized to investigate the effect of the parameters on the drag-reducing performance of superhydrophobic surfaces (SHS). Simulations with a friction Reynolds number of 180 were carried out at solid fraction values of φ s = 1 16 , 1 11 , and 1 4 , and three distribution shapes: aligned, staggered, and random. The top wall is the smooth one, and the bottom wall is a superhydrophobic surface (SHS). Drag reduction and Reynolds stress proﬁles are compared for all cases. The turbulent kinetic energy budget, including production, dissipation, and diffusion, is presented with respect to the solid fraction and type of distribution to investigate the drag reduction mechanism. The sizes of the longitudinal vortices and formation of hairpin vortices are investigated through the observation of coherent structures. The simulation of a post model is a useful method to study the drag reduction for different solid fraction values and distribution geometries. Our study demonstrates that the drag reduction could acquire 42% with the solid fraction value φ s = 116 and an aligned distribution shape for post superhydrophobic surface geometry. Our study also showed the relationship of the Reynolds stress component (R 11 , R 22 , and R 33 ) to the drag reduction with the differences in the solid fraction values and distribution geometry. In which, the R 11 component has the most change between an aligned distribution and a random one. The peak value of R 11 tends to shift away from the SHS wall. In addition, the analysis of the TKE budget over the superhydrophobic surface was performed, which can be adopted as a useful resource in turbulence modeling based on RANS methodology. studied by observing coherent structures. Our results showed there was a signiﬁcant reduction in energy dissipation in the ﬂow process by preventing the production of hairpin vortices induced by the instability of superhydrophobic surfaces and the effects of a solid fraction and distribution geometry on the behavior of the coherent structure.


Introduction
Observations on the ability of natural organisms to reduce drag have long attracted interest due to their influence on everyday life and their potential for enormous influence on society. Superhydrophobic (SHS) surfaces are some of the most attractive methods for artificially reducing drag in pipelines, aircraft, and underwater vehicles. It is speculated that superhydrophobic surface drag reduction occurs in nature in a manner similar to how lotus leaves clean themselves. The microstructural rejection of water prevents it from falling into the hollow created by the nanotextured surfaces, allowing water to roll off effortlessly. As shown in Figure 1, the presence of water on the surface of the boundary layer decreases the velocity gradient, thereby reducing the shear stress. A lot of research has been done recently to understand the drag reduction process on superhydrophobic surfaces by modeling different geometries and roughness properties such as ridges and posts [1,2], grooves [3,4], random textures [5][6][7], and slip length models instead of the no-slip boundary condition on the SHS wall [8,9]. In addition, as Philip's theory [9,10] has produced some profound results for laminar wall flows with alternating nonslip/slip surfaces by modeling different geometries and roughness properties such as ridges and posts [1,2], grooves [3,4], random textures [5−7], and slip length models instead of the noslip boundary condition on the SHS wall [8,9]. In addition, as Philip's theory [9,10] has produced some profound results for laminar wall flows with alternating nonslip/slip boundary conditions and has inspired simulations of superhydrophobic surfaces. Lauga and Brenner [11] extended this methodology in pipeline flow simulations. A comparison was made between longitudinal SHS geometry and transverse SHS geometry at the same solid fraction, s θ , and it was found that the slip lengths for longitudinal geometry were twice that for transverse ones. Im and Lee [12] compared the drag reduction mechanism between the turbulent pipe and channel flows with SHS at the wall. They found that the reduction in drag on the pipe flow was caused by an increase of the streamwise slip velocity and a weakening of Reynolds shear stress. Most studies have been performed using turbulent flows on surfaces with a wide range of generally distinctive geometries corresponding to the respective boundary conditions, such as streamwise slips, spanwise slips, and slips in both directions in a turbulent channel flow. These simulations were conducted using direct numerical simulations (DNS) by Min and Kim [13], which is highly noteworthy. As Min and Kim reported, the drag reduction increased proportionally with the slip length when considering the streamwise slip, whereas it decreased with the spanwise slip. Contrary to popular belief, when the operation is mostly slipping, the drag reduction is less than when the operation is mostly streamwise. Similarly, Fukagata et al. [14] and Mori et al. [15] found a correlation between the drag reduction mechanism and Reynolds number, and they also agree with Min and Kim's [13] conclusions. To widen the scope of Min and Kim's study, Martell et al. [1] utilized DNS to investigate the turbulent channel flow with a superhydrophobic surface for streamwise ridges or posts with varying Reynolds numbers ( Re 180, τ = 365, and 590). In their work, SHS has been proven to minimize the drag in turbulent flows. The average slip velocity is roughly 80% of the bulk velocity, with a 50% decline in wall shear stress. They also provide evidence that the SHS created by a post shape has 30% higher drag reduction than a ridge shape and that microfeature spacing has an influence on drag reduction. Seo and Mani [5−7], on the other hand, were curious about the impact of random textures from the SHS texture distribution and a solid fraction ( s φ = 11~25%) on the hydrodynamic performance when subjected to an overlaying turbulent flow Re 180. τ = Within the same Re , S τ φ and feature spacing, random textures can have a slip length that is 30% shorter than the aligned textures. Even though studies on SHS drag reduction have been conducted for several decades, there is still no consensus on the interactions between the dynamics and distribution geometry that contribute to the presence of drag reduction. Therefore, in the present work, different distributions of SHS (aligned, staggered, and random distributions) were considered to study the effect of distribution. To predict the drag reduction and flow mechanism over different SHS patterns, DNS of the turbulent channel flow with SHS was conducted at the condition of Re 180. τ = The slip boundary condition was applied to the air cavity interface. The interface at the post or ridges is assumed to be a perfectly flat surface, and a no-slip wall was placed at the top of the channel, as well as Seo and Mani [5][6][7], on the other hand, were curious about the impact of random textures from the SHS texture distribution and a solid fraction (φ s = 11~25%) on the hydrodynamic performance when subjected to an overlaying turbulent flow Re τ = 180. Within the same Re τ , φ S and feature spacing, random textures can have a slip length that is 30% shorter than the aligned textures. Even though studies on SHS drag reduction have been conducted for several decades, there is still no consensus on the interactions between the dynamics and distribution geometry that contribute to the presence of drag reduction. Therefore, in the present work, different distributions of SHS (aligned, staggered, and random distributions) were considered to study the effect of distribution. To predict the drag reduction and flow mechanism over different SHS patterns, DNS of the turbulent channel flow with SHS was conducted at the condition of Re τ = 180. The slip boundary condition was applied to the air cavity interface. The interface at the post or ridges is assumed to be a perfectly flat surface, and a no-slip wall was placed at the top of the channel, as well as at the top of each post or ridge. The turbulent kinetic energy budgets were calculated for different distribution textures to investigate the drag reduction mechanism. Figure 2 shows the three different SHS distribution geometries: aligned, staggered, and random roughness models. To generate an idealized surface with randomly distributed roughness, each post was placed using a random algorithm (the random algorithm was Energies 2022, 15, 6645 3 of 21 adopted using two boundary conditions on the SHS (slip and no-slip) as the input variables. They are independent of each other, with a total of 2048 positions on the interface of the superhydrophobic surface, created by the constraint of a slip/no-slip ratio of 0.25 (shown in Figure 2c) with the same solid fraction value as the others.
Energies 2022, 15,6645 3 of 23 at the top of each post or ridge. The turbulent kinetic energy budgets were calculated for different distribution textures to investigate the drag reduction mechanism. Figure 2 shows the three different SHS distribution geometries: aligned, staggered, and random roughness models. To generate an idealized surface with randomly distributed roughness, each post was placed using a random algorithm (the random algorithm was adopted using two boundary conditions on the SHS (slip and no-slip) as the input variables. They are independent of each other, with a total of 2048 positions on the interface of the superhydrophobic surface, created by the constraint of a slip/no-slip ratio of 0.25 (shown in Figure 2c) with the same solid fraction value as the others.

Computational Approach and Numerical Methods
In the present work, the channel flow with a smooth wall at the top and SHS at the bottom was considered. The periodic boundary condition was applied in both the streamwise and spanwise directions, shown in Figures 3 and 4. The dimensionless length of the channel was Lx/H = 6, and the width was LZ/H = 3, where H is the half-height of the channel. This is similar to the former channel flow simulation with 2π and π in the streamwise and spanwise directions, respectively [1,16]. While the simulation does not require dimensions, the simulation results were compared with experimental data, where the working fluid was water at 20 °C, and a half-height of the channel was about 0.15 mm if the posts or ridges were assumed to be 30 µm [17,18]. The size of the computational domain, grid information, and resolution of space and time are presented in Table 1, which are based on the DNS study of fully developed channel flow by Moser et al. [16] and DNS using OpenFoam by Theobald et al. [19]. The simulation cases are shown in Table 2.

Computational Approach and Numerical Methods
In the present work, the channel flow with a smooth wall at the top and SHS at the bottom was considered. The periodic boundary condition was applied in both the streamwise and spanwise directions, shown in Figures 3 and 4. The dimensionless length of the channel was Lx/H = 6, and the width was LZ/H = 3, where H is the half-height of the channel. This is similar to the former channel flow simulation with 2π and π in the streamwise and spanwise directions, respectively [1,16]. While the simulation does not require dimensions, the simulation results were compared with experimental data, where the working fluid was water at 20 • C, and a half-height of the channel was about 0.15 mm if the posts or ridges were assumed to be 30 µm [17,18]. The size of the computational domain, grid information, and resolution of space and time are presented in Table 1, which are based on the DNS study of fully developed channel flow by Moser et al. [16] and DNS using OpenFoam by Theobald et al. [19]. The simulation cases are shown in Table 2.
at the top of each post or ridge. The turbulent kinetic energy budgets were calculated for different distribution textures to investigate the drag reduction mechanism. Figure 2 shows the three different SHS distribution geometries: aligned, staggered, and random roughness models. To generate an idealized surface with randomly distributed roughness, each post was placed using a random algorithm (the random algorithm was adopted using two boundary conditions on the SHS (slip and no-slip) as the input variables. They are independent of each other, with a total of 2048 positions on the interface of the superhydrophobic surface, created by the constraint of a slip/no-slip ratio of 0.25 (shown in Figure 2c) with the same solid fraction value as the others.

Computational Approach and Numerical Methods
In the present work, the channel flow with a smooth wall at the top and SHS at the bottom was considered. The periodic boundary condition was applied in both the streamwise and spanwise directions, shown in Figures 3 and 4. The dimensionless length of the channel was Lx/H = 6, and the width was LZ/H = 3, where H is the half-height of the channel. This is similar to the former channel flow simulation with 2π and π in the streamwise and spanwise directions, respectively [1,16]. While the simulation does not require dimensions, the simulation results were compared with experimental data, where the working fluid was water at 20 °C, and a half-height of the channel was about 0.15 mm if the posts or ridges were assumed to be 30 µm [17,18]. The size of the computational domain, grid information, and resolution of space and time are presented in Table 1, which are based on the DNS study of fully developed channel flow by Moser et al. [16] and DNS using OpenFoam by Theobald et al. [19]. The simulation cases are shown in Table 2.     The open-source code OpenFoam was adopted for the present DNS. The coordinates and flow variables of the simulation were normalized by the channel half-height H and the friction velocity, u , where ω τ is the wall shear, and ρ is the density of the fluid.
The governing equations of the continuity and momentum equations are presented in a nondimensional form below:  The open-source code OpenFoam was adopted for the present DNS. The coordinates and flow variables of the simulation were normalized by the channel half-height H and the friction velocity, u τ = τ ω ρ , where τ ω is the wall shear, and ρ is the density of the fluid. The governing equations of the continuity and momentum equations are presented in a nondimensional form below: where u i , t, and p are the velocity vector, time, and pressure, respectively. The superscript (+) indicates values normalized by the variables mentioned above. Finally, the Reynold number Re τ = u τ H ν appears in the diffusion term. The PISO algorithm, which was proposed by Issa [20], was adopted for an incompressible flow [21,22] to overcome unsteady simulations and become integrated with C ≈ 1. The implicit time advancement employed the second-order Crank-Nicolson scheme for the convection terms. The convection term is discretized with a 2nd-order Gauss linear corrected scheme. A constant pressure gradient was applied as the driving force [23,24].
For the initial flow fields, the near-wall parallel streaks of slow and faster-moving fluid are produced (U + 0 ) as expressed in the equation below [25]: and the spanwise velocity component normal to the streaks is: Energies 2022, 15, 6645

of 21
Time is normalized by the bulk timescale t * , which is equivalent to a flow through time of the domain by the equation: where a + and b + are the constant values for which two relations are chosen to produce a sparse streak spacing (z + ≈200). U m is the mean bulk velocity, and L x is the size of the domain in the streamwise direction. Turbulent kinetic energy budgets such as production, dissipation, and diffusion are expressed as follows [26]: Production: Turbulent diffusion: Vel. P. -grad. corr.: Molecular diffusion: Dissipation: Table 3 shows the mean flow variables of the channel flow in the case of a smooth wall. Based on this condition, the width of the ridge was set to w + = wδ/ν ≈ 0.1875, 0.14062, and 0.09375 H, and the corresponding ratios of the gap/width were 1, 1.6, and 3.0. They are equivalent to sizes of w ≈ 28.125, 21.093 and 14.0625 µm, respectively, in realistic experiments [17,18], as shown in Table 2. In the present work, the solid fraction, φ s , which is calculated by the ratio of the area of the SHS wall on the interface normalized by the projected area of the bottom wall 2 , which is shown in Figure 5. In the experiments [17,18], the value of the solid fraction was about 20% (when the post sizes are assumed to be 30 µm) the range of the solid fraction is 10-20% [27,28] and 6.25% at higher Reynolds numbers). In this study, we considered θ s = 1 16 , 1 11 , and 1 4 , which are presented in Table 4.
, which is shown in Figure 5. In the experiments [17,18], the value of the solid fraction was about 20% (when the post sizes are assumed to be 30 µm, the range of the solid fraction is 10-20% [27,28] and 6.25% at higher Reynolds numbers). In this study, we considered , , and , which are presented in Table 4.  Table 4. Reynolds numbers, post pattern, solid fraction, and case notation in the post type, with name cases from C4 to C11.

Results and Discussion
To validate the present solver, OpenFoam, a fully developed turbulent channel flow with a smooth wall was simulated with the C0 case at the condition of Re 180 τ = . The mean flow variables are shown in Table 3. Figure 6 shows a comparison of Reynolds stress profiles and near-wall velocity profiles with reference DNS results by Moser et al. [16].
Both the Reynolds stress ( ) ' ' i j u u and near-wall velocity profiles (u + ) agree well with the reference DNS results. Turbulent kinetic energy budgets were calculated using Equations (6)-(10) for comparison with the reference data [16]. The production, dissipation, turbulent diffusion, pressure diffusion, and molecular diffusion are plotted in the wall normal direction, and the residuals of all terms are plotted as a black line in Figure 7. Even though there is a small discrepancy with the reference data [16] near the peaks in molecular  Table 4. Reynolds numbers, post pattern, solid fraction, and case notation in the post type, with name cases from C4 to C11.

Results and Discussion
To validate the present solver, OpenFoam, a fully developed turbulent channel flow with a smooth wall was simulated with the C0 case at the condition of Re τ = 180. The mean flow variables are shown in Table 3. Figure 6 shows a comparison of Reynolds stress profiles and near-wall velocity profiles with reference DNS results by Moser et al. [16]. Both the Reynolds stress u i u j and near-wall velocity profiles (u + ) agree well with the reference DNS results. Turbulent kinetic energy budgets were calculated using Equations (6)-(10) for comparison with the reference data [16]. The production, dissipation, turbulent diffusion, pressure diffusion, and molecular diffusion are plotted in the wall normal direction, and the residuals of all terms are plotted as a black line in Figure 7. Even though there is a small discrepancy with the reference data [16] near the peaks in molecular diffusion and turbulent diffusion, all terms are consistent with the reference data. Thus, we confirmed that the present solver can resolve turbulent flow characteristics in the channel geometry.
As justified by the validation discussed above, we employed the present solver to conduct the direct numerical simulation (DNS) of the superhydrophobic channel flow, with the remaining cases indicated in Table 2 as C1-C4. The drag reduction (DR) was calculated as the difference in the skin friction between the smooth wall case and each SHS case as expressed in the equation below: where C f denotes the friction coefficient of the superhydrophobic surface, is the dimensionless streamwise mean velocity, and C D f represents the friction coefficient at the smooth wall, , and Re b Rebis the bulk Reynolds number [20]. diffusion and turbulent diffusion, all terms are consistent with the reference data. Thus, we confirmed that the present solver can resolve turbulent flow characteristics in the channel geometry.   [16] in terms of the production, dissipation, turbulent diffusion, molecular diffusion, and pressure diffusion. Energies 2022, 15, 6645 7 of diffusion and turbulent diffusion, all terms are consistent with the reference data. Thu we confirmed that the present solver can resolve turbulent flow characteristics in the chan nel geometry.    [16] in terms of the production, dissipation, turbulent diffusion, molecular diffusion, and pressure diffusion.

Drag Reduction and Mean Velocity Profile
In this section, the drag reduction and mean velocity profiles of each case were investigated. The instantaneous flow fields were averaged during the nondimensional time t(U m /H)~200 in both streamwise and spanwise directions to obtain the mean flow fields. Based on the smooth wall case (C0), the drag reduction of the ridge type (C1~C3) and post type (C4) cases were calculated using Equation (11), and the results are presented in Table 5.
In the ridge type of SHS, a wider gap results in a larger drag reduction ratio of about 1.5. At the same width-to-gap ratio of the cavity, the post type (C4) case showed 5.5 times more drag reduction than the ridge one (C3). As shown in Figure 8, the mean streamwise velocity profile was compared with the results of Martell et al. [1], which were obtained by the DNS of SHS at the same geometry condition (width/gap = 1) and flow condition (Re τ = 180) as the C1 case. The top part (y/H = 1.0) is at the smooth wall, and the bottom part (y/H = −1.0) is for the ridge type of SHS. The streamwise velocity was normalized using the friction velocity, u τ , which agrees well with Martell et al.'s results [1]. The small discrepancy in the SHS region is not caused by the same boundary conditions at the interface of SHS, but the increase in the mean velocity in this region compared to the smooth wall region is identical to the reference one [1]. Therefore, it is possible to assume that the SHS region has a nonslip boundary condition, and the interface between water and air is a slip boundary condition. Based on a comparison with Martell et al.'s results [1], we confirmed that the present simulation setup for SHS is appropriate, and additional simulations with different solid fractions and post distributions would show reliable results. The mean streamwise velocity profiles for three different solid fractions ( , , and ) in the post types (C4, C5, and C6) are presented with that of the smooth wall case (C0) in Figure 9. The three cases have an aligned distribution of post-type SHS, as mentioned in Table 4. The bulk velocity in the SHS region increases as the solid fraction The mean streamwise velocity profiles for three different solid fractions (φ s = 1 16 , 1 11 , and 1 4 ) in the post types (C4, C5, and C6) are presented with that of the smooth wall case (C0) in Figure 9. The three cases have an aligned distribution of post-type SHS, as mentioned in Table 4. The bulk velocity in the SHS region increases as the solid fraction decreases (which means that the portion of the interface between the liquid and gas increases). Additionally, this increase in the SHS region affects the velocity profile in the smooth wall region, which also increases. The bulk velocity at the condition of φ s = 1 16 is 40% larger than that for a smooth wall. As the solid fraction increases to φ s = 1 11 and φ s = 1 4 , the increase of the bulk velocity decreased to 32% and 17%, respectively. When the decrease in the solid fraction increased the drag reduction, the total mass flux increased with a decrease in the solid fraction, and then, the drag reduction induced an increase in the mass flux in the channel with SHS. As an alternative to calculate the drag reduction, Martell et al. [1] and others employed shear stress reduction to examine the link between parameters and drag reduction. With a constant pressure gradient in the channel, the Reynolds stresses will be indicated by the symbol ij R . The average wall friction is defined as  Figure 10 shows the slip velocity slip U as a percentage of the bulk velocity bulk U for three different distributions of post types: aligned, staggered, and random distribution, with respect to the solid fraction. The solid portion of the superhydrophobic surface was a significant factor influencing the slip velocity of the surface. As expected, the slip velocity normalized by the bulk velocity is inversely related to the solid fraction. Both the As an alternative to calculate the drag reduction, Martell et al. [1] and others employed shear stress reduction to examine the link between parameters and drag reduction. With a constant pressure gradient in the channel, the Reynolds stresses will be indicated by the symbol R ij The average wall friction is defined as u τ = (τ ω ) 2 , in which τ ω = τ T w + τ B w /2, where τ T w denotes the top wall shear stress, and τ B w τwBdenotes the bottom wall shear stress. Since the top wall of our channel is flat, we can compute τ B w by averaging µ ∂u ∂y y=H and τ B w by defining 2τ ω − τ B w . When u B τ = τ B w 1/2 is used to determine the bottom wall friction velocity, it is referred to as the bottom wall friction velocity. Figure 10 shows the slip velocity U slip as a percentage of the bulk velocity U bulk for three different distributions of post types: aligned, staggered, and random distribution, with respect to the solid fraction. The solid portion of the superhydrophobic surface was a Energies 2022, 15, 6645 10 of 21 significant factor influencing the slip velocity of the surface. As expected, the slip velocity normalized by the bulk velocity is inversely related to the solid fraction. Both the aligned and staggered distributions have a negligible impact on the slip velocity for all solid fraction cases in the present work. However, the random distribution has a significant impact on the slip velocity percentage, especially at the low solid fraction. This effect diminishes at the high solid fraction, φ s = 1 4 . The differences in the slip velocity between the aligned and random distributions are 16%, 4%, and 1% at φ s = 1 16 , 1 11 , and 1 4 , respectively. Therefore, we concluded that the solid fraction and distribution geometry are strong factors in the observed dimensionless slip velocity of the SHS. Since the friction Reynolds number, Re τ , and pressure gradient are kept constant throughout the present DNS, the decreased drag on superhydrophobic surfaces would result in a larger fluid mass flow rate. The increase of the mass flow rate induces a larger change in the bulk velocity, which would result in an increased drag on the top wall of the channel. This description can be confirmed in Figure 11, which shows the mean profile for the velocity normalized by the bottom wall friction velocity B u τ in the cases investigated. A nearly linear decrease in the bulk velocity is observed for all investigated cases with the same distribution geometry but different solid fraction values. This trend shows that the relationship between the solid fraction and bulk velocity is a linear function on the SHS wall. Since the friction Reynolds number, Re τ ,Reτ and pressure gradient are kept constant throughout the present DNS, the decreased drag on superhydrophobic surfaces would result in a larger fluid mass flow rate. The increase of the mass flow rate induces a larger change in the bulk velocity, which would result in an increased drag on the top wall of the channel. This description can be confirmed in Figure 11, which shows the mean profile for the velocity normalized by the bottom wall friction velocity u B τ in the cases investigated. A nearly linear decrease in the bulk velocity is observed for all investigated cases with the same distribution geometry but different solid fraction values. This trend shows that the relationship between the solid fraction and bulk velocity is a linear function on the SHS wall.    Figure 12 presents the slip velocity normalized by the bottom wall friction velocity versus the solid fraction value, and Figure 13 shows the percentage of shear stress reduction, ∆τ ω (=τ ω − τ B w ), relative to the average shear stress in the channel, τ ω . These figures show that the drag decreases as the solid fraction increases, and the random distribution of posts showed different drag reductions from the aligned and staggered ones. However, there was only a small difference in drag reduction between the aligned and staggered distributions. In the low solid fraction, φ s = 1 16 , the slip velocity and shear stress difference decreased when compared with the aligned or staggered one. However, in the other two solid fractions, φ s = 1 11 and φ s = 1 4 , the drag reduction increased more than that of the aligned or staggered distributions. This analysis provides further evidence that the featured solid fraction value and distribution geometry play a crucial role in the surface performance.

Reynolds Stresses
In this section, Reynolds stresses, including normal stress (R11, R22, and R33) and the shear one (R12), are plotted in Figures 14 and 15

Reynolds Stresses
In this section, Reynolds stresses, including normal stress (R11, R22, and R33) and the shear one (R12), are plotted in Figures 14 and 15

Reynolds Stresses
In this section, Reynolds stresses, including normal stress (R 11 , R 22 , and R 33 ) and the shear one (R 12 ), are plotted in Figures 14 and 15, which are the cases of φ s = 1 16 and 1 4 , respectively, to investigate the turbulence levels and structure of the near-wall turbulence of SHS. The Reynolds stresses of three different distribution patterns are compared in each figure. Since the differences in drag reduction between the aligned distribution and random one is the largest at φ s = 1 16 and the smallest at φ s = 1 4 , Reynolds normal and shear stress are compared in both solid fraction cases to investigate the different behaviors. and 25% for 1 4 s φ = compared with the peak at the smooth wall. Additionally, the peak of R11 Reynolds stress at the SHS shifted closer to the slip wall. The slip effect can cause a substantial change in the near-wall turbulence. The Reynold stress components R22 and R33 are shown in Figures 14c,d and 15c,d. The peaks of two components, R22 and R33, near the SHS also decreased, similar to the R11 component, but the peak shifted away from the SHS wall, contrary to the R11 component. The peak in the Reynolds shear stress near the SHS wall decreased by two or three times compared to the one near the smooth wall, as shown in Figures 14b and 15b. The biggest difference between an aligned distribution and a random one is shown in the R11 Reynolds stress component. The peak value of this was reduced in both solid fraction cases. Therefore, it is clear that there is a The normal Reynolds stress (R 11 , R 22 , and R 33 ) near the SHS is suppressed compared to that near the smooth wall surface in both cases, and the peak of each component decreased. The degree of decrease is greater in the φ s = 1 16 case than in the φ s = 1 4 case. We further observed that the peak of R 11 Reynolds stress at SHS decreased up to 50% for φ s = 1 16 and 25% for φ s = 1 4 compared with the peak at the smooth wall. Additionally, the peak of R 11 Reynolds stress at the SHS shifted closer to the slip wall. The slip effect can cause a substantial change in the near-wall turbulence. The Reynold stress components R 22 and R 33 are shown in Figure 14c,d and Figure 15c,d. The peaks of two components, R 22 and R 33 , near the SHS also decreased, similar to the R 11 component, but the peak shifted away from the SHS wall, contrary to the R 11 component. The peak in the Reynolds shear stress near the SHS wall decreased by two or three times compared to the one near the smooth wall, as shown in Figures 14b and 15b. The biggest difference between an aligned distribution and a random one is shown in the R 11 Reynolds stress component. The peak value of this was reduced in both solid fraction cases. Therefore, it is clear that there is a relationship between the solid fraction value and distribution geometry for the velocity fluctuations at the SHS surface. relationship between the solid fraction value and distribution geometry for the velocity fluctuations at the SHS surface.

Turbulent Kinetic Energy(TKE) Budget and Flow Structures
For an in-depth analysis of the turbulent characteristics and an energy balance over the SHS, the turbulent kinetic energy budgets in the turbulent channel flow with SHS were calculated using Equations (6)- (10) and were compared with those in the smooth wall cases. Rastegari and Akhavan [29] analyzed turbulent intensities and the mean convective shear stress over superhydrophobic microgrooves and riblet through DNS based on the Lattice Boltzmann method. However, to the best of the authors' knowledge, there are no detailed studies of turbulent kinetic energy budgets such as the production, dissipation, and diffusion in the ridge and post geometries of the superhydrophobic surfaces. Figure  16 shows a comparison of the turbulent kinetic energy budgets between SHS at and with aligned distributions and smooth wall cases.

Turbulent Kinetic Energy(TKE) Budget and Flow Structures
For an in-depth analysis of the turbulent characteristics and an energy balance over the SHS, the turbulent kinetic energy budgets in the turbulent channel flow with SHS were calculated using Equations (6)- (10) and were compared with those in the smooth wall cases. Rastegari and Akhavan [29] analyzed turbulent intensities and the mean convective shear stress over superhydrophobic microgrooves and riblet through DNS based on the Lattice Boltzmann method. However, to the best of the authors' knowledge, there are no detailed studies of turbulent kinetic energy budgets such as the production, dissipation, and diffusion in the ridge and post geometries of the superhydrophobic surfaces. Figure 16 shows a comparison of the turbulent kinetic energy budgets between SHS at φ s = 1 16   The dissipation of TKE is shown in Figure 16a along the wall normal direction. Whether it is a turbulent drag-reducer flow with SHS or simple turbulent flow a with smooth wall, dissipation, ε ij , is always negative, which means that the turbulent energy is dissipated and ε ij has a peak at the wall. This peak in SHS is higher than that at a smooth wall, which is caused by strong energy interactions near the interface between the liquid and gas, and then, there is a rapid decrease to the converged value. Figure 16b shows the behavior in TKE production near SHS and smooth wall. In the smooth wall case, the production increases from zero up to the peak and then decreases. However, in the SHS case, there is some production at the interface, which is not zero, then a gradual decrease at φ s = 1 16 and a peak at a position (y +~1 5) similar to the smooth wall case at φ s = 1 4 . However, the peak values at SHS are only 71% and 48% of that in the smooth wall case. The dramatic drop in SHS production can be attributed to a decrease in Reynolds shear stress within the turbulent boundary layers. Figure 16c indicates the presence of an influencing zone of viscous diffusion focused on the near-wall region in both the SHS and smooth wall cases. In the transition layer, the viscous diffusion achieves a minimum value, which means the energy dissipation induced by viscous diffusion is the maximum. The peak value of viscous diffusion in SHS obviously decreases compared with that in the smooth wall case. It can be observed that the turbulent diffusion in SHS is weaker in the near-wall region than that in the smooth wall case and tends to be stable after the negative peak in Figure 16d. Thus, the turbulent fluctuation kinetic energy decreased gradually to a stable value. The variation in pressure diffusion in Figure 16e shows similar behavior with that of turbulent diffusion.
Based on the analysis of TKE budgets in SHS, we confirmed that the addition of SHS to the turbulent flow obviously changes the production, transport, and dissipation of TKE, especially in the near-wall region. Additionally, the reduction of the turbulent drag is closely related to the dramatic decrease in the production of TKE. For visualization of the coherent structures near the wall in the channel, the Q-criterion (known as the second invariant of the velocity gradient tensor, as established by Hunt et al. [30] and Wu and Moin [31]) is plotted in Figures 17-19. The calculated vortex core isosurface was plotted with the contour of the streamwise velocity with a Q-criterion value of 200. In all cases with SHS, the coherent structures, such as the hairpin or crutch vortices [32][33][34][35], are considerably less dense than those in the smooth wall. When fluids flow over surfaces with solid fraction values varying from 1 16 to 1 4 , the effects of the inhibition become stronger. The same phenomenon occurs when we change the type of distribution geometry from random to aligned. In the staggered distribution, the number of the streamwise vortices decreases compared to the cases of other distributions. On the contrary, the random distributions delayed the reduction of the streamwise vortices. As a result, the energy dissipation in the flow process was significantly reduced.       Figure 20 shows the velocity contour near y + = 15 for the SHS wall for the φ s = 1 16 cases, showing a high drag reduction for the investigated cases. The velocity levels in Figure 20a,b are approximately equal to and greater than the velocity values in Figure 20c at the same position. The effects of the distribution geometry to slip velocity at the SHS wall are apparent.

Conclusions
DNS of the turbulent channel flow at Re 180 τ = were considered, where the bottom wall was a superhydrophobic surface with different solid fraction values and distribution geometries. The effects of the solid component value and distribution geometry on the drag reduction was explored in detail through an analysis of the Reynolds stress profile, turbulent kinetic energy budgets, and coherent structures.
The simulation results showed a drag reduction effect of up to 42% for a solid fraction value 1 16 and an aligned distribution shape. The drag reduction decreased as the scale value increased. We found that there is a variation in the effect of drag reduction at low solid fraction values. The random distribution has less impact on drag reduction than the aligned and staggered distributions. When the solid fraction increased, the effect of drag reduction gradually decreased. In contrast, there was no big difference in drag reduction between the aligned and staggered distributions. According to the simulation results for all geometries, the mean velocity profile was normalized by two quantities: the wall friction velocity u τ and bottom wall friction velocity B u τ . The velocity profiles illustrated a reduction in the mass flux due to a superhydrophobic surface. For the case with a solid fraction value 1 16 s φ = and aligned distribution (which corresponds to case 4), the bulk velocity bulk U was up to 40% higher than the channel with a smooth wall, and the percentage slip velocity slip U normalized by the

Conclusions
DNS of the turbulent channel flow at Re τ = 180 were considered, where the bottom wall was a superhydrophobic surface with different solid fraction values and distribution geometries. The effects of the solid component value and distribution geometry on the drag reduction was explored in detail through an analysis of the Reynolds stress profile, turbulent kinetic energy budgets, and coherent structures.
The simulation results showed a drag reduction effect of up to 42% for a solid fraction value φ s = 1 16 and an aligned distribution shape. The drag reduction decreased as the scale value increased. We found that there is a variation in the effect of drag reduction at low solid fraction values. The random distribution has less impact on drag reduction than the aligned and staggered distributions. When the solid fraction increased, the effect of drag reduction gradually decreased. In contrast, there was no big difference in drag reduction between the aligned and staggered distributions.
According to the simulation results for all geometries, the mean velocity profile was normalized by two quantities: the wall friction velocity u τ and bottom wall friction velocity u B τ . The velocity profiles illustrated a reduction in the mass flux due to a superhydrophobic surface. For the case with a solid fraction value φ s = 1 16 and aligned distribution (which corresponds to case 4), the bulk velocity U bulk was up to 40% higher than the channel with a smooth wall, and the percentage slip velocity U slip normalized by the mean bulk velocity U bulk showed values up to 75%. The latter is a significant parameter, as it is directly related to the drag reduction. Since the pressure gradient in this study is constant, a reduced drag on the SHS will result in an increased drag on the upper surface (smooth wall) and the same total drag in the channel. These results can help us understand the importance of solid fraction and distribution geometry to drag reduction on channels with SHS and correlate these values.
Reynolds stress is a critical quantity to understand turbulent behavior, and it is closely related to the turbulence levels and near-wall turbulence structure. The peak of Reynolds stress tends to shift towards the SHS surface, while it conversely pushes away from the smooth wall. The parameters' streamwise Reynolds stress, R 11 , R11and Reynolds shear stress, R 22 , show a collapse in the profile data when normalized by u B τ and have near-wall turbulence, which is shown in Figures 14 and 15. This led to a deeper investigation into the near-wall turbulence structures.
The kinetic energy budgets were compared at different types of distribution geometries on the superhydrophobic surface. Based on the change of the turbulent component in the wall normal and spanwise directions, the SHS surface changed the production, transport, and dissipation of TKE, and there was a relationship between the solid fraction parameter and change in the total kinetic energy budget. As the solid fraction of the SHS surface increased, the change in the total kinetic energy budget decreased as compared with the smooth surface, especially when approaching the SHS surface. Therefore, the TKE budget through the SHS surface could not reach an equilibrium. The presence of superhydrophobic surfaces near the wall minimized the internal velocity gradients and stress in the near-wall region. As a result, the flow resistance caused by viscosity was significantly reduced. Significant drag reduction occurred, particularly when the solid fraction value was decreased and the distribution geometry was changed. The analysis of the TKE budget over SHS in the present work can be adopted as a useful resource in turbulence modeling based on RANS methodology.
The formation of hairpin vortices and the sizes of the longitudinal vortices were studied by observing coherent structures. Our results showed there was a significant reduction in energy dissipation in the flow process by preventing the production of hairpin vortices induced by the instability of superhydrophobic surfaces and the effects of a solid fraction and distribution geometry on the behavior of the coherent structure.