Development of an Analytical Wall Function for Bypass Transition

: The objective of the present work is to propose an extended analytical wall function that is capable of predicting the bypass transition from laminar to turbulent ﬂow. The algebraic γ transition model, the k − ω turbulence model and the analytical wall function are integrated together in this work to detect the transition onset and start the transition process. The present analytical wall function is validated with the experimental data, the Blasius solution and the law of the wall. With this analytical wall function, the transition onset in the skin friction coefﬁcient is detected and the growth rate of transition is properly generated. The predicted mean velocity proﬁles are found to be in good agreement with the Blasius solution in the laminar ﬂow, the experimental data in the transition zone and the law of the wall in the fully turbulent ﬂow.


Introduction
Theoretically, turbulence and the laminar-to-turbulent flow transition process have been studied, for example in Balonishnikov [1] and Ershkov [2], respectively. For engineering applications, it was reported by Pacciani et al. [3,4] that transition from laminar to turbulent flow plays a key role in modern aero-engines. This complex flow appears for instance in low-pressure turbines where the number of blades per row is reduced to meet the increasing demand of compact and light aircraft engines. Consequently, the load on each blade increases with relatively low Reynolds number operations. Typically, for industrial flow simulations, the wall function is employed when the Reynolds-averaged Navier-Stokes (RANS) equations are solved with the turbulence model. The transition model is used to account for the effects of transition on the mean flow.
The wall function was first proposed by Patankar and Spalding [5] and later on improved by Launder and Spalding [6]. The wall function is useful because the turbulent flow near the wall behaves differently in three different regions: (1) the viscous sublayer, (2) the buffer layer and (3) the log layer, each of which possesses steep variation with complex interfaces among them, leading to the requirement of prohibitively fine mesh, especially for three-dimensional flow. According to the classification by Hanjalić and Launder [7], the wall function can be categorized into five types: standard wall function (SWF), analytical wall function (AWF), simplified analytical wall function (SAWF), blended wall treatment (BWT), and numerical wall function (NWF). Only the AWF is considered here. The AWF was considered as an advanced modeling approach in Saric et al. [8] that had not yet been extensively used in industrial CFD applications. AWF has been developed Fluids 2021, 6, 328 2 of 12 in order to cope with non-equilibrium effects-for example, strong adverse and favorable pressure gradient, separation, transition, compressibility and heat transfer.
The AWF was originally proposed by Craft et al. [9]. With the assumption of the linear eddy viscosity profile beyond the viscous sublayer in the first control volume adjacent to the wall, the AWF formulation is obtained by the analytical integration of the simplified wall-parallel momentum equation in which the wall-parallel convection and the streamwise pressure gradient are balanced with the wall-normal diffusion, leading to an algebraic formulation of the wall shear stress to be prescribed as a wall-normal diffusive flux at the wall boundary. In this first control volume adjacent to the wall, the transport equation of the turbulent kinetic energy (k) is numerically solved using the volume-averaged formulations of its production and dissipation terms. For the dissipation rate of turbulent kinetic energy, its transport equation is not solved in the first control volume adjacent to the wall but its local formulation is prescribed instead. In this way, the non-equilibrium condition can be treated properly.
The AWF has been continually developed and evaluated in various engineering aspects. Craft et al. [10]) applied it to a 2D downward-directed buoyancy-modified turbulent wall jet. Suga et al. [11] extended its capability to account for the effects of fine-grain surface roughness on turbulent flow and heat transfer. Suga [12] extended its capability to turbulent flow and heat transfer of a high-Prandtl-number fluid with and without roughness. Suga and Nishiguchi [13] applied it in the interface region between a porous wall and a clear fluid in order to bridge the flows inside and outside the porous medium. Suga and Kubo [14] proposed an extended version of AWF to model the mass transfer and the concentration field across undeformable air-water interfaces with and without shear at the interface over a range of Schmidt number 1 ≤ Sc ≤ 1000. Suga et al. [15] modified the coefficient α of the AWF eddy viscosity from a constant to a function of mean strain rate and tested it on heat transfer of backward-facing step flow and an impinging jet. Amano et al. [16] evaluated its performance on turbulent flow and heat transfer in a 3D two-pass cooling channel. Omranian et al. [17] evaluated its performance on the turbulent natural convection flows in various cavity configurations with differentially heated walls. Wang et al. [18] applied it to 2D shock wave/turbulent boundary layer interaction. Chedevergne [19] modified the AWF of Suga et al. [11] using the roughness corrections of Aupoix [20,21] to improve the prediction of turbulent flow and heat transfer over rough walls in the unstructured mesh framework. Saric et al. [8] implemented the AWF formulation of the energy equation (AWF-e) of Suga [22] into the AVL software and tested it in high-Prandtl-number turbulent flows in an IC engine and an E-motor cooling jacket. However, from the literature review, there has been no AWF proposed for transitional flows. The main objective of the present research work is to extend the AWF capability to transitional flow prediction.
From the literature review in Fu and Wang [23], Durbin [24,25] and Dick and Kubacki [26], transition models can be categorized into three classes: two-equation, one-equation and algebraic (zero-equation) models. For two-equation transition models, Langtry and Menter [27] proposed the γ − Re θ transition model, where γ is the intermittency factor and Re θ is the momentum-thickness Reynolds number, while Juntasaro and Ngiamsoongnirn [28] developed the γ − k L transition model, where k L is the kinetic energy of laminar fluctuations. Juntasaro and Narejo [29] further improved the γ − k L transition model to account for pressure gradient effects. Xu et al. [30] proposed another γ − k L transition model. The γ − ν LF transition model of Xu et al. [31] was an extension of the work of Xu et al. [30]; however, the k L -equation was replaced by the ν LF -equation, where ν LF is the eddy viscosity of the laminar fluctuations. For one-equation transition models, a variety of γ transition models were proposed: Lodefier et al. [32], Wang and Fu [33], Durbin [34], Ge et al. [35], Menter et al. [36] and Juntasaro et al. [37]. Besides the intermittency factor γ, Mayle and Schulz [38], Walters and Leylek [39], Walters and Cokljat [40] and Medina et al. [41] developed the k L transition model. As an alternative to the kinetic energy of laminar fluctuations k L , Lopez and Walters [42] proposed the ν 2 transition model, where ν 2 is the wall-normal velocity fluctuation made approximately equal to the difference between the turbulent kinetic energy k and the kinetic energy of laminar fluctuations k L , i.e., ν 2 ≈ k − k L . For zero-equation transition models, Kubacki and Dick [43,44] proposed an algebraic γ transition model that requires no transport equation. Sandhu and Ghosh [45] modified the k-equation of the SST k − ω turbulence model, where k is the turbulent kinetic energy and ω is the specific dissipation rate of k, by (1) multiplying its production and destruction terms by the algebraic γ expression and (2) adding three extra terms (diffusion, source and sink) to account for the effects of transition on the mean flow.
Since the AWF has been constructed to save the computing time for industrial flow simulations, the algebraic γ transition model of Kubacki and Dick [43,44] is thus selected here for its simplicity, compared to other transition models, in order to essentially serve the same economical purpose as the AWF.

Transition and Turbulence Models
The k − ω turbulence model of Wilcox [46] was modified in combination with the algebraic γ transition model by Dick (2016a, 2016b) for the predictive capability of capturing the effect of transition on the mean flow. In this work, some constants of the algebraic γ transition model of Kubacki and Dick [43,44] are re-calibrated for its compatibility with AWF.

k − ω Turbulence Model
The k − ω turbulence model of Wilcox [46] was modified by Kubacki and Dick [43,44] for the prediction of laminar-to-turbulent flow transition as follows: where k is the turbulent kinetic energy, ω is the specific dissipation rate of k, ν is the kinematic viscosity, (σ * , σ, α * , β * , σ d ) are the model constants and β is the model function.
The intermittency factor γ in Equation (1) was used by Kubacki and Dick [43,44] to detect the transition onset and start the transition process by controlling the production term of the k-equation with its formulation, described later in Section 2.2. The production term P k was modified by Kubacki and Dick [43,44] using the small-scale eddy viscosity ν s as: where S is the magnitude of the mean strain rate, whose definition will be provided after Equation (7), and ν s was defined by Kubacki and Dick [43,44] as: The small-scale specific dissipation rate of turbulent kinetic energy ω s in Equation (4) retains the original form of Wilcox [46] as: where C lim and a 1 are the model constants with their standard values. The small-scale turbulent kinetic energy k s in Equation (4) was modeled by Kubacki and Dick [43,44] to take into account a shear-sheltering mechanism with the following definition: where the shear-sheltering function f ss was adopted from Walters and Cokljat [40] by Kubacki and Dick [43] with the following definition: where C ss = 4.5 is the model constant, re-calibrated in this work for compatibility with AWF. In Equation (3) and Equation (5), S = 2S ij S ij is the magnitude of the mean strain rate with S ij = 0.5 ∂U i /∂x j + ∂U j /∂x i where U i and U j are the velocity components. In Equation (7), Ω = 2Ω ij Ω ij is the magnitude of the mean rotation rate with Based on Walters and Cokljat [40], the turbulent kinetic energy k was divided by Kubacki and Dick [43,44]) into two parts: k s (small scale) and k l (large scale). The smallscale part k s was already defined in Equation (6) while the large-scale part k l is calculated by: Similarly, the eddy viscosity ν T used in the momentum equations was also divided by Kubacki and Dick [43,44]) into two parts: ν s (small scale) and ν l (large scale). The small-scale part ν s was already defined in Equation (4) whereas the large-scale part ν l is calculated by: where the large-scale specific dissipation rate of turbulent kinetic energy ω l was proposed by Kubacki and Dick [43,44]) to have the same functional form as ω s as follows: where a 2 is the model constant that was calibrated by Kubacki and Dick [43,44] for the transition prediction. According to Wilcox [46] and Kubacki and Dick [43,44], the model constants and functions used can be summarized as follows: Wilcox [46]: Kubacki and Dick [43,44]:

Algebraic γ Transition Model
The algebraic γ transition model of Kubacki and Dick [43,44] is expressed as: where y * = y √ k/ν, with y being the wall-normal distance and A γ = 45 is the model constant, re-calibrated in this work for compatibility with AWF.

Analytical Wall Function
The analytical wall function (AWF) is used only for the first layer of all near-wall cells (hereafter the near-wall cell P) adjacent to the surface. According to Craft et al. [9], Suga [47] and Hanjalić and Launder [7], the AWF formulations can be summarized as seen below.
3.1. Eddy Viscosity of the Near-Wall Cell P As proposed in Craft et al. [9], the eddy viscosity of the near-wall cell P was assumed to have the following variation: where α = C µ C , C µ = 0.09, C = 2.55, y * ν = 10.7 and y * P = y P k s,P /ν at the center of the near-wall cell P. In Equation (12), k s,P is proposed in this work to calculate y * p .

Wall Shear Stress at Face s (South) of the Near-Wall Cell P
At face s of the near-wall cell P, the wall shear stress that is proposed in this work to account for laminar, transitional and turbulent flows is modeled in the framework of AWF as follows: where ρ is the fluid density, U ∞ is the free-stream velocity, C f ,Blasius = 0.664/ √ Re x is the local skin friction coefficient of Blasius solution with Re x = U ∞ x/ν, k s,P is the small-scale part of k given in Equation (6) at the center of the near-wall cell P, and A U is expressed as: where µ is the dynamic viscosity, U n is the mean velocity at face n (north) of the near-wall cell P, y * n = y n √ k P /ν is the dimensionless wall-normal distance at face n of the near-wall cell P, and C U is expressed as: where P is the pressure, and the first and second terms in the bracket are the wall-parallel convection and the streamwise pressure gradient, respectively. In Equation (14), A U is the constant of integration that appears when the simplified x-momentum equation is integrated in the wall-normal direction y over the near-wall cell P using the eddy viscosity in Equation (12). The formulation of A U in Equation (14) is obtained by satisfying the wall-parallel velocity profiles U of both the viscous sublayer and the turbulent-flow layer with boundary conditions at the wall, at the north face of the near-wall cell P and at the interface between the viscous sublayer and the turbulent-flow layer.

k-Equation of the Near-Wall Cell P
The k-equation of Craft et al.
[9] is directly adopted in this work only in the near-wall cell P. Its volume-averaged sink term is ε, which is calculable as a function of (k, ν, y) as shown in Equation (17). The k-equation of the near-wall cell P is numerically solved using the following volume-averaged production and dissipation terms, respectively: where y * n = y n k s,P /ν and k s,P are proposed in this work to calculate P k , y * ε = 5.1 and C M is expressed as: The volume-averaged intermittency factor γ of the near-wall cell P is proposed in this work as a multiplier of the volume-averaged production term of k-equation P k in Equation (16) to detect the transition onset and start the transition process in AWF. γ is obtained by integrating the intermittency factor γ in Equation (11) over the near-wall cell P as: where the model constant A γ = 30 is proposed in this work.

ω-Equation of the Near-Wall Cell P
The ω-equation of the near-wall cell P is not numerically solved but ω P at the center of the near-wall cell P is prescribed as: where β 1 = 0.0708 and y * ω = 6α/β 1 = 19.4491.

Results and Discussion
The present AWF in combination with the currently modified versions of the k − ω turbulence model of Wilcox [46] and the algebraic γ transition model of Kubacki and Dick [43,44], hereafter AWF-transition, is implemented into an in-house CFD code which is based on the cell-centered finite-volume method. The convection terms are discretized by the QUICK scheme while the second-order central difference scheme is used for the diffusion terms. The SIMPLE algorithm is employed to couple the pressure and the velocity. The Rhie-Chow interpolation is used to handle the collocated grid arrangement. The convergence criteria is set to be equal to 10 −4 .
The standard T3A test case of ERCOFTAC by Coupland [48] is used to validate the proposed AWF-transition formulation for predicting bypass transition. The computational domain and boundary conditions of T3A are illustrated in Figure 1. The structured mesh used is composed of 300 (clustering near the leading edge) × 100 (uniform with ∆y + ∼ = 40 at outlet) in the x-and y-directions, respectively, which is selected after performing gridindependent study. The mesh distribution is displayed in Figure 2. The flow condition of T3A is summarized in Table 1. The standard T3A test case of ERCOFTAC by Coupland [48] is used to validate the proposed AWF-transition formulation for predicting bypass transition. The computational domain and boundary conditions of T3A are illustrated in Figure 1. The structured mesh used is composed of 300 (clustering near the leading edge) × 100 (uniform with ∆ + ≅ 40 at outlet) in the x-and y-directions, respectively, which is selected after performing grid-independent study. The mesh distribution is displayed in Figure 2. The flow condition of T3A is summarized in Table 1.   Using the inlet condition in Table 1, the decay of the free-stream turbulence intensity is well matched between the prediction and the experimental data in the case of T3A, as shown in Figure 3. proposed AWF-transition formulation for predicting bypass transition. The computational domain and boundary conditions of T3A are illustrated in Figure 1. The structured mesh used is composed of 300 (clustering near the leading edge) × 100 (uniform with ∆ + ≅ 40 at outlet) in the x-and y-directions, respectively, which is selected after performing grid-independent study. The mesh distribution is displayed in Figure 2. The flow condition of T3A is summarized in Table 1.   Using the inlet condition in Table 1, the decay of the free-stream turbulence intensity is well matched between the prediction and the experimental data in the case of T3A, as shown in Figure 3.  Using the inlet condition in Table 1, the decay of the free-stream turbulence intensity is well matched between the prediction and the experimental data in the case of T3A, as shown in Figure 3. Figure 4 shows the skin friction coefficient distribution along a flat plate where the AWF-transition result is compared with the experimental data. It reveals that the present AWF-transition can capture the laminar-to-turbulent flow transition, although a coarse mesh is employed.
The mean velocity profile at x = 0.695 m in the transition zone is plotted in wall units (u + = U/u τ , y + = yu τ /ν) and in the same dimensionless form as Blasius solution U/U ∞ , η = y √ U ∞ /νx in Figure 5a,b, respectively, where the AWF-transition result is compared with the experimental data in both scaling forms, including the law of the wall in Figure 5a and the Blasius solution in Figure 5b. It reveals that the present AWF-transition is capable of predicting the mean velocity profile in the transition zone. It can be noticed from the experimental data that the flow is neither fully turbulent nor laminar because the log layer is not yet established and the Blasius solution is no longer valid at this location in the transition zone.   Figure 4 shows the skin friction coefficient distribution along a flat plate where the AWF-transition result is compared with the experimental data. It reveals that the present AWF-transition can capture the laminar-to-turbulent flow transition, although a coarse mesh is employed. The mean velocity profile at = 0.695 in the transition zone is plotted in wall units ( + = / , + = / ) and in the same dimensionless form as Blasius solution ( / ∞ , = √ ∞ / ) in Figure 5a,b, respectively, where the AWF-transition result is compared with the experimental data in both scaling forms, including the law of the wall in Figure 5a and the Blasius solution in Figure 5b. It reveals that the present AWF-transition is capable of predicting the mean velocity profile in the transition zone. It can be noticed from the experimental data that the flow is neither fully turbulent nor laminar because the log layer is not yet established and the Blasius solution is no longer valid at this location in the transition zone.    The mean velocity profile at = 0.695 in the transition zone is plotted in wall units ( + = / , + = / ) and in the same dimensionless form as Blasius solution ( / ∞ , = √ ∞ / ) in Figure 5a,b, respectively, where the AWF-transition result is compared with the experimental data in both scaling forms, including the law of the wall in Figure 5a and the Blasius solution in Figure 5b. It reveals that the present AWF-transition is capable of predicting the mean velocity profile in the transition zone. It can be noticed from the experimental data that the flow is neither fully turbulent nor laminar because the log layer is not yet established and the Blasius solution is no longer valid at this location in the transition zone.  In Figure 6, the velocity profile at = 0.195 in the laminar zone is shown where the AWF-transition result is compared with both the experimental data and the Blasius solution. Even though there are very few grid points in the wall-normal direction inside the laminar boundary layer, the AWF-transition result agrees well with both the experimental data and the Blasius solution. In Figure 6, the velocity profile at x = 0.195 m in the laminar zone is shown where the AWF-transition result is compared with both the experimental data and the Blasius solution. Even though there are very few grid points in the wall-normal direction inside the laminar boundary layer, the AWF-transition result agrees well with both the experimental data and the Blasius solution. In Figure 7, the mean velocity profile at = 0.895 in the fully turbulent zone is plotted in wall units, where the AWF-transition result is compared with both the experimental data and the law of the wall. The AWF-transition result is in good agreement with the experimental data and the law of the wall. At this location in the fully turbulent zone, the log layer is well established as demonstrated by the experimental data and also the AWF-transition result. However, the AWF-transition result is slightly over-predicted compared to the experimental data. In Figure 7, the mean velocity profile at x = 0.895 m in the fully turbulent zone is plotted in wall units, where the AWF-transition result is compared with both the experimental data and the law of the wall. The AWF-transition result is in good agreement with the experimental data and the law of the wall. At this location in the fully turbulent zone, the log layer is well established as demonstrated by the experimental data and also the AWF-transition result. However, the AWF-transition result is slightly over-predicted compared to the experimental data.
In Figure 7, the mean velocity profile at = 0.895 in the fully turbulent zone is plotted in wall units, where the AWF-transition result is compared with both the experimental data and the law of the wall. The AWF-transition result is in good agreement with the experimental data and the law of the wall. At this location in the fully turbulent zone, the log layer is well established as demonstrated by the experimental data and also the AWF-transition result. However, the AWF-transition result is slightly over-predicted compared to the experimental data.

Conclusions
In this work, the capability of AWF was extended to the prediction of the bypasstransition effect on the mean flow. To capture the transition onset and the transition process, the presently modified versions of the algebraic transition model of Kubacki and Dick [43,44] and the − turbulence model of Wilcox [46] were integrated into the proposed AWF-transition of an in-house CFD code. For validation, the T3A test case of bypass transition was used. This AWF-transition showed the predictive capability of detecting the transition onset and the transition process in the distribution of over a flat plate. Moreover, the mean velocity profiles were well predicted in laminar, transitional and fully turbulent regions.

Conclusions
In this work, the capability of AWF was extended to the prediction of the bypasstransition effect on the mean flow. To capture the transition onset and the transition process, the presently modified versions of the algebraic γ transition model of Kubacki and Dick [43,44] and the k − ω turbulence model of Wilcox [46] were integrated into the proposed AWF-transition of an in-house CFD code. For validation, the T3A test case of bypass transition was used. This AWF-transition showed the predictive capability of detecting the transition onset and the transition process in the distribution of C f over a flat plate. Moreover, the mean velocity profiles were well predicted in laminar, transitional and fully turbulent regions.