Influence of small wavy roughness on flatplate boundary layer natural transition

One of the effective ways to reduce viscous drag around an airfoil is by delaying the boundary layer transition. In this study, we analyzed the influence of a small wavy roughness on a two-dimensional, natural boundary layer transition, using direct numerical simulation that resolved each small roughness. A parametric study was conducted on the wavy roughness wavelength. Our results show that in some cases the transition delays whose characteristics depend on the roughness wavelength. In a detailed analysis, we found that the wavy roughness firstly affects the process of primary vortex growth, Tollmien–Schlichting (TS) instability. In addition, we found that the secondary vortex pairing also depended on it. In the most transition-delayed cases, the roughness wavelength was different far from the TS instability one, and the vortex pairing occurred firstly in upstream however not much in downstream, keeping the vortex size is kept small.


Introduction
In recent years, the demand for air transportation has continued to increase. According to a survey by the Ministry of Land, Infrastructure and Transport, the world's number of flights worldwide has increased by approximately 80% between 1990and 2014(The Ministry of Land, Infrastructure and Transport, 2007. Consequently, CO2 and NOx emissions have increased by approximately 80% and 200%, respectively in the same period. This trend is expected to continue. Therefore, it is necessary to develop an environmentally friendly high-performance aircraft. Friction drag is one of the highest contributors to the total drag acting on an aircraft. Reducing the friction drag acting on an aircraft significantly improves its fuel efficiency and environmental friendliness. There have been several ways to reduce friction drag; controlling drag-reduction in the full turbulent state, turbulence control, otherwise keeping the flow laminar or quasi-stable state, and delaying the transition to the turbulent state, laminarization. One of the turbulence control strategies is the use of a riblet. A riblet is a grooved surface regularly fabricated in the streamwise direction (Caram et al., 1991). The study by Walsh (1992) showed that the use of the riblet produces a friction drag reduction effect. The study suggested that the use of the riblet to reduce friction drag can be achieved by changing the shape of the vortex near the wall surface. Currently, JAXA is developing a riblet for aircraft (Japan Aerospace Exploration Agency, 2018). However, the practical use of riblets in an aircraft has not yet been implemented because of problems associated with their maintainability and durability. In addition, there is turbulence control that can be used to decrease turbulent energy production and friction drag by changing the surface shape and the fluid itself (Yakeno et al., 2014 etc.). On the other hand, laminarization is one of the transition control technologies to reduce viscous friction drag significantly. Many friction drag reduction methods have been proposed in existing studies. Boundary layer transition control for laminarization is indispensable in the development of a high-performance aircraft and has been used in many aircrafts for years. Recently, attention has been focused on new laminarization technique, with improvements in surface processing technology. This study focused on the laminarization technology of the boundary layer transition.
A microroughness structure is one of a potential way to delay the transition and also to reduce friction drag. The microroughness surface structure is even smaller than the riblet. Tani (1988) reported that the use of a microroughness surface structure reduced friction drag in the turbulent state. This effect was confirmed in a wind tunnel experiment using a flat plate with randomly distributed microroughness surface structures (Oguri et al., 1996). Surprisingly, the roughness height of this microroughness surface structure is a scale previously considered aerodynamically smooth. The roughness Reynolds number & (= ℎ % / ) is 0.23-1.5. The microroughness shape is completely randomly distributed and it is clear from the above that the friction drag reduction mechanism of the microroughness surface structure is different from that of the riblet. Besides, a boundary layer transition delayed owing to the presence of a fibrous micro-surface structure known as the fiber surface (Oguri et al., 1998). The fiber surface delayed the bypass transition. Oguri et al. (1998) predicted the suppression of a streak structure that occurs at the initial stage of the boundary layer transition. Furthermore, Kikuchi et al. (2004) performed several predictions on the transition delay mechanism of the fiber surface.
Previous studies have demonstrated that the use of a distributed microroughness surface and a fiber-structure surface resulted in delaying transition and friction drag reduction. However, detailed flow-structure modification and the control mechanisms have not been clarified. A good understanding of these mechanisms will improve the development of a surface structure to achieve friction drag reduction. Additionally, there have been no references for the control performance on a natural transition without any artificial disturbance in upstream, while Ma'mum et al. (2014) investigated the effects of a surface corrugation on the artificial Tollmien-Schlichting (T-S) wave experimentally. Therefore, in this study, we analyzed the effect of a wavy rough surface smaller than the boundary layer displacement thickness of different wavelengths on a natural transition using direct numerical simulation (DNS). As a first step, we investigated effects on a two-dimensional transition due to the computational constraints.

Governing Equation
The governing equation in this study is the compressible Navier-Stokes equation given below.
Here, the variables Q, E, F, G, ( , ( , ( are defined as: The Navier-Stokes equation is made dimensionless by introducing the density $ * , the speed of sound $ * , and the dynamic viscosity $ * of the free-stream, as well as the reference length . * . In Eq. (1), is the Reynolds number defined with the free-stream velocity $ * ; In this study, we consider the input laminar boundary layer thickness !" * as the reference length . * and define the following Reynolds number # by setting it as 1000, with the free-stream Mach number $ of 0.2. All the variables without any subscription in the present paper are non-dimensionalized with using !" * , $ * , $ * and $ * . In this Mach number condition, we confirm that both the temperature and the density fluctuation vary by less than approximately 1% over the roughness, indicating a quasi-incompressible flow.
The dynamic viscosity ( ) is given by the Sutherland's law for air, where = 198.6 and the Prandtl number Pt and the specific heat ratio are 0.72 and 1.4, respectively.

Numerical Scheme
In this study, the equations are solved in the generalized curvilinear coordinate system, where the spatial derivatives for the convective terms, viscous terms, metrics, and Jacobian are evaluated using the following sixth-order compact difference scheme (Lele, 1992). The first spatial derivatives are obtained as given in the following formula: Tameike, Yakeno and Obayashi, Journal of Fluid Science and Technology, Vol.16, No.1 (2021) is the number of grids, ! is a scalar quantity， ! 6 is the first spatial derivative of ! , and h is the width between the grid points. The relationship between the coefficients , , and , is derived from various order of Taylor series expansions. In the case of the sixth-order compact scheme, the coefficients are evaluated as follows: Substituting these coefficients into Eq. (6) gives A filtering procedure for the compact difference scheme is needed to suppress numerical instabilities (Gaitonde et al., 2000). In this study, the following tenth-order filtering scheme was used.
where is the solution vector, and ‚ is the filtered quantity. In Eq. (11), the tenth-order (N=6) filter is obtained with The term > is a free parameter satisfying the inequality -0.5 < > ≤ 0.5; > is set as 0.495 with the 6 th order compact scheme in the boundary-layer transition study by Kawai et al. (2008). A high accuracy time integration method is required to analyze the unsteady flow phenomena. In this study, we used the third-order three-step total variation diminishing (TVD) Runge-Kutta time integration method (Shu et al., 1988).

Boundary Conditions
In this study, we performed a two-dimensional simulation of the flow over a flat plate by resolving each small wallroughness. The boundary conditions were defined at four points: streamwise inlet boundary, streamwise outflow boundary, wall-normal outflow boundary, and wall surface. A laminar Blasius boundary layer profile was used at the streamwise inlet boundary. The pressure p and the density ρ at the streamwise inlet boundary were fixed as $ and $ , respectively by setting the total energy e using the following equation, where the subscript Blasius indicates the value of the Blasius boundary layer profile. Figure 3 shows the cases verified in this study. In each case, the wavy surface roughness height h was set to be 0.5, which is smaller than the boundary layer displacement thickness. And the wavelength was varied from 3.12 to 10.0. These wavelengths are under a size of the TS instability roll vortex at the present boundary layer. The effect of the wavy rough surface wavelength on boundary-layer transition was verified by comparing Cases 1 to 5. Cases 2 and 6 have a 4 At the streamwise outflow boundary, the velocities and the pressure were set as their gradients are zero. At the wallnormal outflow boundary, the slip boundary condition was set as the streamwise direction velocity and the pressure were fixed as $ and $ , respectively, and the wall velocity w was set as the gradient is zero. We applied the non-slip boundary and adiabatic condition to the wall surface.

Computational Setup 3.1 Computational Condition
The computational domain is shown in Fig. 1. The test section size is I (= * I / * !" ) = 1200 and J (= * J / * !" ) = 100 in the streamwise and the wall-normal directions, non-dimensional with the input boundary layer thickness * !" . The surface is smooth from the inlet boundary = 0 to = 120, and it changes to wavy roughness from = 120 to = 1200. A sponge layer is placed near the outflow boundary with size is I = 3000 and J = 1000.
In this study, we used a structured grid that resolves each small roughness. Figure 2 shows the shape of the wavy rough surface and the computational grid. The wavy rough surface was a sine curve. The boundary from the smooth surface to the wavy rough surface is connected smoothly using the following sigmoid function. similar wavelength . However, the entrance shape in the two cases is opposite in phase, and the starting position of the first roughness slope is different. In this study, the effect of the entrance shape of the wavy rough surface on boundarylayer transition was verified by comparing Cases 2 and 6.

Computational Grid
The present computational grid resolution was selected based on a grid-convergence study. Table 1 presents the medium and fine grid resolutions in the test section, grid number, and the minimum grid width in the streamwise and the wall-normal directions, respectively. The total numbers of grids in the test section were 1,645,880 and 4,937,640. We compared the root mean square values of wall-normal velocity using two different resolution grids for Case 2 (Fig. 4), i.e., medium and fine. In the fine grid case, the increase of the wall-normal velocity RMS was slightly behind compared to the medium grid one. However, for the comparison of control performance among tested cases, its difference was not very considerable and we decided to use the medium grid for this study.

Inflow boundary layer
In this study, the wavy rough surface starts at = 120. Figure 5 shows the streamwise velocity u profiles at = 100 in the Cases 1-6 and it is confirmed the inflow profiles are similar. Table 2 presents the 99% thickness of the boundary layer, displacement and momentum thicknesses. In Cases 1-6, the wavy rough surface height ℎ (= 0.5) is 23.7% of the inflow boundary layer thickness KK and 78.6% of * .

Effect of wavy rough entrance shape
As mentioned above, the wavy rough surface starts from = 120. That is, the shape of the flat plate changes from a smooth surface to a wavy rough surface. The effect of the wavy rough entrance shape on boundary layer transition has not been clarified in previous studies. In this study, the effect was verified by comparing two cases, namely Case 2 and    Case 6. Figure 6 shows the entrance shape of each case. In Case 2, the wavy rough surface starts from a hill. In contrast, the wavy rough surface starts from the valley in Case 6. The effect of the differences in the entrance shape on the boundary-layer transition was verified. Figure 7 shows the contour map of the instantaneous value of the wall-normal velocity w in Cases 2 and 6 simultaneously. The figure shows that the two-dimensional roll vortices are generated on the wavy rough surface in each case. We assume that the vortex generation occurs as a transition, and then measure the tendency as the root mean square value of the wall-normal velocity w, obtained from sufficient period data of a saturated flow state. Figure 8 shows the contour map of the root mean square (RMS) value of the wall-normal velocity w. It can be observed from the figure that   the RMS values of w in each case are almost the same. That is, the boundary layer transition occurs at the same position. This result indicates that the wavy rough surface entrance shape has minimal effect on the boundary-layer transition in this study. Therefore, in the other cases except Case 2 and Case 6, the entrance shape of the wavy rough surface was the same as in Case 2.

Effect of wavy rough surface wavelength
The effect of wavy rough surface wavelength on the boundary-layer transition was verified. In this study, we compared Cases 1 to 5 with different wavelengths. The wavelength of each case is as shown in Figure 3. Case 1 has the smallest wavelength of = 3.12, whereas Case 5 has the largest wavelength of = 10. Figure 9 shows the contour map of the instantaneous wall-normal velocity w, and Fig. 10 shows the contour map of the RMS value of w for each case. As shown in Fig. 9, roll vortices were generated in each test case. In Case 5 with the largest wavelength, vortices were generated from the entrance of the roughness surface and were found to be the strongest. As the wavelength decreases, the size and the intensity of a vortex decrease. In Cases 1 and 2 in which the roughness wavelength was relatively smaller, the vortices size and intensity increased intermittently. Therefore, the vortex growth process changes according to the wavelength . Figure 10 shows that increasing the wavelength tends to advance the vortex generation starting position. The position recedes as the wavy surface wavelength decreases. However, although Case 1 has a smaller wavelength than Case 2, the boundary layer transition position of Case 1 has advanced compared with that of Case 2. That is, the boundary layer transition position moved backward when the wavy rough surface wavelength decreased, whereas it moved forward when the wavelength decreased to less than a certain value.
In this study, two-dimensional roll vortices are generated and boundary layer transition occurs on the small wavy rough surface. It is considered that this is due to the two-dimensional instability of a boundary layer. Furthermore, we examined this phenomenon in detail using a neutral stability curve. Figure 11 shows the neutral stability curve of Blasius flow (Kobayashi et al., 2015). The vertical axis is the dimensionless frequency, L (= * * / $ * @ × 10 M ), whereas the horizontal axis is the Reynolds number L , based on the boundary layer displacement thickness (= * / !" * ). Here, we calculated the disturbance frequency L from the wavy surface wavelength by assuming that the roughness produces periodic disturbances in the boundary layer with the convection velocity * / $ * (= 1.0). In addition, the value of the horizontal axis L is calculated using the displacement thickness immediately after the start of the wavy rough surface at x = 120. In this study, the value of the Reynolds number L is 800. Figure 11 shows a plot of each case we calculated. The figure shows that Case 5 is in the "Unstable" region of the neutral stability curve at x = 120. Moreover, as the wavelength decreases, it gradually moves to a stable region. Therefore, we considered that the reason for the early boundary layer transition in Case 5 is that the disturbance generated from the wavy rough surface coincided with Tollmien-Schlichting (TS) instability. Hence, as the first condition to delay a transition, a wavy rough surface wavelength that avoids the frequency that induces TS instability delays the boundary-layer transition.
The above result indicates that as the roughness wavelength decreases in Case 5 and the disturbance generated from the wavy rough surface deviates from the frequency range that induces TS instability, i.e., the "Unstable" region shown in Fig. 11. Further, the TS transition delayed. However, Fig. 10 shows that in Case 2 the wall-normal velocity RMS becomes smaller compared to that of Case 1. Moreover, the instantaneous flow velocity state in Fig. 9 shows that the vortex grows intermittently. It is observed that vortex pairing occurred in downstream in every case. We considered that the intermittent vortex growth observed is due to pairing. The instantaneous snapshots in Fig. 12 shows the vortex pairing characteristics in Case 5. In the figure, the interval between the images is 20,000 steps ( = ∆ * * / !" * = 100). The vortex pairing is considered to be one of the transition processes and causes a rapid vortex-size enlargement in each case. The frequency of the vortices that are advected on the wavy rough surface decreases when vortex pairing occurs. The characteristics of intermittently growing vortices differ between Case 1 and Case 2, as shown in Fig. 9. The number of vortex pairings also seems to be different. Fast Fourier transform (FFT) analysis was performed for Case 1 and Case 2 to compare the occurrence of vortex-pairing. A sufficient time series of wall-normal w data was used in the analysis. The vortices generated on the wavy rough surface grow as they are advected in the streamwise direction. Therefore, the vortex disturbance frequency was calculated at the position where the RMS value of w is the largest, and at all points in the streamwise direction in = 120~1200. Fig. 13 shows the results of the FFT analysis in Case 1 and Case 2. In Fig. 13,

Short
Long the vertical axis F is the dimensionless frequency, which is defined as F (= * !" * / $ * ) and is the most dominant detected.
In Case 1, one of the first vortices of a frequency F of 0.028 is observed at approximately = 200, as shown in Fig.  13 (a). Then, the dominant frequency gradually decreases in the streamwise direction. In contrast, the same vortices are generated around = 200 in Case 2, however, the frequency decreases rapidly in downstream. That is, the vortex pairing occurs at the position where the frequency decreases discontinuously. In Case2, further vortex pairing does not occur thereafter, and the dominant frequency is almost constant. From the above results, we considered that the reason the wall-normal velocity RMS is suppressed in Case 2 compared with Case 1 is because of vortex pairing. In Case 2, the pairing mostly occurs in the upstream first half and very less in the downstream latter half of the wavy rough surface; the vortex size is kept small. That decreases the dominant frequency and the RMS value of w. Therefore, we consider that the boundary layer transition position is delayed as shown in Fig. 10 (b), that is, a certain wavelength surface roughness affects the secondary nonlinear vortex paring and keeps the vortex size small, resulting in decreasing RMS value.

Conclusion
In this study, we performed direct numerical simulations to verify the effect of a wavy rough surface on a twodimensional natural boundary layer transition. First, we examined the effect of the entrance shape of the wavy rough surface. We compared two cases with similar wavelengths and roughness heights. The shape of the entrance of the wavy rough surface was changed by setting the phases such that they are opposite to each other. In each case, the wall-normal velocity RMS value almost coincided, that is, the entrance shape has almost no effect on the natural transition when the roughness height is 23.7% of the thickness of the inflow boundary layer, in this study.
Second, we also verified the effect of the wavelength of the surface roughness. We found that the first twodimensional roll vortices of a boundary layer transition were facilitated when the wavy rough surface wavelength is in the frequency range that induces TS instability. In addition, the vortex pairing was observed in downstream. In the most w-velocity-fluctuation suppressed case, the vortex paring occurred however additional pairing did not occur significantly in downstream and the vortex size on the wavy rough surface was kept small. In conclusion, a small wavy surface affects not only the primary but also the secondary two-dimensional boundary layer transition. The TS instability is encouraged by that coincide roughness wavelength, on the other hand, there is a case in which the w-velocity fluctuation increases when the roughness wavelength is much smaller than that. The roughness wavelength should be avoided a frequency of the TS instability in order to delay the natural boundary layer transition. Additionally, in the optimal wavelength case, the w-velocity fluctuation was significantly suppressed; the vortex pairing occurred in upstream while the vortex size was kept small in downstream, and as a result, the intensity of vortices tended to decrease. The effect of roughness on the vortex pairing will need further investigation. We confirmed that the roughness affected not only the first TS instability but also the second vortex paring in the transition delay, and the roughness shape should be decided based on those characteristics.