A new boundary layer integral method based on the universal velocity profile

Cantwell, Brian J.; Bilgin, Eylul; Needels, Jacob T.

doi:10.1063/5.0100367

A recently developed mixing length model of the turbulent shearing stress in wall bounded flows has been used to formulate a universal velocity profile (UVP) that provides an effective replacement for the widely used Coles wall-wake formulation. Comparisons with both direct numerical simulation and experimental data demonstrate the ability of the profile to approximate a wide variety of wall-bounded flows. The UVP is uniformly valid from the wall to the boundary layer edge and for all Reynolds numbers from zero to infinity. There is no presumption of logarithmic dependence of the velocity profile outside the viscous wall layer so the profile can accurately approximate low Reynolds number turbulent boundary layers. The effect of a pressure gradient is included in the UVP through the introduction of a modified Clauser parameter that correlates well with the parameters that determine the wake portion of the velocity profile. The inherent dependence of the UVP on Reynolds number, extended to include the effect of pressure gradient, enables it to be used as the basis of a new method for integrating the von Kármán boundary layer integral equation for a wide variety of attached wall bounded flows. To illustrate its application, the UVP is used to determine the zero-lift drag coefficient of the Joukowsky 0012 and NACA (National Advisory Committee for Aeronautics) 0012 airfoils over a wide range of chord Reynolds numbers.

I. INTRODUCTION

In volume III of the Durand series on Aerodynamic Theory, Prandtl¹ points out that the boundary layer approximation derived on a flat wall can be extended to a curved surface. He writes: “In this case, a coordinate system can be chosen in which the arc length along the fixed surface can be introduced as the x-coordinate, and the perpendicular distance from the surface as the y-coordinate.” The idea is illustrated in Fig. 1.

FIG. 1.

View large Download slide

Curved surface boundary layer coordinate notation. The wing shown is a Joukowsky 0012.

Using this construct, a number of simplified approaches were developed to produce approximate solutions suitable for estimating the drag of planar aerodynamic shapes. For laminar boundary layers, the most well-known integral methods are due to Thwaites² and Pohlhausen.³ The latter was extended recently by Majdalani and Xuan.⁴ There are also a limited number of methods that integrate the Kármán equation for turbulent boundary layers. One of the most widely used is the method due to Head.⁵ He assumes the existence of two universal functions, F and G, for the turbulent boundary layer. The first function relates a modified boundary layer shape factor to the conventional shape factor, and the second relates the conventional shape factor to the entrainment velocity at the outer edge of the boundary layer. When these functions are combined with the friction law of Ludwieg and Tillman,⁶ all of the basic characteristics of the boundary layer can be determined. The main weakness is that the functions F and G are determined from experimental correlations that exhibit very limited universality and the Ludwieg–Tillman friction law is valid only over a narrow range of Reynolds numbers.

Low Mach number turbulent flow over a flat plate at high Reynolds number is described by the incompressible boundary layer approximation (1)

\begin{matrix} \frac{\partial}{\partial x} (u u) + \frac{\partial}{\partial y} (u v) + \frac{\partial}{\partial y} (\bar{u' v'}) - u_{e} (x) \frac{d u_{e}}{d x} - ν \frac{\partial^{2} u}{\partial y^{2}} = 0, \\ \frac{\partial u}{\partial x} + \frac{\partial v}{\partial y} = 0, \end{matrix}

(1)

subject to the no-slip and free stream conditions

u (0) = v (0) = 0, u (δ_{h}) = u_{e},

(2)

where u, u_e, and v are mean velocities.

In the region near the stagnation point, the boundary layer approximation breaks down (Fig. 2). Outside of this region, if the surface radius of curvature is large compared to the boundary layer thickness, no significant pressure variation will occur in the wall-normal direction and wall curvature plays no role in the boundary layer approximation.⁷ Theoretical underpinning of Prandtl's idea can be found in the invariance of the boundary layer approximation, (1), under the Lie translation group with an arbitrary function of x added to the y coordinate (Cantwell,⁸ Sec. 10.5, p. 301).

FIG. 2.

View large Download slide

Coordinate notation near the wing leading edge.

The free stream velocity $u_{e} (x)$ is determined by solving Laplace's equation for the potential flow about the body with the Kutta condition applied at the trailing edge. With $u_{e} (x)$ known, along with a suitable model of the turbulent shear stress $\bar{u' v'}$ ⁠, the calculation of boundary layer characteristics can proceed.

Boundary layer integral methods all use the von Kármán⁹ equation (3) derived by integrating the boundary layer equations (1) over the height of the boundary layer

\frac{d δ_{2}}{d x} + (2 δ_{2} + δ_{1}) \frac{1}{u_{e}} \frac{d u_{e}}{d x} - {(\frac{u_{τ}}{u_{e}})}^{2} = 0.

(3)

The friction velocity is defined as

u_{τ} (x) \equiv {(\frac{τ_{w}}{ρ})}^{1 / 2} .

(4)

The displacement thickness is defined as

δ_{1} (x) = \int_{0}^{δ_{h}} (1 - \frac{u}{u_{e}}) d y .

(5)

and the momentum thickness is

δ_{2} (x) = \int_{0}^{δ_{h}} \frac{u}{u_{e}} (1 - \frac{u}{u_{e}}) d y,

(6)

where δ_h indicated in Fig. 1 is a boundary layer thickness defined in Cantwell¹⁰ as the equivalent channel half height thickness. The thickness, δ_h, will be discussed further in Sec. II B.

Wing drag is calculated using the construction in Fig. 3. The differential force acting in the $\tilde{x}$ direction on, say, the wing upper surface is

d F_{\tilde{x}} = τ_{w} (\tilde{x}) d x cos (θ) = τ_{w} (\tilde{x}) d \tilde{x} .

(7)

FIG. 3.

View large Download slide

Drag force on a differential section of the wing surface.

The viscous drag coefficient based on the chord of a symmetric wing of the type shown in Fig. 1 is

C_{d_{v}} = \frac{D_{v}}{\frac{1}{2} ρ u_{\infty}^{2} c} = 2 \int_{0}^{1} {(\frac{u_{e}}{u_{\infty}})}^{2} C_{f} d (\frac{\tilde{x}}{c}),

(8)

where the local (based on u_e) friction coefficient is

C_{f} = \frac{τ_{w}}{\frac{1}{2} ρ u_{e}^{2}} = 2 {(\frac{u_{τ}}{u_{e}})}^{2} .

(9)

II. THE UNIVERSAL VELOCITY PROFILE (UVP)

The universal velocity profile comes from the solution of the fully developed channel (or pipe) flow equation

\frac{\partial}{\partial y} \bar{u' v'} + \frac{1}{ρ} \frac{d p_{e} (x)}{d x} - ν \frac{\partial^{2} u}{\partial y^{2}} = 0,

(10)

where the pressure gradient is a constant, Cantwell.¹¹

Express Eq. (10) in wall normalized coordinates

y^{+} = \frac{y u_{τ}}{ν}; u^{+} = \frac{u}{u_{τ}}; τ^{+} = - \frac{\bar{u' v'}}{u_{τ}^{2}} .

(11)

Integrating (10) with respect to y⁺ from the lower wall of the channel, the result is the stress balance

τ^{+} + \frac{d u^{+}}{d y^{+}} - (1 - \frac{y^{+}}{R_{τ}}) = 0,

(12)

where $R_{τ}$ is the friction Reynolds number

R_{τ} = \frac{u_{τ} δ_{h}}{ν} .

(13)

The turbulent shear stress is modeled using classical mixing length theory (von Kármán;¹² see also Prandtl¹³ and van Driest¹⁴). Let

τ^{+} = {(λ (y^{+}) \frac{d u^{+}}{d y^{+}})}^{2} .

(14)

When (14) is substituted into (12), the result is a quadratic equation for $d u^{+} / d y^{+}$ ⁠. We use the positive root

\frac{d u^{+}}{d y^{+}} = - \frac{1}{2 λ {(y^{+})}^{2}} + \frac{1}{2 λ {(y^{+})}^{2}} {(1 + 4 λ {(y^{+})}^{2} (1 - \frac{y^{+}}{R_{τ}}))}^{1 / 2} .

(15)

The singularity at the wall where the mixing length, $λ$ ⁠, goes to zero can be easily removed as follows.^15,16 Multiply (15) by one,

\begin{matrix} \frac{d u^{+}}{d y^{+}} = (- \frac{1}{2 λ {(y^{+})}^{2}} + \frac{1}{2 λ {(y^{+})}^{2}} {(1 + 4 λ {(y^{+})}^{2} (1 - \frac{y^{+}}{R_{τ}}))}^{1 / 2}) \\ \times \frac{\frac{1}{λ {(y^{+})}^{2}} + \frac{1}{λ {(y^{+})}^{2}} {(1 + 4 λ {(y^{+})}^{2} (1 - \frac{y^{+}}{R_{τ}}))}^{1 / 2}}{\frac{1}{λ {(y^{+})}^{2}} + \frac{1}{λ {(y^{+})}^{2}} {(1 + 4 λ {(y^{+})}^{2} (1 - \frac{y^{+}}{R_{τ}}))}^{1 / 2}} . \end{matrix}

(16)

Carry out the multiplication in (16) and canceling terms, the result, equivalent to (15), but with the singularity removed is

\frac{d u^{+}}{d y^{+}} = \frac{2 (1 - \frac{y^{+}}{R_{τ}})}{1 + {(1 + 4 λ {(y^{+})}^{2} (1 - \frac{y^{+}}{R_{τ}}))}^{1 / 2}},

(17)

where $\frac{d u^{+}}{d y^{+}} (0) = 1$ and $\frac{d u^{+}}{d y^{+}} (R_{τ}) = 0$ ⁠. The limiting velocity gradients at the wall and centerline are

\begin{matrix} \lim_{\begin{matrix} y^{+} \to 0 \\ λ \to 0 \end{matrix}} \frac{d u^{+}}{d y^{+}} = 1 - \frac{y^{+}}{R_{τ}}, \\ \lim_{y^{+} \to R_{τ}} \frac{d u^{+}}{d y^{+}} = 1 - \frac{y^{+}}{R_{τ}}, \end{matrix}

(18)

independent of the choice of $λ (y_{+})$ ⁠.

Equation (17) is integrated from the wall to y⁺ to obtain the velocity profile in the form of an integral dependent on the non-dimensional mixing length function $λ (y^{+})$ at a given $R_{τ}$

\begin{matrix} u^{+} (y^{+}) = \int_{0}^{y^{+}} \frac{2 (1 - \frac{s}{R_{τ}})}{1 + {(1 + 4 λ {(s)}^{2} (1 - \frac{s}{R_{τ}}))}^{1 / 2}} d s \end{matrix},

(19)

where $u^{+} (0) = 0$ and $u^{+} (R_{τ}) = u_{e} / u_{τ}$ ⁠. At low Reynolds number, $λ \to 0$ and Eq. (19) approaches the laminar pipe/channel flow solution

\lim_{R_{τ} \to 0} u^{+} = y^{+} (1 - \frac{y^{+}}{2 R_{τ}}),

(20)

where, in the laminar limit, $R_{τ} = {(2 u_{e} δ_{h} / ν)}^{1 / 2}$ ⁠.

The mixing length model introduced by Cantwell¹¹ to approximate pipe data is

λ (y^{+}) = \frac{k y^{+} (1 - e^{- {(\frac{y^{+}}{a})}^{m}})}{{(1 + {(\frac{y^{+}}{b R_{τ}})}^{n})}^{1 / n}} .

(21)

The model (21) contains five free parameters. The constant k is closely related to the Kármán constant. The parameter a constitutes a wall damping length scale. The wall model is similar to the exponential decay proposed by van Driest¹⁴ except for the exponent m that helps to determine the shape and thickness of the near wall velocity profile. Near the wall

\begin{matrix} \lim_{y^{+} \to 0} λ & = k y^{+} (1 - e^{- {(\frac{y^{+}}{a})}^{m}}) . \end{matrix}

(22)

The outer flow term in the denominator of (21) includes a length scale b proportional to the fraction of the boundary layer thickness where wake-like behavior begins as well as an exponent n that helps shape the outer part of the profile.

The reason for making a distinction between k and the Kármán constant is that the UVP accurately approximates low Reynolds number turbulent velocity profiles that have no distinct logarithmic layer and a fairly complicated dependence on k. Whereas at high Reynolds number (⁠ $R_{τ} > 2000 / k$ ⁠),¹¹ k becomes simply a scale factor on u⁺, y⁺ and the damping length scale, a. In that regime, a distinct logarithmic dependence in the velocity profile begins to emerge and k can be viewed as equivalent to the usual definition of the Kármán constant, κ.

The UVP provides an effective replacement for the classical wall-wake formulation.¹⁷ The velocity profile (19) with the mixing length model (21) is uniformly valid over $0 \leq y \leq δ_{h}$ and $0 \leq R_{τ} < \infty$ ⁠. Therefore, there is no need for a buffer layer function and there is no discontinuity in the velocity derivative at the outer edge of the boundary layer. At very low Reynolds number, the velocity profile reverts to the laminar solution. The profile is directly connected to a model of the turbulent shear stress that can be used in computations based on the full Reynolds Averaged Navier–Stokes equations. Finally, the UVP and the shear stress model can be used to study turbulence characteristics such as turbulent kinetic energy (TKE) production.

A. The high Reynolds number explicit form of the UVP

The classical wall-wake formulation is silent as to the Reynolds number at which the formulation actually applies. Whereas for $R_{τ} \geq 2000 / k$ ⁠, the integral form of the UVP, Eq. (19) approaches the explicit form^10,11

\begin{matrix} \frac{u}{u_{τ}} = \frac{1}{k} ln (k y^{+}) + \frac{1}{k} ϕ (k a, m, b, n, \frac{y}{δ_{h}}) . \end{matrix}

(23)

Importantly, the shape function, $ϕ$ ⁠, is independent of $R_{τ}$ ⁠. As long as $R_{τ}$ is above about $2000 / k$ ⁠, the explicit form of the UVP can be used to replace the integral form (19) outside the viscous sublayer and buffer layer, $y^{+} > 132$ (see Table II in Subrahmanyam et al.¹⁷).

The shape function for the zero pressure gradient (ZPG) boundary layer using the parameter values in Table I is shown in Fig. 4. A best fit eighth-order polynomial is overlaid on $ϕ$ ⁠. In Sec. IV, the shape function will be extended to include pressure gradient effects and the eighth-order fit will be extended to a bivariate polynomial. In practice, the shape function is determined by simply evaluating the UVP at any $R_{τ} \geq 2000 / k$ ⁠, say $R_{τ} = 10^{5}$ or 10⁶, to produce $ϕ$ ⁠. Any velocity profile with $R_{τ} \geq 2000 / k$ can then be generated by integrating Eq. (19) only out to about $y^{+} = 132$ ⁠.¹⁷ The rest of the profile out to $y^{+} = R_{τ}$ is generated using Eq. (23). Complete details are given in Cantwell^10,11 and Subrahmanyam et al.¹⁷

TABLE I.

Average model parameters with standard deviation for basic wall flows. Ranges of $R_{τ}$ for each flow are as follows: Pipe (3327–530 023), Channel (550–8016), ZPG boundary layer (1343–17 207).

Flow	$\bar{k}$	σ_k	$\bar{a}$	σ_a	$\bar{m}$	σ_m	$\bar{b}$	σ_b	$\bar{n}$	σ_n
Pipe (21 profiles)	0.4092	0.0057	20.0950	0.381	1.6210	0.0379	0.3195	0.0157	1.6190	0.1204
Channel (7 profiles)	0.4086	0.0179	22.8673	1.599	1.2569	0.0292	0.4649	0.0485	1.3972	0.1213
ZPG boundary layer (11 profiles)	0.4233	0.0068	24.9583	0.663	1.1473	0.0373	0.1752	0.0060	2.1707	0.2238

Flow	$\bar{k}$	σ_k	$\bar{a}$	σ_a	$\bar{m}$	σ_m	$\bar{b}$	σ_b	$\bar{n}$	σ_n
Pipe (21 profiles)	0.4092	0.0057	20.0950	0.381	1.6210	0.0379	0.3195	0.0157	1.6190	0.1204
Channel (7 profiles)	0.4086	0.0179	22.8673	1.599	1.2569	0.0292	0.4649	0.0485	1.3972	0.1213
ZPG boundary layer (11 profiles)	0.4233	0.0068	24.9583	0.663	1.1473	0.0373	0.1752	0.0060	2.1707	0.2238

FIG. 4.

View large Download slide

The universal velocity profile shape function using the parameter values in Table I evaluated at $R_{τ} = 10^{6}$ overlaid by an eighth-order polynomial fit (dashed line).

In this way, the UVP can be evaluated up to arbitrary large Reynolds numbers enabling the structure of wall bounded flows to be explored in the limit of infinite Reynolds number.^10,11,15,18

The high Reynolds number friction law is generated by evaluating (23) at the boundary layer edge

\begin{matrix} \frac{u_{e}}{u_{τ}} = \frac{1}{k} l n (k R_{τ}) + \frac{1}{k} ϕ (k a, m, b, n, 1), \end{matrix}

(24)

shown as the dashed line in Fig. 5. The blue curve in Fig. 5 is the exact friction law: Eq. (19) evaluated at $y^{+} = R_{τ}$ ⁠. Note that the two curves merge around $R_{τ} = 2000 / k$ ⁠.

FIG. 5.

View large Download slide

The universal velocity profile friction law for a zero pressure gradient boundary layer. The dashed line is the high Reynolds number limit generated by evaluating Eq. (23) at $y^{+} = R_{τ}, y / δ_{h} = 1$ ⁠.

The integral method introduced here relies on the UVP to define the velocity profile over the surface in question. In Subrahmanyam et al.,^17,19 the UVP is shown to provide an accurate approximation to zero pressure gradient (ZPG) boundary layer computational data from Simens et al.,²⁰ Borrell et al.,²¹ Sillero et al.,²² and Eitel-Amor et al.²³ The UVP accurately approximates the ZPG experimental data of Klebanoff,²⁴ DeGraaff and Eaton,²⁵ and Baidya et al.²⁶ as well as the adverse pressure gradient data of Perry and Marusic.²⁷ In the present paper, the UVP will be shown to accurately approximate the favorable pressure gradient (sink flow) experimental data of Jones.²⁸

In fact, as noted earlier, the velocity profile that accurately approximates all these flows is fundamentally the same channel flow solution. One flow is distinguished from another only by the values of the five empirical constants $(k, a, m, b, n)$ that characterize it. Since $\partial u / \partial y = 0$ at the edge of the boundary layer, the UVP defines an unambiguous overall thickness for the boundary layer; the equivalent channel half height, δ_h, referred to earlier. For a given flow geometry, the Reynolds number dependence of the model parameters $(k, a, m, b, n)$ tends to be quite weak and so average values of the model parameters can provide a good approximation to the velocity profile over a wide range of Reynolds numbers. In Cantwell,¹⁰ a boundary layer integral method utilizing the universal velocity profile was used to generate the integral properties of the zero pressure gradient boundary layer on a wall over the full Reynolds number range $0 \leq R_{τ} < \infty$ ⁠. In the present paper, this method is generalized to boundary layer flows with pressure gradients and the parameters of the mixing length model used to approximate the wake portion of the velocity profile in the UVP are allowed to vary depending on the pressure gradient. Importantly, Subrahmanyam et al.¹⁷ showed that the wall parameters (k, a, m) exhibit very little change in the presence of an adverse pressure gradient.

To facilitate the later discussion, four functions will be defined here; F₀, F₁, F₂, and F₃. As noted above, the boundary layer friction law is generated by evaluating the UVP at $y^{+} = R_{τ}$

\begin{matrix} \frac{u_{e}}{u_{τ}} = u^{+} (R_{τ}) = \int_{0}^{R_{τ}} \frac{2 (1 - \frac{y^{+}}{R_{τ}})}{1 + {(1 + 4 λ {(y^{+})}^{2} (1 - \frac{y^{+}}{R_{τ}}))}^{1 / 2}} d y^{+} . \equiv F_{0} (R_{τ}) . \end{matrix}

(25)

This is the blue curve in Fig. 5 just discussed. The local displacement thickness Reynolds number expressed in wall units is

R_{δ_{1}} = \frac{u_{e} δ_{1}}{ν} = \frac{u_{e}}{u_{τ}} \int_{0}^{R_{τ}} (1 - \frac{u_{τ}}{u_{e}} u^{+}) d y^{+} \equiv F_{1} (R_{τ}),

(26)

and the local momentum thickness Reynolds number is written as

R_{δ_{2}} = \frac{u_{e} δ_{2}}{ν} = \int_{0}^{R_{τ}} u^{+} (1 - \frac{u_{τ}}{u_{e}} u^{+}) d y^{+} \equiv F_{2} (R_{τ}) .

(27)

For the UVP integral method, we will also need the derivative of the momentum thickness Reynolds number

\frac{d F_{2}}{d R_{τ}} \equiv F_{3} (R_{τ}) .

(28)

B. Optimal parameters

The model parameters $(k, a, m, b, n)$ for the mixing length function (21) are selected by minimizing the sum of total squared error between a given data profile and the universal velocity profile (19) using the cost function

G = \sum_{i = 1}^{N} {(u^{+} (k, a, m, b, n, R_{τ}, y_{i}^{+}) - u_{i}^{+} (y_{i}^{+}))}^{2} .

(29)

Parameter values $(\bar{k}, \bar{a}, \bar{m}, \bar{b}, \bar{n})$ averaged over a wide range of Reynolds numbers for smooth wall channel, pipe and zero-pressure-gradient boundary layer flow are shown in Table I along with the variances of each parameter over the related dataset. The corresponding ranges of $R_{τ}$ are given in the caption. Throughout the present paper, the ZPG boundary layer model constants are assumed to be the boundary layer average values shown in Table I, $(\bar{k}, \bar{a}, \bar{m}, \bar{b}, \bar{n}) = (0.4233, 24.9583, 1.1473, 0.1752, 2.1707)$ from Table I in Subrahmanyam et al.¹⁷ Generally speaking the parameters $(k, a, m, b, n)$ depend, at most weakly on the Reynolds number for a given flow geometry and at high Reynolds number the dependence appears to be especially weak. The evidence for this is the analysis of $(k, a, m, b, n)$ for the Princeton Superpipe data¹¹ that spans three orders of magnitude in the Reynolds number.²⁹ Table I includes the variance in the optimal parameters over the set of profiles for each flow giving a sense of the variation over the ranges of $R_{τ}$ for each flow geometry. For the zero pressure gradient boundary layer, the variation in k, a, m, and b is less than $4 %$ while the variation of n is about $10 %$ ⁠.

The statistical variation in parameter values indicated in Table II implies a similar variation in the wall friction derived from the UVP. Figures 6 and 7 provide a relatively simple measure of the extent of the variation in friction associated with variations in the parameter values by simply showing the friction law with parameter values increased by one standard deviation and decreased by one standard deviation with no attempt to account for any correlation that might exist between parameters. The variation in friction coefficient shown in Fig. 7 translates directly to a comparable variation in calculated drag. For later reference, at $R_{τ} = 5000$ ⁠, the friction coefficient based on the mean parameters in Table I is $C_{f} = 0.002 378$ ⁠. With the mean parameters all increased by one standard deviation $C_{f^{+}} = 0.002 463$ and with the parameters decreased by one standard deviation $C_{f^{-}} = 0.002 277$ for an overall variation of approximately $7.5 %$ ⁠. In Sec. VIII where airfoil drag comparisons are discussed, this uncertainty in C_f will be discussed again.

TABLE II.

Run data, Reynolds number, optimal model parameters, and RMS error for adverse pressure gradient boundary layer datasets from Perry and Marusic.²⁷ Initial free stream values are $u_{\infty} = 10$ and $u_{\infty} = 30$ m/s. Channel half height thickness for these data are at $u = 0.998 u_{e}$ ⁠; $δ_{h} = δ_{0.998}$ ⁠. Reprinted with permission from Subrahmanyam et al., J. Fluid Mech. 933, A16 (2022). Copyright 2022 Author(s), licensed under a Creative Commons Attribution (CC BY) License.¹⁷

x (m)	$u_{e} (\frac{m}{s})$	$R_{δ_{1}}$	$R_{δ_{2}}$	β	β_c	$δ_{998} (m)$	$R_{δ_{998}}$	$R_{τ}$	$(\frac{u_{e}}{u_{τ}})$	k	a	m	b	n	$u^{+}_{rms}$
1.20	10.361	3 165	2 282	0.0	0.0	0.031 79	21 439	912	23.51	0.4287	25.18	1.1528	0.212 2	2.111	0.120
1.80	9.976	5 226	3 734	0.65	1.115	0.050 19	32 606	1285	25.37	0.4239	24.59	1.1583	0.183 9	1.705	0.253
2.24	9.256	6 410	4 342	1.45	2.432	0.055 43	33 456	1195	28.00	0.4223	24.70	1.1460	0.123 7	2.461	0.163
2.64	8.588	8 606	5 517	2.90	4.760	0.070 55	39 406	1252	31.47	0.4255	24.79	1.1121	0.093 31	2.399	0.152
2.88	8.155	11 235	6 879	4.48	7.223	0.086 34	46 043	1337	34.44	0.4296	25.54	1.1173	0.073 83	3.178	0.220
3.08	7.896	12 397	7 213	7.16	11.326	0.092 63	47 598	1248	38.13	0.4302	25.20	1.0938	0.061 51	4.578	0.231
1.20	30.704	8 772	6 564	0.0	0.0	0.033 53	64 807	2461	26.34	0.4095	23.21	1.1161	0.228 5	1.743	0.0705
1.80	29.054	12 401	9 073	0.71	1.230	0.044 15	80 849	2870	28.17	0.4088	23.25	1.1468	0.172 2	1.922	0.0895
2.24	27.035	16 307	11 587	1.39	2.378	0.055 26	94 275	3137	30.05	0.4112	23.42	1.1710	0.133 1	2.335	0.0942
2.64	25.150	21 634	14 736	2.74	4.606	0.069 68	110 700	3373	32.82	0.4035	22.87	1.1352	0.106 6	2.663	0.0984
2.88	23.885	25 854	17 020	3.96	6.567	0.080 54	121 760	3471	35.08	0.4018	22.80	1.1213	0.090 7	3.164	0.1183
3.08	22.908	31 767	20 052	6.07	9.901	0.093 73	136 290	3587	37.99	0.4083	23.53	1.1339	0.074 2	4.069	0.1673

x (m)	$u_{e} (\frac{m}{s})$	$R_{δ_{1}}$	$R_{δ_{2}}$	β	β_c	$δ_{998} (m)$	$R_{δ_{998}}$	$R_{τ}$	$(\frac{u_{e}}{u_{τ}})$	k	a	m	b	n	$u^{+}_{rms}$
1.20	10.361	3 165	2 282	0.0	0.0	0.031 79	21 439	912	23.51	0.4287	25.18	1.1528	0.212 2	2.111	0.120
1.80	9.976	5 226	3 734	0.65	1.115	0.050 19	32 606	1285	25.37	0.4239	24.59	1.1583	0.183 9	1.705	0.253
2.24	9.256	6 410	4 342	1.45	2.432	0.055 43	33 456	1195	28.00	0.4223	24.70	1.1460	0.123 7	2.461	0.163
2.64	8.588	8 606	5 517	2.90	4.760	0.070 55	39 406	1252	31.47	0.4255	24.79	1.1121	0.093 31	2.399	0.152
2.88	8.155	11 235	6 879	4.48	7.223	0.086 34	46 043	1337	34.44	0.4296	25.54	1.1173	0.073 83	3.178	0.220
3.08	7.896	12 397	7 213	7.16	11.326	0.092 63	47 598	1248	38.13	0.4302	25.20	1.0938	0.061 51	4.578	0.231
1.20	30.704	8 772	6 564	0.0	0.0	0.033 53	64 807	2461	26.34	0.4095	23.21	1.1161	0.228 5	1.743	0.0705
1.80	29.054	12 401	9 073	0.71	1.230	0.044 15	80 849	2870	28.17	0.4088	23.25	1.1468	0.172 2	1.922	0.0895
2.24	27.035	16 307	11 587	1.39	2.378	0.055 26	94 275	3137	30.05	0.4112	23.42	1.1710	0.133 1	2.335	0.0942
2.64	25.150	21 634	14 736	2.74	4.606	0.069 68	110 700	3373	32.82	0.4035	22.87	1.1352	0.106 6	2.663	0.0984
2.88	23.885	25 854	17 020	3.96	6.567	0.080 54	121 760	3471	35.08	0.4018	22.80	1.1213	0.090 7	3.164	0.1183
3.08	22.908	31 767	20 052	6.07	9.901	0.093 73	136 290	3587	37.99	0.4083	23.53	1.1339	0.074 2	4.069	0.1673

FIG. 6.

View large Download slide

The blue line is the UVP friction law for a zero pressure gradient boundary layer using mean parameters $(\bar{k}, \bar{a}, \bar{m}, \bar{b}, \bar{n})$ from Table II. The dashed line is the UVP friction law with each parameter increased by one standard deviation, $(\bar{k} + σ_{k}, \bar{a} + σ_{a}, \bar{m} + σ_{m}, \bar{b} + σ_{b}, \bar{n} + σ_{n})$ ⁠. The dotted line is the UVP friction law with each parameter decreased by one standard deviation, $(\bar{k} - σ_{k}, \bar{a} - σ_{a}, \bar{m} - σ_{m}, \bar{b} - σ_{b}, \bar{n} - σ_{n})$ ⁠.

FIG. 7.

View large Download slide

Same as Fig. 6 but expressed in terms of the friction coefficient, $C_{f} = 2 / {(u_{e} / u_{τ})}^{2}$ ⁠.

It should be noted that Fig. 5 and Table II can be compared to Fig. 6 and Table I in Cantwell.¹⁰ The average values of $(k, a, m, b, n)$ in Cantwell¹⁰ differ slightly from the ZPG parameter values used in the present paper but are the same as in the earlier conference paper by Subrahmanyam et al.¹⁹ The reason is that boundary layer data from Baidya et al.²⁶ is included in the average values quoted in Subrahmanyam et al.¹⁷ but was not used in determining average optimal parameters in the 2021 papers where it was decided to use direct numerical simulation (DNS) data only, since it reaches all the way to the wall. In retrospect that condition is not necessary, although experimental data do tend to have higher uncertainty. As new data become available, the average optimal parameters will continue to be updated although changes can be expected to be relatively small as time goes on and the sample size grows.

In the presence of a pressure gradient, the parameters $(b, n)$ will be allowed to vary while (k, a, m) will be held fixed at the average values in Table I. The key question of how $(b, n)$ are related to the pressure gradient will be addressed in Sec. IV.

1. The equivalent channel half height thickness

In Fig. 1, the boundary layer thickness, δ_h, is identified as the equivalent channel half height thickness. The idea behind δ_h originates in the fact that the UVP accurately approximates not only pipe flow and channel flow but also boundary layer flows with a zero, positive and negative pressure gradient.^10,11,17

The UVP has a well-defined outer edge, where $u / u_{e} = 1$ and $\partial u / \partial y = 0$ ⁠. However, these properties of the boundary layer profile are only approached asymptotically. The usual practice is to choose an arbitrary point where $u / u_{e}$ is slightly less than one and $\partial u / \partial y$ is slightly greater than zero. Probably the most common choice is the boundary layer thickness corresponding to $u_{e} / U_{freestream} = 0.99$ ⁠, where U_freestream is the free stream velocity reported with the data for a given profile.

The UVP can be used to remove some of this arbitrariness. When optimal parameters for the UVP are determined for boundary layer data, the choice of cutoff point for the boundary layer edge has an impact on the rms error. If δ_h is too large, data points approaching the free stream where the velocity is virtually constant will be included in the profile. The optimization procedure will try to fit these points and this will tend to reduce the accuracy over the whole profile leading to increased error. If δ_h is too small, the data will be cutoff too far short of the boundary layer edge. The condition $\partial u / \partial y = 0$ will be imposed where the derivative is not quite zero. This also leads to increased error. There is a choice of δ_h that minimizes the error between the UVP and the data. See pp. 12 and 13 and Figs. 17–23 in Cantwell¹⁰ where δ_h is first defined.

At present, techniques for measuring and reporting data vary widely among researchers and so the boundary layer edge velocity corresponding to δ_h also varies. Generally, the best fit occurs at $u_{e} = (0.993 - 0.998) U_{freestream}$ ⁠.

C. Evolution of the UVP from the laminar to the turbulent profile as $R_{τ}$ increases from 0 to $\infty$

Figure 8 shows how the shape of the universal velocity profile becomes fuller as $R_{τ}$ is increased by 7 orders of magnitude. At the lowest Reynolds numbers, the profile is essentially the laminar channel/pipe velocity profile. At the highest Reynolds number, large increases produce very small changes in the velocity profile. At extreme Reynolds numbers, the velocity profile approaches plug flow but astronomically large values are needed to approach this state.¹¹

FIG. 8.

View large Download slide

The universal velocity profile at $R_{τ} = 1, 10, 100, 1000, 10^{4}, 10^{5}, 10^{6}$ ⁠, and $10^{7}$ ⁠. Parameters $(k, a, m, b, n)$ are the ZPG boundary layer values from Table I.

Figure 5 shows the UVP friction law along with low Reynolds number and high Reynolds number limits for a zero pressure gradient boundary layer. Based on Fig. 5, the profile appears to begin becoming fuller at about $R_{τ} = 30$ corresponding to $R_{δ_{1}} = 122, R_{δ_{2}} = 51$ and $R_{x} =$ 15 700 with a shape factor H = 2.39. The profile approaches a fully turbulent shape at about $R_{τ} = 500$ corresponding to $R_{δ_{1}} = 2030, R_{δ_{2}} = 1373$ and $R_{x} = 645 000$ with a shape factor H = 1.48. These changes in the profile inferred from Fig. 5 approximately correspond to the first and third profiles in Fig. 8. The evolution of the profile encompasses the nominal transition range from the critical Reynolds number of the Blasius profile, $R_{τ} = 121, R_{δ_{1}} = 520$ ⁠, and $R_{x} = 104 000$ ⁠,⁷ to a nominal transition Reynolds number, $R_{δ_{1}} = 1600, R_{τ} = 392$ and $R_{x} = 474 000$ ⁠. In other words, the UVP begins to evolve at a Reynolds number well below the instability limit.

The critical Reynolds number for a smooth plate is highly dependent on the free stream turbulence level usually denoted T_u. The question of transition in the UVP boils down to determining what transition assumption one is making when the UVP is used throughout the flow from a laminar leading edge to a fully turbulent boundary layer. To explore this question, the UVP friction law is plotted again in Fig. 9 along with transition data from Schubauer and Klebanoff³⁰ and Coupland.³¹ Note that large $u_{e} / u_{τ}$ implies low C_f. Extremely low values of T_u can lead to critical Reynolds numbers exceeding $R_{τ} = 2000$ corresponding to approximately $R_{e x} = 2 700 000$ while the UVP implies turbulent values of the friction coefficient beginning at about $R_{τ} \approx 200$ corresponding to $R_{e x} \approx 200 000$ ⁠.

FIG. 9.

View large Download slide

The universal velocity profile friction law for a zero pressure gradient boundary layer along with transitional friction data from Schubauer and Klebanoff³⁰ (natural transition $T_{u} = 0.03 %$ —black) and Coupland³¹ (ERCOFTAC case T3A $T_{u} = 0.9 %$ —red, case T3A $T_{u} = 3.0 %$ —magenta).

Methods for modeling the transition region on a smooth airfoil, for example, Kaynak et al.,³² available in the aerodynamic simulation suite SU2,³³ correlate the transition distance with T_u. In order to accurately determine the drag of a given smooth airfoil, the free stream turbulence needs to be well characterized. In the present paper, in the absence of such data, we will assume the UVP holds beginning at the forward stagnation point while recognizing that that assumes that the critical Reynolds number is close to the minimum required to sustain turbulent flow, $R_{τ_{transition}} \approx 200$ ⁠. In other words, the UVP corresponds to a tripped boundary layer.

Figure 9 can be compared to Fig. 14 in Cantwell¹¹ where the UVP friction law for pipe flow is compared with data over 7 orders of magnitude in the pipe Reynolds number. Based on the discussion here, one can conclude that that comparison is for a UVP pipe profile corresponding to a disturbed pipe entry flow.

D. Friction law comparisons

In view of the accurate approximation of the UVP to DNS and experimental data, it should be expected that the UVP generates an accurate friction law. This is crucial if the UVP integral method presented here is to be useful and it needs to be confirmed.

Figure 10 shows the UVP friction law from Figs. 5 and 9. Above $R_{τ} = 2000 / k$ ⁠, the UVP friction law is

\frac{u_{e}}{u_{τ}} = \frac{1}{0.423 3} ln (R_{τ}) + 8.907 74.

(30)

FIG. 10.

View large Download slide

Five friction laws compared to the UVP (blue line). The orange line is the Ludwieg–Tillman law⁶ used in Head's method. The magenta line is the friction law developed by Nash.³⁴ The red and black lines are two versions of the Coles–Fernholz law with $(k, C) = (0.41, 5.0)$ (red line) and $(k, C) = (0.384, 4.1)$ (black line). The brown line is the Spalart–Allmaras law, Eq. (34), deduced from data in Polewski and Cizmas.³⁶

In Fig. 10, the UVP law is compared with five other laws. The Ludwieg–Tillman relation⁶ used in Head's integral method is

\frac{u_{e}}{u_{τ}} = (\frac{2}{0.246}) e^{1.561 H {(R_{δ_{2}})}^{0.268}},

(31)

where H is the shape factor, $H = δ_{1} / δ_{2}$ ⁠.

Also shown in Fig. 10 is the relation developed by Nash.³⁴ Nash based his friction law on a study of several classical approaches including Ludwieg–Tillman. Nash's approach involves the following three equations:

\begin{matrix} \frac{u_{e}}{u_{τ}} = \frac{1}{0.400} ln (R_{δ_{1}}) + 3.7 + K' (G), \\ K' (G) = \frac{3}{2} G + \frac{2110}{G^{2} + 200} - 18.5, \\ G = (\frac{u_{e}}{u_{τ}}) \frac{H - 1}{H} . \end{matrix}

(32)

Nash's law agrees quite closely with the UVP over a wide range of Reynolds numbers.

Also included in Fig. 10 are two versions of what Nagib et al.³⁵ call the Coles–Fernholz law.

\frac{u_{e}}{u_{τ}} = \frac{1}{k} ln (R_{δ_{2}}) + C .

(33)

Finally, Polewski and Cizmas³⁶ carried out a detailed comparison of several turbulence models. In their Fig. 9(a), the skin friction coefficient on a flat plate generated by three codes using the Spalart–Allmaras³⁷ model is plotted up to $R_{x} = 10^{7}$ ⁠. The curves agree very closely. The Spalart–Allmaras friction law at this Reynolds number deduced from this plot is

\frac{u_{e}}{u_{τ}} = \frac{1}{0.4292} ln (R_{τ}) + 8.9414.

(34)

Note that for the comparisons in Figs. 10 and 11, the $R_{δ_{1}}$ and $R_{δ_{2}}$ that appear in Eqs. (31)–(33) are related to $R_{τ}$ using the integrals (5) and (6) expressed in terms of $R_{τ}$ ⁠, Eqs. (26) and (27). See also the explicit high Reynolds number relations, Eqs. (40) and (43), in Cantwell.¹⁰

FIG. 11.

View large Download slide

Four friction laws compared to the UVP (blue line) and plotted as C_f in linear coordinates. The magenta line is the friction law developed by Nash.³⁴ The red and black lines are the two versions of the Coles–Fernholz law. The brown line is the Spalart–Allmaras model deduced from data in Polewski and Cizmas.³⁶

The closest to the UVP law in Fig. 10 are the Coles–Fernholz law with Coles' recommended values $(k, C) = (0.41, 5.0)$ and the Spalart–Allmaras relation. Outside of the range $1000 < R_{τ} < 10 000$ ⁠, the Ludwieg–Tillman law leads to unrealistically high $u_{e} / u_{τ}$ (low C_f) values and will not be considered further.

Figure 11 shows the UVP, Nash, Spalart–Allmaras, and the two Coles–Fernholz laws in linear coordinates over the range $0 < R_{τ} < 25 000$ which covers the vast bulk of available boundary layer data (⁠ $0 < R_{δ_{1}} < 1.009 \times 10^{5}$ corresponding to $0 < R_{δ_{2}} < 7.93 \times 10^{4}$ corresponding to $0 < R_{x} < 7.38 \times 10^{7}$ ⁠). All five curves in Fig. 10 agree within about 5% over the range above $R_{τ} > 1000$ ⁠. The closest agreement with the UVP law in this Reynolds number range are the laws due to Nash³⁴ and Spalart and Allmaras.³⁷ Note that only the UVP merges with the laminar limit, $R_{τ} \to 0$ as shown in Fig. 5.

Nagib et al.³⁵ carried out measurements of the skin friction on a flat plate in the Reynolds number range $10 000 < R_{δ_{2}} < 70 000$ ⁠. They used that data, along with older data corrected to form a common experimental dataset, in a comprehensive study of several widely used classical friction laws. The result is shown in Fig. 12 where their Fig. 7 is overlaid by the friction law generated by the UVP (blue circles) and the friction law of Nash³⁴ (magenta dashed line). Note that the $R_{τ}$ and $R_{δ_{2}}$ domains in Figs. 11 and 12 represent almost identical Reynolds number ranges.

FIG. 12.

View large Download slide

Friction coefficient from the UVP compared with Fig. 7 from Nagib et al.³⁵ The friction law of Nash³⁴ is also included (magenta dashed line).

Figures 10–12 show that the UVP friction law is consistent with the vast body of friction data in the literature.

III. THE VON KÁRMÁN INTEGRAL EQUATION EXPRESSED IN TERMS OF $R_{τ}$

Define the overall flow Reynolds number in terms of the free stream velocity $u_{\infty}$ and an appropriate body length scale which we will take to be the airfoil leading edge radius of curvature r shown in Fig. 2

R_{e} = \frac{u_{\infty} r}{ν} .

(35)

The main drag results are expressed in terms of the chord Reynolds number, R_chord, obtained by multiplying Eq. (35) by $c / r = 91.3737$ for the J0012 airfoil or $c / r = 63.5870$ for the NACA0012 airfoil. Introduce the length scale r and let $ξ = x / r$ ⁠. The dimensionless potential flow velocity at the surface of the body is

\frac{u_{e}}{u_{\infty}} = U (ξ) .

(36)

Now re-express the Kármán equation (3) in terms of $R_{τ}$ ⁠. First, express the momentum thickness, in terms of the momentum thickness Reynolds number (27),

δ_{2} = \frac{ν}{u_{e}} \int_{0}^{R_{τ}} u^{+} (1 - \frac{u_{τ}}{u_{e}} u^{+}) d y^{+} = \frac{ν}{u_{e}} F_{2} (R_{τ}) .

(37)

Differentiate (37) with respect to x

\begin{matrix} \frac{d δ_{2}}{d x} = - \frac{ν}{u_{e}^{2}} \frac{d u_{e}}{d x} F_{2} (R_{τ}) + \frac{ν}{u_{e}} \frac{d F_{2} (R_{τ})}{d R_{τ}} \frac{d R_{τ}}{d x} \\ = - \frac{δ_{2}}{u_{e}} \frac{d u_{e}}{d x} + \frac{ν}{u_{e}} \frac{d F_{2} (R_{τ})}{d R_{τ}} \frac{d R_{τ}}{d x} . \end{matrix}

(38)

Using (38) in (3), the Kármán equation is now written

\frac{ν}{u_{e}} \frac{d F_{2}}{d R_{τ}} \frac{d R_{τ}}{d x} + (δ_{2} + δ_{1}) \frac{1}{u_{e}} \frac{d u_{e}}{d x} - {(\frac{u_{τ}}{u_{e}})}^{2} = 0.

(39)

Normalize x by the nose radius, r, and use (26) and (27) to express δ₁ and δ₂ in (39) in terms of the displacement and momentum thickness Reynolds numbers. Use the definition (25) for the friction law and the definition (36) for the inviscid velocity.

\frac{d F_{2}}{d R_{τ}} \frac{d R_{τ}}{d ξ} + (F_{2} + F_{1}) \frac{1}{U} \frac{d U}{d ξ} - \frac{u_{\infty} r}{ν} \frac{U}{F_{0}^{2}} = 0.

(40)

The form of the Kármán equation that we need to solve numerically is

\frac{d R_{τ}}{d ξ} = \frac{U R_{e}}{F_{0}^{2} F_{3}} (1 - \frac{F_{0}^{2}}{R_{e}} (F_{2} + F_{1}) \frac{1}{U^{2}} \frac{d U}{d ξ}),

(41)

where (28) has been used. The zero pressure gradient form of (41) was used by Cantwell¹⁰ to determine integral measures of the ZPG boundary layer over $0 \leq R_{τ} < \infty$ ⁠.

IV. THE EFFECT OF PRESSURE GRADIENT ON THE UVP

Subrahmanyam et al.¹⁷ used the UVP to approximate the adverse pressure gradient data of Perry and Marusic²⁷ for the two upstream velocity cases, $u_{\infty} = 10$ and $u_{\infty} = 30$ m/s. The near wall parameters $(k, a, m)$ did not vary significantly with the Reynolds number, whereas the outer flow wake parameters (b, n) varied quite a lot. The results from Subrahmanyam et al.¹⁷ are shown in Table II along with the fit of the UVP to the Ref. 27 data in Figs. 13 and 14.

FIG. 13.

View large Download slide

Reprinted with permission from Subrahmanyam et al., J. Fluid Mech. 933, A16 (2022). Copyright 2022 Author(s) licensed under a Creative Commons Attribution (CC BY) License.¹⁷ Comparison between the universal velocity profile and the $u_{\infty} = 10$ m/s adverse pressure gradient data of Perry and Marusic²⁷ (open red circles).

FIG. 14.

View large Download slide

Reprinted with permission from Subrahmanyam et al. J. Fluid Mech. 933, A16 (2022). Copyright 2022 Author(s) licensed under a Creative Commons Attribution (CC BY) License.¹⁷ Comparison between the universal velocity profile and the $u_{\infty} = 30$ m/s adverse pressure gradient data of Perry and Marusic²⁷ (open red circles).

The sink flow data from Jones et al.²⁸ will be used to determine how the wake parameters of the UVP, (b, n) vary for favorable pressure gradients. These results are provided in Table III and the comparison between the UVP and the data is shown in Figs. 15–17. The quantity K in Table III is essentially the inverse of the sink flow Reynolds number,

K = \frac{ν}{u_{e}^{2}} \frac{d u_{e}}{d x} = \frac{ν}{u_{\infty} L},

(42)

where $u_{\infty}$ is the converging channel entry velocity, x is the coordinate along the lower wall and L is the effective length of the sink flow channel. The optimal value of the UVP parameter k for these data tends to be relatively high compared to ZPG and adverse pressure gradient data, particularly near the entrance to the converging channel and at the lowest Reynolds number. This may be a consequence of the relatively low values of $R_{τ}$ that characterize the data. Comparable values of k were observed at low $R_{τ}$ in the pipe data (Cantwell¹¹ Fig. 9, PSP cases 1–5). The mean parameter values in Table I for pipe flow were averaged over PSP cases 6–26.

TABLE III.

Run data, Reynolds number, optimal model parameters, and RMS error for favorable pressure gradient boundary layer datasets from Jones et al.²⁸ Converging channel entry velocities are $u_{0} = 5.0, 7.5$ and 10.0 m/s. Channel half height thickness for these data are at $u = 0.996 u_{e}$ ⁠; $δ_{h} = δ_{0.996}$ ⁠. Kinematic viscosity is $ν = 1.51 \times 10^{- 5} m / s^{2}$ ⁠. Note, according to Jones et al.,²⁸ $u_{e} / u_{0} = 1 / (1 - x / L)$ and the calibrated sink length is $L = 5.60 m$ ⁠.

$K \times 10^{7}$	x (mm)	$u_{0} (\frac{m}{s})$	$R_{δ_{1}}$	$R_{δ_{2}}$	$- β$	$- β_{c}$	$R_{δ_{996}}$	$R_{τ}$	$(\frac{u_{e}}{u_{τ}})$	k	a	m	b	n	$u^{+}_{rms}$
5.39	800	5.0	1112	780	0.2436	0.4145	8656	429	20.16	0.4686	26.88	1.2136	0.3931	1.3295	0.098
5.39	1600	5.0	1629	1192	0.3780	0.6548	14 130	681	20.75	0.4630	25.58	1.2408	0.5284	1.1644	0.068
5.39	2200	5.0	1648	1209	0.3880	0.6726	14 806	709	20.90	0.4670	27.89	1.0869	0.5806	1.1343	0.105
5.39	2800	5.0	1946	1449	0.4713	0.8222	18 244	861	21.20	0.4555	26.89	1.1196	0.6934	1.0908	0.089
5.39	3280	5.0	2138	1606	0.5276	0.9241	20 549	960	21.40	0.4445	26.04	1.1563	0.6831	1.1818	0.067
5.39	3580	5.0	2226	1687	0.5500	0.9671	22 468	1049	21.41	0.4448	25.96	1.1920	0.8189	1.1514	0.068
3.59	800	7.5	1496	1069	0.3613	0.4125	11 752	555	21.17	0.4602	26.62	1.1906	0.3738	1.2986	0.092
3.59	1600	7.5	2000	1470	0.3382	0.5868	17 325	798	21.70	0.4644	28.91	1.0413	0.4647	1.1745	0.085
3.59	2200	7.5	2350	1755	0.4030	0.7040	21 366	977	21.86	0.4557	27.30	1.0822	0.5752	1.0988	0.060
3.59	2800	7.5	2827	2155	0.4850	0.8748	27 899	1262	22.12	0.4431	26.19	1.1527	0.7463	1.0917	0.065
3.59	3280	7.5	2928	2245	0.5202	0.9191	29 762	1338	22.25	0.4374	25.93	1.1726	0.7271	1.1931	0.055
3.59	3580	7.5	3027	2339	0.5344	0.9473	32 204	1452	22.17	0.4403	25.78	1.2037	0.8700	1.1579	0.062
2.70	800	10.0	1862	1343	0.3161	0.4092	15 001	690	21.74	0.4564	28.35	1.0045	0.4324	1.1316	0.097
2.70	1600	10.0	2555	1904	0.3421	0.5969	22 546	1013	22.27	0.4529	27.56	1.0538	0.5090	1.1120	0.074
2.70	2200	10.0	2873	2160	0.3908	0.6846	26 114	1164	22.44	0.4525	26.90	1.0669	0.5491	1.0864	0.047
2.70	2800	10.0	3372	2577	0.4689	0.8273	32 769	1444	22.70	0.4402	25.79	1.1307	0.6833	1.0766	0.047
2.70	3280	10.0	3725	2851	0.5297	0.9351	35 899	1564	22.95	0.4364	25.45	1.1227	0.6653	1.0704	0.073
2.70	3580	10.0	3936	3044	0.5575	0.9888	39 990	1746	22.91	0.4368	25.31	1.1607	0.7504	1.0923	0.095

$K \times 10^{7}$	x (mm)	$u_{0} (\frac{m}{s})$	$R_{δ_{1}}$	$R_{δ_{2}}$	$- β$	$- β_{c}$	$R_{δ_{996}}$	$R_{τ}$	$(\frac{u_{e}}{u_{τ}})$	k	a	m	b	n	$u^{+}_{rms}$
5.39	800	5.0	1112	780	0.2436	0.4145	8656	429	20.16	0.4686	26.88	1.2136	0.3931	1.3295	0.098
5.39	1600	5.0	1629	1192	0.3780	0.6548	14 130	681	20.75	0.4630	25.58	1.2408	0.5284	1.1644	0.068
5.39	2200	5.0	1648	1209	0.3880	0.6726	14 806	709	20.90	0.4670	27.89	1.0869	0.5806	1.1343	0.105
5.39	2800	5.0	1946	1449	0.4713	0.8222	18 244	861	21.20	0.4555	26.89	1.1196	0.6934	1.0908	0.089
5.39	3280	5.0	2138	1606	0.5276	0.9241	20 549	960	21.40	0.4445	26.04	1.1563	0.6831	1.1818	0.067
5.39	3580	5.0	2226	1687	0.5500	0.9671	22 468	1049	21.41	0.4448	25.96	1.1920	0.8189	1.1514	0.068
3.59	800	7.5	1496	1069	0.3613	0.4125	11 752	555	21.17	0.4602	26.62	1.1906	0.3738	1.2986	0.092
3.59	1600	7.5	2000	1470	0.3382	0.5868	17 325	798	21.70	0.4644	28.91	1.0413	0.4647	1.1745	0.085
3.59	2200	7.5	2350	1755	0.4030	0.7040	21 366	977	21.86	0.4557	27.30	1.0822	0.5752	1.0988	0.060
3.59	2800	7.5	2827	2155	0.4850	0.8748	27 899	1262	22.12	0.4431	26.19	1.1527	0.7463	1.0917	0.065
3.59	3280	7.5	2928	2245	0.5202	0.9191	29 762	1338	22.25	0.4374	25.93	1.1726	0.7271	1.1931	0.055
3.59	3580	7.5	3027	2339	0.5344	0.9473	32 204	1452	22.17	0.4403	25.78	1.2037	0.8700	1.1579	0.062
2.70	800	10.0	1862	1343	0.3161	0.4092	15 001	690	21.74	0.4564	28.35	1.0045	0.4324	1.1316	0.097
2.70	1600	10.0	2555	1904	0.3421	0.5969	22 546	1013	22.27	0.4529	27.56	1.0538	0.5090	1.1120	0.074
2.70	2200	10.0	2873	2160	0.3908	0.6846	26 114	1164	22.44	0.4525	26.90	1.0669	0.5491	1.0864	0.047
2.70	2800	10.0	3372	2577	0.4689	0.8273	32 769	1444	22.70	0.4402	25.79	1.1307	0.6833	1.0766	0.047
2.70	3280	10.0	3725	2851	0.5297	0.9351	35 899	1564	22.95	0.4364	25.45	1.1227	0.6653	1.0704	0.073
2.70	3580	10.0	3936	3044	0.5575	0.9888	39 990	1746	22.91	0.4368	25.31	1.1607	0.7504	1.0923	0.095

FIG. 15.

View large Download slide

Comparison between the universal velocity profile and the $K = 5.39 \times 10^{- 7}$ sink flow data of Jones et al.²⁸ (open red circles).

FIG. 16.

View large Download slide

Comparison between the universal velocity profile and the $K = 3.59 \times 10^{- 7}$ sink flow data of Jones et al.²⁸ (open red circles).

FIG. 17.

View large Download slide

Comparison between the universal velocity profile and the $K = 2.70 \times 10^{- 7}$ sink flow data of Jones et al.²⁸ (open red circles).

A. A modified Clauser parameter

The pressure gradient term in (41) can be arrived at from another direction. The strength of the pressure gradient in a turbulent boundary layer is often expressed in terms of the Clauser³⁸ parameter $β = (δ_{1} / τ_{w}) (d p / d x)$ ⁠. See, for example, Perry and Marusic.²⁷ Here, we will introduce a modified Clauser parameter

β_{c} = ((δ_{1} + δ_{2}) / τ_{w}) (d p_{e} / d x) .

(43)

The Bernoulli equation is used to express β_c in terms of the local free stream velocity, and the thicknesses are expressed in terms of the local thickness Reynolds numbers. The result is identical to the pressure gradient term in the Kármán equation (41)

β_{c} = - \frac{F_{0}^{2}}{R_{e}} (F_{2} + F_{1}) \frac{1}{U^{2}} \frac{d U}{d ξ} .

(44)

The UVP parameters n and b are plotted in Figs. 18 and 19 against β_c. The solid line

n (β_{c}) = 1.419 350 + 0.271 499 β_{c}

(45)

approximates the data for $n (β_{c})$ in Fig. 18 reasonably well with root mean square (rms) error 0.197. Two approximations to the data for $b (β_{c})$ are shown in Fig. 19. The dashed line is the simpler of the two and is given by

b (β_{c}) = 0.018 193 8 + \frac{0.286 852}{1 + 0.654 161 β_{c}}

(46)

with rms error 0.0440. Equation (46) gives a good approximation to most of the data but is not accurate near $β_{c} = 0$ ⁠. This is a problem because over most of the airfoil surface the pressure gradient is quite small as can be seen in the lower plot in Fig. 22. This means that β_c is also quite small over a substantial part of the airfoil surface. The solid line is given by

b (β_{c}) = 0.018 193 8 + \frac{0.286 852}{1 + 0.654 161 β_{c}} - 0.14 \frac{e^{- 2 β_{c}^{2}}}{{(2.2 + β_{c})}^{2 / 3}}

(47)

with rms error 0.052 5. The exponential term provides the correction needed to better approximate the data in the neighborhood of $β_{c} = 0$ without significantly increasing the error away from $β_{c} = 0$ ⁠.

FIG. 18.

View large Download slide

Correlation between the model parameter n and β_c for the data in Tables II and III.

FIG. 19.

View large Download slide

Correlation between the model parameter b and β_c for the data in Tables II and III. The solid line provides a more accurate approximation near $β_{c} = 0$ ⁠.

The effect of this change can be seen by evaluating C_f at $β_{c} = 0$ for $R_{τ} = 10^{4}$ ⁠. If the parameters b and n from Table I, $(b, n) = (0.1752, 2.1707)$ ⁠, are used $C_{f} = 0.002 13$ ⁠. If Eqs. (45) and (46) are used, $(b, n) = (0.3050, 1.4194)$ and $C_{f} = 0.002 38$ ⁠, which is $12 %$ higher than the Table I value and will give almost a $12 %$ higher airfoil drag. If Eqs. (45) and (47) are used, $(b, n) = (0.2223, 1.4194)$ and the friction coefficient is $C_{f} = 0.002 15$ which is less than $1 %$ higher than the Table I value and within the $7.5 %$ error in C_f generated when the optimal parameter values are increased or decreased by one standard deviation. Given the relatively limited data near $β_{c} = 0$ ⁠, the best approximation to the friction is generated by Eqs. (45) and (47). More data for weak pressure gradient conditions is needed.

In the presence of a pressure gradient, the wall parameters $(k, a, m)$ will be kept constant at the ZPG boundary layer values in Table I while the wake parameters $(b, n)$ will be treated as functions of β_c using the correlations (45) and (47). This will be the approach used for the drag calculations presented in Secs. VI–VIII.

The important conclusion of this is that the length scale function (21) no longer depends on the parameters b and n independently, but depends only on the pressure gradient parameter β_c.

λ (k, a, m, b, n, y^{+}) \to λ (k, a, m, β_{c}, y^{+})

(48)

and

λ (y^{+}) = \frac{k y^{+} (1 - e^{- {(\frac{y^{+}}{a})}^{m}})}{{(1 + {(\frac{y^{+}}{b (β_{c}) R_{τ}})}^{n (β_{c})})}^{1 / n (β_{c})}} .

(49)

It should be noted that the number of significant figures provided in this paper as well as in previous work on the UVP^10,11,17 is intended to allow an interested reader to reproduce our results with accurate comparisons. It is not intended to reflect experimental or computed accuracy of the data.

V. LEADING EDGE TREATMENT

In Sec. VI, the J0012 airfoil boundary layer calculation will be considered in detail but first it is necessary to analyze the flow at the airfoil leading edge. To model the leading edge flow depicted in Fig. 2, we will use the inviscid solution about an elliptic wing presented in Van Dyke.³⁹ Let

α \equiv 1 + {(2 \frac{r}{c})}^{1 / 2},

(50)

where r is the leading edge radius of curvature, and c is the chord of the elliptic wing. The coordinate of the upper surface of the ellipse with the leading edge located at $\tilde{x} / c = 0$ is

\frac{\tilde{y}}{c} = \frac{1}{2} (α - 1) {(1 - {(2 \frac{\tilde{x}}{c} - 1)}^{2})}^{1 / 2}

(51)

and the surface velocity is

U (\frac{\tilde{x}}{c}) = \frac{α}{{(1 + {(α - 1)}^{2} \frac{{(2 \frac{\tilde{x}}{c} - 1)}^{2}}{1 - {(2 \frac{\tilde{x}}{c} - 1)}^{2}})}^{1 / 2}} .

(52)

The coordinate $x (\tilde{x})$ along the wing surface is

\frac{x}{c} = \frac{1}{2} (E (\frac{π}{2}, 1 - α) - E (Arcsin (1 - 2 \frac{\tilde{x}}{c}), 1 - α)),

(53)

where the function $E (θ, 1 - α)$ in Eq. (53) is the incomplete elliptic integral of the second kind.

Using Eqs. (51)–(53), the following limits are found to hold near the leading edge. Recall $ξ = x / r$

\begin{matrix} \lim_{\frac{\tilde{x}}{c} \to 0} \frac{\tilde{y}}{c} = \frac{x}{c}, \\ \lim_{\frac{\tilde{x}}{c} \to 0} \frac{x}{c} = (α - 1) {(\frac{\tilde{x}}{c})}^{1 / 2}, \\ \lim_{\frac{\tilde{x}}{c} \to 0} U = \frac{2 α}{α - 1} {(\frac{\tilde{x}}{c})}^{1 / 2}, \\ \lim_{\frac{\tilde{x}}{c} \to 0} U = \frac{2 α}{{(α - 1)}^{2}} (\frac{x}{c}), \\ \lim_{\frac{\tilde{x}}{c} \to 0} U = α ξ . \end{matrix}

(54)

A. The Hiemenz solution near the leading edge

With the local inviscid solution (54) known, the viscous laminar solution near the leading edge can be generated from the classical Hiemenz solution at a stagnation point,

\begin{matrix} η = α^{1 / 2} {(\frac{u_{\infty} r}{ν})}^{1 / 2} (\frac{y}{r}) \\ \frac{v (x, y)}{u_{\infty}} = - α^{1 / 2} (\frac{ν}{u_{\infty} r}) f (η) \\ \frac{u (x, y)}{u_{\infty}} = α ξ f' (η) \\ \frac{Ω (x, y) r}{u_{\infty}} = - \frac{r}{u_{\infty}} \frac{\partial u}{\partial y} = - α^{3 / 2} {(\frac{u_{\infty} r}{ν})}^{1 / 2} ξ f ″ (η) . \end{matrix}

(55)

Substitute the variables in (55) into the vorticity equation

\begin{matrix} u \frac{\partial Ω}{\partial x} + v \frac{\partial Ω}{\partial y} - ν \frac{\partial^{2} Ω}{\partial x^{2}} - ν \frac{\partial^{2} Ω}{\partial y^{2}} = 0. \end{matrix}

(56)

The result is

α^{1 / 2} {(ν r u_{\infty})}^{1 / 2} (f^{″ ″} (η) + f (η) f ″' (η) - f' (η) f ″ (η)) = 0.

(57)

The fourth-order ordinary differential equation (ODE) in (57) can be integrated once leading to the third-order equation governing the Hiemenz solution.

\begin{matrix} f^{‴} (η) + f (η) f ″ (η) + 1 - f' {(η)}^{2} = 0, \\ f (0) = 0, f' (0) = 0, f' (\infty) = 1. \end{matrix}

(58)

The solution of (58) is shown in Fig. 20. Near the leading edge, the Hiemenz solution leads to the following limits:

\begin{matrix} \lim_{ξ \to 0} {(\frac{u_{τ}}{u_{e}})}^{2} = \frac{f ″ (0)}{α^{1 / 2} R_{e}^{1 / 2} ξ}, \\ \lim_{ξ \to 0} \frac{δ_{1} u_{\infty}}{ν} = \frac{c_{1}}{α^{1 / 2} R_{e}^{1 / 2}}, \\ \lim_{ξ \to 0} \frac{δ_{2} u_{\infty}}{ν} = \frac{c_{2}}{α^{1 / 2} R_{e}^{1 / 2}}, \\ \lim_{ξ \to 0} R_{τ} = η_{h} f ″ {(0)}^{1 / 2} α^{1 / 4} R_{e}^{1 / 4} ξ^{1 / 2}, \end{matrix}

(59)

where $f ″ (0) = 1.232 59, c_{1} = 0.647 836$ ⁠, and $c_{2} = 0.292 282$ ⁠. The shape factor of the Hiemenz flow at the forward stagnation point is $H = δ_{1} / δ_{2} = 2.216$ ⁠, which can be compared to $H = 2.5$ for channel flow and $H = 2.604$ for the ZPG Blasius boundary layer.

FIG. 20.

View large Download slide

Hiemenz streamwise velocity and vorticity.

B. The UVP near the leading edge

In the UVP integral boundary layer method presented in this paper, the simplest approach will be used where the UVP is assumed to apply all the way from the forward stagnation point of the airfoil. Therefore, it is important to consider the degree of error associated with this assumption by comparing the UVP with the “exact” Hiemenz solution.

The UVP approaches the laminar profile, (20), at the wing leading edge. From this profile, the following limits can be determined:

\begin{matrix} \lim_{R_{τ} \to 0} F_{0} = \frac{R_{τ}}{2}, \\ \lim_{R_{τ} \to 0} F_{1} = \frac{R_{τ}^{2}}{6}, \\ \lim_{R_{τ} \to 0} F_{2} = \frac{R_{τ}^{2}}{15}, \\ \lim_{R_{τ} \to 0} F_{3} = \frac{2 R_{τ}}{15} . \end{matrix}

(60)

Substitute the limits in (60) into the Kármán equation (41). The result is

\frac{d}{d ξ} (\frac{R_{τ}^{4}}{120}) = R_{e} U - (\frac{R_{τ}^{4}}{120}) \frac{1}{U} \frac{d U}{d ξ} .

(61)

Equation (61) is a first-order ODE with the following integral:

\lim_{R_{τ} \to 0} R_{τ} = {(\frac{120 R_{e}}{U^{7}} \int_{0}^{ξ} U^{8} d ξ)}^{\frac{1}{4}},

(62)

where the initial condition $U (0) = 0$ has been used. Substitute the last relation in (54), $U = α ξ$ ⁠, into (62). The result is

\lim_{R_{τ} \to 0} R_{τ} = {(\frac{40}{3})}^{1 / 4} α^{1 / 4} R_{e}^{1 / 4} ξ^{1 / 2} .

(63)

Compare the UVP limit (63) to the Hiemenz limit in Eq. (59). In the last relation in (59), let $η_{h} = 2.381$ corresponding to $u / u_{e} = 0.99$ ⁠,

\lim_{R_{τ} \to 0} \frac{R_{τ}_{UVP}}{R_{τ}_{Hiemenz}} = \frac{{(\frac{40}{3})}^{1 / 4}}{2.381 f ″ {(0)}^{1 / 2}} = 0.723 .

(64)

The laminar limiting profile of the UVP, Eq. (20), can be expressed in terms of η by noting that $y / δ_{h} = η / η_{h}$

\frac{u (x, y)}{u_{\infty}} = α ξ (\frac{2 η}{η_{h}}) (1 - \frac{1}{2} \frac{η}{η_{h}}) .

(65)

The Hiemenz and UVP velocity profiles are compared in Fig. 21.

FIG. 21.

View large Download slide

Comparison between the Hiemenz velocity profile (black) and UVP (blue) near the forward stagnation point.

The conclusion from this discussion is that near the leading edge, the UVP limit (63) has the same dependence on α, R_e, and ξ as the “exact” Hiemenz limit in (59) but the magnitude of $R_{τ}$ ⁠, and therefore, the friction is about 28% below the Hiemenz limit. Keeping this in mind, the UVP will be assumed to apply beginning at the leading edge in the airfoil examples discussed later.

This is not the only approach that could be used. The laminar boundary layer equations could be solved beginning at the leading edge up to a transition region modeled with an appropriate transition criterion. This could then be followed by a turbulent boundary layer approximated by the UVP. For airfoil flows with a substantial laminar region, this would be more accurate but would require more information on free stream turbulence, and other factors that could influence transition. The viscous drag coefficient of low Reynolds number airfoil flows would be roughly $10 %$ lower using this approach compared to using the UVP over the whole airfoil surface. The reason is that, even though the UVP friction might be lower very close to the leading edge as just discussed, the friction associated with the UVP channel flow limit is $10 %$ or more higher than the Blasius-like friction of the rest of the laminar part of the flow,¹⁰ depending on the pressure gradient.

In this paper our focus is on high Reynolds number and so we will use the UVP all the way to the leading edge and comparisons will be with tripped airfoil data.

VI. ITERATIVE PROCEDURE AND SOLUTION FOR THE J0012 AIRFOIL AT $R_{chord} = 10^{7}$ USING THE INTEGRAL FORM OF THE UVP

The J0012 airfoil geometry and inviscid surface velocity is shown in Fig. 22. Figure 23 shows the relation between the surface and chordwise coordinates including an inset showing the leading edge. The infinite slope of the wing near the leading edge and the singularity in the inviscid velocity derivative are both alleviated when the flow is expressed in terms of the surface coordinate $ξ = x / r$ rather than the chordwise coordinate $\tilde{ξ} = \tilde{x} / r$ ⁠. The potential flow distributions of $U (\tilde{ξ})$ for the J0012 and NACA0012 airfoils are shown in Fig. 24.

FIG. 22.

View large Download slide

Joukowsky 0012 airfoil with free stream velocity and velocity derivative. Maximum surface velocity (⁠ $U = 1.158 9$ ⁠) is at $\tilde{x} / c = 0.174 45$ ⁠. Zero airfoil slope is at $\tilde{x} / c = 0.368 47$ ⁠.

FIG. 23.

View large Download slide

Joukowsky 0012 surface coordinate vs chordwise coordinate. Inset shows shows the relation near the leading edge. Leading edge radius of curvature is $r / c = 0.010 944$ ⁠. Airfoil trailing edge is at $\tilde{x_{c}} / r = 91.3737$ corresponding to $x_{c} / r = 92.9941$ ⁠.

FIG. 24.

View large Download slide

Potential flow velocity distributions about the J0012 (blue) and NACA0012 (magenta) airfoils.

Figure 25 shows eight successive solutions of the Kármán equation (41) over the J0012 airfoil at an airfoil chord Reynolds number, $R_{chord} = 10^{7}$ ⁠. The integral form of the UVP, Eq. (19), is used in this example. In Sec. VII, an alternative approach using the high Reynolds number explicit form of the UVP (23) will be described. The calculations were carried out in Mathematica^TM on a 2017 MacBook Pro (OS 10.12.6, 2.9 GHz Intel core i7 processor). The iterative procedure converges relatively rapidly.

FIG. 25.

View large Download slide

Successive iterations of $R_{τ} (x / r)$ for the J0012 airfoil. Inset shows the covergence of the trailing edge value of $R_{τ_{T E}}$ with iteration number. The final iteration (blue) has a trailing edge value of $R_{τ_{T E}} = 4673$ ⁠.

A. Step 1—Initial data, initial calculation of $R_{δ_{1}} (R_{τ})$ and $R_{δ_{2}} (R_{τ})$ ⁠, first integration of the Kármán equation

Initial data required by the method begins with the kinematic viscosity ν and free stream velocity, $u_{\infty}$ ⁠. The wing profile is defined, particularly the chord c and the leading edge radius r, which are used to determine α in Eqs. (50) and (63). The Reynolds number is $R_{e} = u_{\infty} r / ν$ ⁠. The wing profile is specified as

\tilde{y} / c = f (\tilde{x} / c) .

(66)

From the wing profile is calculated the wing surface coordinate x used as the independent variable in the UVP integral method,

\frac{x}{c} = \int_{0}^{\frac{\tilde{x}}{c}} {(1 + f^{'} {(\frac{\tilde{x}}{c})}^{2})}^{1 / 2} d (\frac{\tilde{x}}{c}) .

(67)

See Figs. 23 and 24, where x and $\tilde{x}$ are normalized by the wing leading edge radius by multiplication by c/r. The velocity distribution $U (ξ)$ ⁠, along with $d U / d ξ$ are determined from a potential solution about the body in question. For the UVP integral method to work effectively, it is important that the free stream velocity and its derivative are accurately known and smooth.

Prior to integrating the Kármán equation or focusing on the airfoil, the integral form of the UVP, Eq. (19), is used to determine the functions $F_{0} (R_{τ}), F_{1} (R_{τ}), F_{2} (R_{τ})$ ⁠, and $F_{3} (R_{τ})$ over the $R_{τ}$ range of interest. For the calculations presented here $0 \leq R_{τ} < 10^{10}$ ⁠. This is done with the constants $(k, a, b, m, n)$ fixed at the values given in Table I, i.e., for $β_{c} = 0$ ⁠. Since nested integration is required to produce the boundary layer thicknesses, this step can be quite slow when the range of interest includes large values of $R_{τ}$ ⁠. However, these constant-b-constant-n thicknesses only need to be determined once and stored for use in the first iteration each time the drag for a new R_chord is to be computed.

In a bootstrap procedure, the Kármán equation (68) is solved numerically for the initial distribution $R_{τ_{1}} (ξ)$ using the Mathematica^TM function NDSolve,

\begin{matrix} \frac{d R_{τ_{1}} (ξ)}{d ξ} = \frac{U (ξ) R_{e}}{F_{0} {(R_{τ_{1}} (ξ))}^{2} F_{3} (R_{τ_{1}} (ξ))} \\ \times (1 - \frac{F_{0} {(R_{τ_{1}} (ξ))}^{2}}{R_{e}} (F_{2} (R_{τ_{1}} (ξ)) + F_{1} (R_{τ_{1}} (ξ))) \frac{1}{U {(ξ)}^{2}} \frac{d U}{d ξ}) . \end{matrix}

(68)

The initial value of ξ is chosen to be close to the airfoil leading edge. For the J0012 case, with $U (ξ)$ known analytically, ξ_initial was chosen to be $10^{- 7}$ ⁠. For the NACA 0012 case discussed later, with $U (ξ)$ known approximately, ξ_initial was chosen to be 0.1. The initial distribution, $R_{τ_{1}} (ξ)$ ⁠, is the black curve in Fig. 25 with $R_{τ_{1}_{T E}} = 9784$ ⁠.

The thickness functions $F_{0} (R_{τ_{1}}), F_{1} (R_{τ_{1}}), F_{2} (R_{τ_{1}})$ ⁠, and $F_{3} (R_{τ_{1}})$ are computed over the range $0 < R_{τ_{1}} < R_{τ_{1}_{T E}}$ and are used to calculate the first nonzero distribution of the modified Clauser parameter, $β_{c} (R_{τ_{1}})$ ⁠. After the first iteration, the wall parameters $(k, a, m)$ are kept constant at the ZPG boundary layer values in Table I while $n (β_{c})$ and $b (β_{c})$ are allowed to vary over the airfoil surface according to the correlations (45) and (47).

B. Step 2

The monotonically increasing friction Reynolds number distribution over the airfoil computed in the first iteration, $R_{τ_{1}} (ξ)$ together with $β_{c} (R_{τ_{1}})$ are then used to prepare for iteration 2. New thickness functions $F_{0_{1}} (R_{τ_{1}}), F_{1_{1}} (R_{τ_{1}}), F_{2_{1}} (R_{τ_{1}})$ ⁠, and $F_{3_{1}} (R_{τ_{1}})$ are computed over the range $0 < R_{τ_{1}} < R_{τ_{1}_{T E}}$ ⁠. The Kármán equation (69)

\begin{matrix} \frac{d R_{τ_{2}} (ξ)}{d ξ} = \frac{U (ξ) R_{e}}{F_{0_{1}} {(R_{τ_{2}} (ξ))}^{2} F_{3_{1}} (R_{τ_{2}} (ξ))} \\ \times (1 - \frac{F_{0_{1}} {(R_{τ_{2}} (ξ))}^{2}}{R_{e}} (F_{2_{1}} (R_{τ_{2}} (ξ)) + F_{1_{1}} (R_{τ_{2}} (ξ))) \frac{1}{U {(ξ)}^{2}} \frac{d U}{d ξ}) \end{matrix}

(69)

is integrated to generate $R_{τ_{2}} (ξ)$ ⁠. This is the red curve in Fig. 25.

C. Step 3—Finish

The process is repeated until the nth iterate when further changes in $R_{τ_{n}_{T E}}$ are sufficiently small (⁠ $< 2 %$ ⁠). This is the blue curve in Fig. 25. The iterative procedure is shown schematically in Appendix B.

The time consuming steps in this procedure are the computations of the thickness functions, F₁ and F₂ which require nested integrations at npoints over the airfoil surface where npoints may vary from 500 for $R_{chord} = 10^{5}$ up to 2500 for $R_{chord} = 10^{12}$ ⁠. At $R_{chord} = 10^{12}$ each iteration may take 3–4 h on the laptop described above. In this example, eight iterations are needed to reach convergence with $R_{τ_{1}_{T E}}$ changing by less than $2 %$ over the last two iterations. The corresponding change in $C_{d_{v}}$ over the last two iterations is $0.3 %$ ⁠. The inset shows the convergence of $R_{τ_{i}_{T E}} (\frac{c}{r})$ at the end of which, $R_{τ_{8}_{T E}} = 4673$ ⁠. At this point, all of the boundary layer characteristics are known along with the velocity profile at any point.

The procedure could be repeated with the airfoil surface modified by the displacement thickness however the close correspondence between the NACA0012 potential flow and the boundary layer edge velocity from the SU2 converged solution at $R_{chord} = 10^{7}$ shown in Fig. 35 indicates that further iterations would produce little change in the UVP drag.

Notice that the trailing edge values of β_c in Fig. 29 never exceed $β_{c} = 18$ ⁠. However, the user needs to be aware that during each iteration, the $β_{c_{i}}$ calculated from the current $R_{τ_{i}} (ξ)$ distribution becomes very large, far larger than 18, as the trailing edge is approached where $u_{τ}$ becomes relatively small (but not 0 because of the cusped trailing edge) and the displacement and momentum thickness Reynolds numbers in Eq. (43) become large. However note the decreasing values of $R_{τ_{T E}}$ with each iteration. It turns out that when $β_{c_{i}}$ is used to calculate $R_{τ_{i + 1}}$ ⁠, the calculation reaches the trailing edge at a value of $R_{τ_{i + 1_{T E}}}$ that is considerably less than $R_{τ_{i_{T E}}}$ with a corresponding $β_{c_{i}}$ that is less than 18.

Figures 26–32 show the convergence of the various functions required to integrate the Kármán equation (41) and determine the airfoil drag. Two other boundary layer properties of importance are the channel half height thickness (overall boundary layer thickness), δ_h, and the shape factor, $H = δ_{1} / δ_{2}$ ⁠, shown in Figs. 32 and 33. The overall thickness is determined from

\frac{δ_{h}}{c} = \frac{R_{τ} (ξ)}{\frac{u_{\infty} c}{ν} {(\frac{C_{f} (ξ)}{2})}^{1 / 2} U (ξ)} .

(70)

FIG. 26.

View large Download slide

Iterative distributions of the wall friction coefficient C_f for the J0012 airfoil. Inset shows the covergence of the trailing edge value of $δ_{h} / c$ with iteration number. Final trailing edge value is $C_{f} = 0.001 055$ ⁠.

FIG. 27.

View large Download slide

Iterative distributions of the displacement thickness $δ_{1} / c$ for the J0012 airfoil. Inset shows the covergence of the trailing edge value of $δ_{1} / c$ with iteration number. Final trailing edge value is $δ_{1} / c = 0.006 201$ ⁠.

FIG. 28.

View large Download slide

Iterative distributions of the momentum thickness $δ_{2} / c$ for the J0012 airfoil. Inset shows the convergence of the trailing edge value of $δ_{2} / c$ with iteration number. Final trailing edge value is $δ_{2} / c = 0.003 759$ ⁠. All eight iterations overlay almost exactly.

FIG. 29.

View large Download slide

Iterative distributions of the modified Clauser parameter, β_c, for the J0012 airfoil. Inset shows the covergence of the trailing edge value of β_c with iteration number. Final trailing edge value is $β_{c} = 17.238$ ⁠. Note that $β_{c} = 0$ over the first iteration.

FIG. 30.

View large Download slide

Iterative distributions of the UVP parameter b for the J0012 airfoil. Inset shows the covergence of the trailing edge value of b with iteration number. Final trailing edge value is $b_{T E} = 0.041 56$ ⁠.

FIG. 31.

View large Download slide

Iterative distributions of the UVP parameter n for the J0012 airfoil. Inset shows the covergence of the trailing edge value of n with iteration number. Final trailing edge value is $n_{T E} = 6.0994$ ⁠.

FIG. 32.

View large Download slide

Iterative distributions of boundary layer thickness $δ_{h} / c$ for the J0012 airfoil. Inset shows the covergence of the trailing edge value of $δ_{h} / c$ with iteration number. Final trailing edge value is $δ_{h} / c = 0.024 98$ ⁠.

FIG. 33.

View large Download slide

Iterative distributions of the boundary layer shape factor $H = δ_{1} / δ_{2}$ for the J0012 airfoil. Inset shows the convergence of the trailing edge value of H with iteration number. Final trailing edge value is $H = 1.6496$ ⁠.

With the friction coefficient known over the airfoil surface, the viscous drag coefficient is determined using Eq. (8). The convergence of the J0012 and NACA0012 viscous drag coefficients at $R_{chord} = 10^{7}$ are shown in Fig. 34. Figure 35 shows the potential flow over the NACA0012 wing.

FIG. 34.

View large Download slide

Convergence of the J0012 (blue open circles) and NACA0012 (magenta open circles) viscous drag coefficients at $R_{chord} = 10^{7}$ ⁠. Converged values are $C_{d v_{J 0012}} = 0.006 30$ and $C_{d v_{NACA 0012}} = 0.006 39$ ⁠.

FIG. 35.

View large Download slide

Potential flow velocity distribution about the NACA0012 airfoil (magenta) compared with the inviscid, boundary layer edge velocity (black) from SU2 at $R_{chord} = 10^{7}$ ⁠.

VII. HIGH REYNOLDS NUMBER

In the presence of a pressure gradient, the wall parameters $(k, a, m)$ are kept constant at the ZPG boundary layer values in Table I while the wake parameters $(b, n)$ are treated as functions of β_c using the correlations (45) and (47). Figure 36 shows the shape function, $ϕ$ Eq. (23), for a range of values of β_c. Changes in $ϕ$ with β_c are smooth and monotonic. A polynomial of high order is a safe, reliable fit in this situation since there is no question of $y / δ_{h}$ falling outside of the range $0 \leq y / δ_{h} \leq 1$ ⁠. The family of shape functions shown in Fig. 36 is accurately approximated by the eighth-order polynomial

\begin{matrix} ϕ (β_{c}, \frac{y}{δ_{h}}) = c_{0} (β_{c}) + c_{1} (β_{c}) (\frac{y}{δ_{h}}) + c_{2} (β_{c}) {(\frac{y}{δ_{h}})}^{2} \\ + c_{3} (β_{c}) {(\frac{y}{δ_{h}})}^{3} + c_{4} (β_{c}) {(\frac{y}{δ_{h}})}^{4} + c_{5} (β_{c}) {(\frac{y}{δ_{h}})}^{5} \\ + c_{6} (β_{c}) {(\frac{y}{δ_{h}})}^{6} + c_{7} (β_{c}) {(\frac{y}{δ_{h}})}^{7} + c_{8} (β_{c}) {(\frac{y}{δ_{h}})}^{8} . \end{matrix}

(71)

FIG. 36.

View large Download slide

The shape function, $ϕ$ ⁠, for various β_c in the range $- 1.0 < β_{c} < 18.0$ evaluated at $R_{τ} = 10^{6}$ ⁠. Successive curves are in increments of 0.5.

The coefficients $c_{0} (β_{c})$ to $c_{8} (β_{c})$ are approximated by tenth-order polynomial functions of the pressure gradient parameter β_c and are strictly limited to the range $- 1.0 < β_{c} < 18.0$ indicated in Fig. 36. The expressions for c₀ to c₈ are provided in Appendix A. It is important to note that $ϕ (β_{c}, y / δ_{h})$ in Eq. (71) is universal in the sense that it applies to any pressure gradient wall flow that falls in the range $- 1.0 < β_{c} < 18.0$ ⁠. Higher order fits could be used to increase the range of β_c beyond 18.0.

The explicit high Reynolds number, variable pressure gradient form of the UVP reduces to the simple form

\begin{matrix} \lim_{k R_{τ} > 2000} k u^{+} = ln (k y^{+}) + ϕ (k a, m, β_{c}, \frac{y}{δ_{h}}) . \end{matrix}

(72)

At high Reynolds number the parameter k, essentially the Kármán constant, is simply a scale factor on y⁺, u⁺ and the damping length scale a. Since b and n are correlated with the modified Clauser parameter, β_c, the remaining optimal parameters governing the k-scaled UVP are only ka and m.

With $ϕ (k a, m, β_{c}, y / δ_{h})$ known explicitly as a bivariate polynomial, the integrals (26) and (27) can be carried out. This leads to polynomial expressions for the friction and thickness functions $F_{0} (β_{c}, R_{τ}), F_{1} (β_{c}, R_{τ}), F_{2} (β_{c}, R_{τ})$ ⁠, and $F_{3} (β_{c}, R_{τ})$ ⁠. Once these functions are known, the UVP integral method no longer requires the computation of nested integrals. As a result, the calculation time is independent of the Reynolds number. The viscous drag coefficient for any chord Reynolds number can be determined in a few seconds. Above about $R_{chord} = 10^{6}$ ⁠, this high Reynolds number explicit UVP method matches the results using the integral form of the UVP almost exactly. See Fig. 38 and Table IV.

TABLE IV.

Chord Reynolds number, $C_{d_{v}}$ for the J0012 airfoil, $C_{d_{v}}$ for the NACA0012 airfoil using both the integral and explicit forms of the UVP, $C_{d_{v}}$ and $C_{d_{p}}$ computed using SU2,³³ tripped $C_{d_{v}}$ data from Ladson.⁴⁰

R_chord	$C_{d_{v}_{J 0012}}$	$C_{d_{v}_{NACA 0012_{UVPintegral}}}$	$C_{d_{v}_{NACA 0012_{UVPexplicit}}}$	$C_{d_{v}_{NACA 0012_{S U 2}}}$	$C_{d_{p}_{NACA 0012_{S U 2}}}$	$C_{d_{v}_{Ladso n_{80 - 180 grit}}}$
10⁵	0.014 797 0	0.014 817 4	0.013 499 8	0.012 746 0	0.004 738 7	⋯
$5 \times 10^{5}$	0.010 352 9	0.010 397 7	0.010 032 5	0.010 108 7	0.002 664 7	⋯
10⁶	0.009 062 1	0.009 147 5	0.009 003 4	0.008 930 0	0.002 252 9	⋯
$2 \times 10^{6}$	0.008 040 3	0.008 147 7	0.008 057 8	⋯	⋯	0.008 53
$4 \times 10^{6}$	0.007 214 1	0.007 295 5	0.007 244 0	⋯	⋯	0.007 33
$5 \times 10^{6}$	0.006 977 5	0.007 050 9	0.007 005 5	0.007 039 7	0.001 608 7	⋯
$6 \times 10^{6}$	0.006 787 1	0.006 862 6	0.006 820 3	⋯	⋯	0.006 82
$8.95 \times 10^{6}$	0.006 397 9	0.006 488 3	0.006 435 3	⋯	⋯	0.006 51
10⁷	0.006 296 5	0.006 394 3	0.006 294 3	0.006 431 9	0.001 495 7	⋯
$1.2 \times 10^{7}$	0.006 137 3	0.006 228 17	0.006 170 8	⋯	⋯	0.006 53
$5 \times 10^{7}$	0.005 112 8	0.005 102 1	0.005 082 9	0.005 128 0	0.001 092 0	⋯
10⁸	0.004 716 9	0.004 716 8	0.004 6502	0.004 696 2	0.000 933 7	⋯
10⁹	0.003 512 5	0.003 547 7	0.003 535 7	0.003 400 6	0.000 560 9	⋯
10¹⁰	0.002 755 5	0.002 814 7	0.002 767 3	⋯	⋯	⋯
10¹¹	0.002 204 9	0.002 147 2	0.002 217 3	⋯	⋯	⋯
10¹²	0.001 802 8	0.001 764 5	0.001 812 6	⋯	⋯	⋯

R_chord	$C_{d_{v}_{J 0012}}$	$C_{d_{v}_{NACA 0012_{UVPintegral}}}$	$C_{d_{v}_{NACA 0012_{UVPexplicit}}}$	$C_{d_{v}_{NACA 0012_{S U 2}}}$	$C_{d_{p}_{NACA 0012_{S U 2}}}$	$C_{d_{v}_{Ladso n_{80 - 180 grit}}}$
10⁵	0.014 797 0	0.014 817 4	0.013 499 8	0.012 746 0	0.004 738 7	⋯
$5 \times 10^{5}$	0.010 352 9	0.010 397 7	0.010 032 5	0.010 108 7	0.002 664 7	⋯
10⁶	0.009 062 1	0.009 147 5	0.009 003 4	0.008 930 0	0.002 252 9	⋯
$2 \times 10^{6}$	0.008 040 3	0.008 147 7	0.008 057 8	⋯	⋯	0.008 53
$4 \times 10^{6}$	0.007 214 1	0.007 295 5	0.007 244 0	⋯	⋯	0.007 33
$5 \times 10^{6}$	0.006 977 5	0.007 050 9	0.007 005 5	0.007 039 7	0.001 608 7	⋯
$6 \times 10^{6}$	0.006 787 1	0.006 862 6	0.006 820 3	⋯	⋯	0.006 82
$8.95 \times 10^{6}$	0.006 397 9	0.006 488 3	0.006 435 3	⋯	⋯	0.006 51
10⁷	0.006 296 5	0.006 394 3	0.006 294 3	0.006 431 9	0.001 495 7	⋯
$1.2 \times 10^{7}$	0.006 137 3	0.006 228 17	0.006 170 8	⋯	⋯	0.006 53
$5 \times 10^{7}$	0.005 112 8	0.005 102 1	0.005 082 9	0.005 128 0	0.001 092 0	⋯
10⁸	0.004 716 9	0.004 716 8	0.004 6502	0.004 696 2	0.000 933 7	⋯
10⁹	0.003 512 5	0.003 547 7	0.003 535 7	0.003 400 6	0.000 560 9	⋯
10¹⁰	0.002 755 5	0.002 814 7	0.002 767 3	⋯	⋯	⋯
10¹¹	0.002 204 9	0.002 147 2	0.002 217 3	⋯	⋯	⋯
10¹²	0.001 802 8	0.001 764 5	0.001 812 6	⋯	⋯	⋯

Referring back to Sec. II A, it should be pointed out that when Eq. (72) is used to generate the thicknesses, Eqs. (26) and (27), at high Reynolds number there is an error in neglecting the integration from the wall to the beginning of the log region. The log region for a ZPG boundary layer begins at approximately $y^{+} = 132$ ⁠.¹⁷ Figure 37 shows the ratio of the integration using the integral form of the UVP from the wall to $y^{+} = 132$ divided by integration using the explicit form of the UVP from 0 to $R_{τ}$ ⁠, as a function of $R_{τ}$ ⁠. Above $R_{τ} = 10^{5}$ the error is less than $0.5 %$ ⁠.

FIG. 37.

View large Download slide

Estimate of the error when the explicit form of the UVP, Eq. (72), is used to determine boundary layer thickness Reynolds numbers F₁ (dotted line) and F₂ (dashed line) neglecting the near wall region. Log region for a ZPG boundary layer begins at $y^{+} = 132$ ⁠.¹⁷

VIII. VISCOUS DRAG COEFFICIENT OF THE JOUKOWSKY AND NACA 0012 AIRFOILS

The J0012 airfoil was chosen for this study because it has an exact, and, therefore, highly accurate, potential flow solution and a well-defined cusped trailing edge. The NACA 0012 (with closed trailing edge) was chosen because this airfoil has been used as a standard for testing computational results, for assessing wall interference effects, and for developing correction procedures for wind tunnel testing. There is a variety of formulas in the literature for the NACA 0012 airfoil surface, consequently there is a variety of profiles used to study the flow, some with unclosed trailing edges and some with closed trailing edges. The profile used in this study is defined by

\frac{\tilde{y}}{c} = a_{1} {(\frac{\tilde{x}}{c})}^{\frac{1}{2}} + a_{2} (\frac{\tilde{x}}{c}) + a_{3} {(\frac{\tilde{x}}{c})}^{2} + a_{4} {(\frac{\tilde{x}}{c})}^{3} + a_{5} {(\frac{\tilde{x}}{c})}^{4}

(73)

with constants, $a_{1} = 0.177 349 856, a_{2} = - 0.075 600 000, a_{3} = - 0.212 843 959 1, a_{4} = 0.173 640 303 0$ ⁠, and $a_{5} = - 0.062 546 200 2$ ⁠, which produce a closed trailing edge.

The NACA 0012 airfoil drag was determined using the UVP integral method and compared with computations of the drag using the aerodynamic design software suite, SU2.³³ For the presentation, the boundary layer thicknesses, drag coefficient, and Reynolds number are normalized by the airfoil chord. Figure 35 shows the NACA0012 potential flow $U (\tilde{ξ})$ determined using conformal mapping, compared with the free stream $U (\tilde{ξ})$ computed from the SU2 pressure distribution at $R_{chord} = 10^{7}$ ⁠. The pressure drag from SU2 at this Reynolds number was $C_{d_{p_{S U 2}}} = 0.001 496$ ⁠.

The final C_dv vs R_chord results are shown in Fig. 38 and the data used to produce this figure are provided in Table IV. The nearly overlapping blue (J0012) and magenta (NACA0012) curves in Fig. 38 are generated using the UVP integral method with the mean wall parameter values, $(\bar{k}, \bar{a}, \bar{m})$ ⁠, given in Table I. Recalling the discussion of Figs. 6 and 7 in Sec. II B, increasing or decreasing the optimal parameter values by one standard deviation would lead to a $7.5 %$ change in C_f at $R_{τ} = 5000$ ⁠. The magenta error bar in Fig. 38 is approximately the variation of $C_{d_{v}}$ that would occur if optimal parameters in the UVP were increased or decreased by the σ's given in Table I. The position of the magenta error bar at $R_{chord} = 7 \times 10^{6}$ is approximately at $R_{τ} = 5000$ ⁠. It should be noted that, while the error associated with changing the optimal parameters this way actually depends on $R_{τ}$ ⁠, the magenta error bar provides a reasonable approximation to the error in the calculated $C_{d_{v}}$ based on the set of boundary layer cases evaluated to date.

FIG. 38.

View large Download slide

$C_{d_{v}}$ for the J0012 airfoil (blue curve, see inset) and NACA 0012 airfoil (magenta curve) determined using the integral form of the UVP. The magenta vertical error bar is the variation of $C_{d_{v}}$ that would occur if optimal parameters in the UVP were increased or decreased by the variances given in Table I. The green curve is $C_{d_{v}}$ calculated using the high Reynolds number explicit form of the UVP. Red filled circles are $C_{d_{v}}$ computed using SU2³³ with the Spalart–Allmaras model.³⁷ Open black circles are NACA 0012 tripped data from Ladson⁴⁰ averaged over grit sizes 80–180. The black dashed curve is the correlation of tripped NACA 0012 data from McCroskey.⁴¹ The black vertical error bar represents the variation (⁠ $Δ C_{d v} = 0.001 2$ ⁠) in the tripped data evaluated by McCroskey.

The green (NACA0012) curve was generated using the explicit form of the UVP which allows high Reynolds number cases to be computed in a few seconds. The three UVP curves are almost identical at Reynolds numbers above $R_{chord} > 10^{6}$ ⁠.

The open black circles in Fig. 38 are from Ladson⁴⁰ who used carborundum strips placed at $\tilde{x} / c = 0.05$ to generate a tripped boundary layer over an NACA0012 wing. The data points shown are averages over grit sizes 80, 120, and 180 and Mach numbers 0.15 and 0.30. The black dashed curve in Fig. 38 is a fit to high quality tripped data from several sources (including Ladson) collected and evaluated by McCroskey⁴¹

C_{d} = 0.001 7 + \frac{0.91}{{(L o g_{10} (R_{chord}))}^{2.58}} .

(74)

Note that an interpolation of the pressure drag coefficient, $C_{d p_{S U 2}}$ is subtracted from the C_d data in Ladson⁴⁰ and McCroskey⁴¹ to generate the black open circles and dashed line in Fig. 38.

The McCroskey correlation (74) lies slightly below the UVP results and has a slightly smaller slope and intersects the blue, magenta, and green curves at about $R_{chord} = 10^{7}$ ⁠. The black error bar in Fig. 38 at $R_{chord} = 10^{7}$ indicates the variation in the data in Fig. 4 of McCroskey.⁴¹

The filled red circles in Fig. 38 are results from SU2³³ using the Spalart–Allmaras model.³⁷ No transition model was used in the computations. Between $R_{chord} = 10^{6}$ and 10⁹, the difference between the UVP drag and the SU2 data are generally less than 0.0002. The results for $R_{chord} = 10^{9}$ were obtained using a procedure consistent with the lower R_chord cases. The finest grid available from the Langley Turbulence Modeling Resource webpage⁴² gave a wall spacing of $10^{- 7}$ m corresponding to a minimum $y^{+} ≳ 1$ over significant portions of the airfoil instead of the preferred $y^{+} \leq 1$ ⁠. This inability to accurately capture the viscous sublayer over much of the airfoil may explain the slight under-prediction of SU2 relative to the UVP at this high Re. While further validation in high Re flows is still needed, this case does demonstrate close agreement of the UVP method with established turbulence models at high Re.

The good agreement between the UVP and SU2 results is to be expected given the fact that the flat plate wall friction coefficients generated by both Eqs. (30) and (34) are in good agreement.

The results in Fig. 38 are all consistent with the point made in Sec. II C that the UVP represents a tripped boundary layer with a critical Reynolds number close to the minimum needed to sustain turbulent flow.

A. Drag force comparison between the UVP and SU2 NACA0012 results at $R_{chord} = 10^{7}$

To investigate the UVP and SU2 drag coefficient results in a little more detail, the viscous drag force distribution along the NACA 0012 airfoil is shown in Figs. 39 and 40.

FIG. 39.

View large Download slide

Viscous drag force on the NACA0012 airfoil at $R_{c} = 10^{7}$ as a function of position. Comparison between the UVP first iteration (magenta, dashed) and the UVP final iteration (magenta, solid).

FIG. 40.

View large Download slide

Viscous drag force on the NACA0012 airfoil at $R_{c} = 10^{7}$ as a function of position. Comparison between the UVP final iteration (magenta) and the viscous drag computed using SU2 with the Spalart–Allmaras model (black).

Figure 39 shows the drag force distribution generated by the first and last iterations of the UVP integral method. The first iteration is remarkably good considering that there is no adjustment of b or n. Most of the variation in drag is accounted for by the relationship between the free stream velocity and wall friction inherent to the UVP. The largest difference between the first and last iterations occurs near the airfoil trailing edge where the adverse pressure gradient adjustments of b and n are needed to properly capture the drag force as well as other boundary layer characteristics in that region.

In Fig. 40 the SU2 distribution, shown in black, is compared to the final iteration of the UVP distribution, shown in magenta. The two curves almost overlap over most of the wing. Over the last 2% of the chord the UVP drag drops rapidly to almost zero near the finite angle trailing edge, while the SU2 curve rises slightly in response to the beginning of the airfoil wake.

IX. CONCLUDING REMARKS

The universal velocity profile, developed originally for pipe flow,¹¹ provides a very accurate approximation to channel flow and boundary layer flow with zero, favorable, and adverse pressure gradient. A modified Clauser parameter, β_c, is defined for the first time and shown in Figs. 18 and 19 to correlate well with the parameters b and n that characterize the wake portion of the UVP. This result provides a complete characterization of the Reynolds number and pressure gradient dependence of the UVP enabling it to be used to create the new integral method presented in this paper. Interestingly, the first iteration on C_dv with b and n fixed at the zero pressure gradient values is remarkably close to the final value, although the distributions of $R_{τ}$ and C_f generated by the modified Clauser parameter are essential to properly capture the trailing edge flow.

For all calculations, the parameters that characterize the wall flow, $(\bar{k}, \bar{a}, \bar{m})$ were held constant at the values in Table I independent of the pressure gradient.¹⁷ As additional boundary layer data become available, the average optimal parameter values will continue to change slightly although the impact on calculated values of $C_{d_{v}}$ will be decreasingly small as the sample size grows.

The $C_{d_{v}}$ generated by the UVP in Fig. 38 and Table IV is generally in very good agreement with the McCroskey/Ladson data and the $C_{d_{v}}$ data generated by SU2. The conclusion that can be drawn from the results in Table I and Fig. 38 is that the UVP integral method is an effective method for determining viscous drag on smooth aerodynamic shapes. The conclusion from the discussion of Figs. 5 and 9 and the comparison between the UVP and Ladson/McCroskey data in Fig. 38 is that when the UVP is used beginning at the airfoil leading edge, the drag should be regarded as that of a tripped boundary layer with a critical Reynolds number close to the minimum required to sustain turbulent flow.

To handle an untripped case, the laminar and transitional part of the flow will need to be determined using a transition criterion for the boundary layer along with additional information about the free stream. One might use the well-known model of boundary layer transition where transition is considered complete after the amplification of unstable initial disturbances by a factor e^N.^43,44 For a boundary layer developing in a very low free stream turbulence environment, for example, Schubauer and Klebanoff,³⁰ N would be selected to be in the range 9–12, whereas for a boundary layer developing in a highly turbulent free stream one might choose $N \leq 1$ ⁠. The problem of how best to incorporate a transition criterion into the universal velocity profile is a subject of current research.

Large eddy simulation of wall flows is faced with the conflict between the need to reduce grid resolution near the wall in order to reduce computational cost and the requirement to accurately determine the wall friction. Since the large eddy simulation (LES) grid under-resolves the fine scale motions near the wall that are responsible for generating the wall stress, empirical wall functions need to be used between the first grid point and the wall in order to determine the friction. These functions inevitably encounter the log-layer mismatch problem.^45–48 The UVP provides an accurate value of the mean stream-wise and wall-normal velocity components as well as the turbulent shear stress throughout the viscous wall layer without assuming a log or power law profile and may provide a useful alternative to the current approach. That the UVP can provide an accurate connection to the wall is evidenced by the fact that virtually all the experimental non-DNS datasets that we have approximated thus far only contain data that is well outside the wall layer.

The inherent dependence of the UVP on Reynolds number, extended to include the effect of pressure gradient enables it to be used as the basis of a new method for integrating the Kármán equation for a wide variety of attached, wall bounded flows. There is really no practical limit to the Reynolds number that can be computed suggesting that, with modification for the effects of compressibility and/or roughness, it can be applied to very large-scale aerodynamic, hydrodynamic, and geophysical flows.

Ultimately, fundamental questions about the optimal parameters will need to be answered as to their numerical values, dependence on flow geometry (seen clearly in Table I), and possible weak dependence on Reynolds number, free stream turbulence, surface roughness, Mach number, and so on. It is not clear when the high Reynolds number boundary layer data required to address these questions will become available.

ACKNOWLEDGMENTS

The authors would like to express appreciation to Stanford University for continued support. The computational resources made available through the Stanford Sherlock High Performance Computing Cluster are gratefully acknowledged. All grids used for the computations came from the Langley Turbulence Modeling Resource webpage.

AUTHOR DECLARATIONS

Conflict of Interest

The authors have no conflicts to disclose.

Author Contributions

Brian Joseph Cantwell: Conceptualization (lead); Data curation (lead); Formal analysis (lead); Funding acquisition (lead); Investigation (equal); Methodology (lead); Resources (equal); Software (equal); Supervision (lead); Writing – original draft (lead); Writing – review & editing (lead). Eylul Bilgin: Conceptualization (supporting); Data curation (supporting); Formal analysis (supporting); Investigation (equal); Methodology (supporting); Resources (equal). Jacob Troy Needels: Investigation (equal); Methodology (supporting); Resources (equal); Software (equal).

DATA AVAILABILITY

The data that support the findings of this study are available within the article.

APPENDIX A: SHAPE FUNCTION COEFFICIENTS

The polynomial approximations to the coefficients in Eq. (71) are listed below. These approximations are only valid in the range $- 1.0 < β_{c} < 18.0$ ⁠. Outside that range the higher order terms will produce very large errors,

\begin{matrix} c_{0} (β_{c}) = 3.185 213 759 059 279 5 - 0.000 097 351 124 875 088 01 β_{c} \\ + 0.000 857 641 541 537 847 5 β_{c}^{2} \\ - 0.000 134 299 418 144 835 46 β_{c}^{3} \\ - 0.000 020 654 101 916 920 754 β_{c}^{4} \\ + 0.000 015 020 911 822 593 57 β_{c}^{5} \\ - 2.968 623 256 456 977 \times 10^{- 6} β_{c}^{6} \\ + 2.963 701 061 289 852 4 \times 10^{- 7} β_{c}^{7} \\ - 1.631 850 736 220 117 2 \times 10^{- 8} β_{c}^{8} \\ + 4.735 140 769 133 537 \times 10^{- 10} β_{c}^{9} \\ - 5.670 359 008 088 255 5 \times 10^{- 12} β_{c}^{10}, \end{matrix}

(A1)

\begin{matrix} c_{1} (β_{c}) = 0.177 875 635 158 609 52 - 0.243 068 542 472 120 7 β_{c} \\ - 0.182 315 590 796 708 98 β_{c}^{2} + 0.114 167 281 428 819 74 β_{c}^{3} \\ - 0.032 508 461 93 233 937 β_{c}^{4} + 0.004 881 950 187 888 647 β_{c}^{5} \\ - 0.000 419 972 498 600 616 4 β_{c}^{6} \\ + 0.000 020 756 027 245 732 183 β_{c}^{7} \\ - 5.360 810 884 925 711 \times 10^{- 7} β_{c}^{8} \\ + 4.864 234 131 520 179 6 \times 10^{- 9} β_{c}^{9} \\ + 2.745 541 114 505 071 \times 10^{- 11} β_{c}^{10}, \end{matrix}

(A2)

\begin{matrix} c_{2} (β_{c}) = 8.752 630 174 726 951 + 5.383 284 010 214 323 β_{c} \\ - 1.802 896 664 601 566 7 β_{c}^{2} - 0.077 034 046 540 169 49 β_{c}^{3} \\ + 0.335 381 660 158 488 7 β_{c}^{4} - 0.103 409 008 148 911 44 β_{c}^{5} \\ + 0.015 699 101 594 675 673 β_{c}^{6} \\ - 0.001 372 392 667 304 845 8 β_{c}^{7} \\ + 0.000 070 043 102 706 461 β_{c}^{8} \\ - 1.940 086 907 438 102 7 \times 10^{- 6} β_{c}^{9} \\ + 2.254 603 431 206 475 2 \times 10^{- 8} β_{c}^{10}, \end{matrix}

(A3)

\begin{matrix} c_{3} (β_{c}) = - 29.320 201 478 527 7 - 13.009 083 141 290 633 β_{c} \\ + 9.891 072 103 799 035 β_{c}^{2} + 0.286 430 751 940 286 4 β_{c}^{3} \\ - 1.264 366 997 570 824 β_{c}^{4} + 0.357 951 024 377 768 16 β_{c}^{5} \\ - 0.051 953 118 683 825 3 β_{c}^{6} \\ + 0.004 449 565 462 995 644 β_{c}^{7} \\ - 0.000 225 332 996 880 537 36 β_{c}^{8} \\ + 6.231 250 458 034 746 \times 10^{- 6} β_{c}^{9} \\ - 7.250 246 433 619 231 \times 10^{- 8} β_{c}^{10}, \end{matrix}

(A4)

\begin{matrix} c_{4} (β_{c}) = 63.370 268 858 329 155 + 21.206 924 878 802 89 β_{c} \\ - 28.490 362 428 461 17 β_{c}^{2} - 1.383 788 687 404 143 2 β_{c}^{3} \\ + 2.999 048 905 796 227 5 β_{c}^{4} - 0.761 144 546 069 162 7 β_{c}^{5} \\ + 0.104 860 072 527 180 98 β_{c}^{6} \\ - 0.008 803 651 704 295 447 β_{c}^{7} \\ + 0.000 444 105 810 851 922 5 β_{c}^{8} \\ - 0.000 012 319 542 850 362 865 β_{c}^{9} \\ + 1.441 663 408 041 264 3 \times 10^{- 7} β_{c}^{10}, \end{matrix}

(A5)

\begin{matrix} c_{5} (β_{c}) = - 87.693 764 826 178 6 - 22.352 131 224 850 766 β_{c} \\ + 47.942 643 158 529 59 β_{c}^{2} + 3.144 148 516 656 517 β_{c}^{3} \\ - 4.414 688 172 210 721 β_{c}^{4} + 1.011 217 929 468 684 7 β_{c}^{5} \\ - 0.132 607 179 383 751 05 β_{c}^{6} \\ + 0.010 971 698 041 512 124 β_{c}^{7} \\ - 0.000 554 916 304 796 180 4 β_{c}^{8} \\ + 0.000 015 534 289 801 032 056 β_{c}^{9} \\ - 1.837 286 959 416 604 7 \times 10^{- 7} β_{c}^{10}, \end{matrix}

(A6)

\begin{matrix} c_{6} (β_{c}) = 70.968 812 707 064 7 + 12.583 124 149 174 417 β_{c} \\ - 46.888 434 832 647 61 β_{c}^{2} - 3.734 105 839 184 008 5 β_{c}^{3} \\ + 3.884 301 374 624 679 β_{c}^{4} - 0.802 627 524 234 329 7 β_{c}^{5} \\ + 0.099 762 336 519 354 65 β_{c}^{6} \\ - 0.008 144 173 595 693 59 β_{c}^{7} \\ + 0.000 414 872 495 566123 67 β_{c}^{8} \\ - 0.000 011 778 008 804 226 057 β_{c}^{9} \\ + 1.413 661 707 265 926 8 \times 10^{- 7} β_{c}^{10}, \end{matrix}

(A7)

\begin{matrix} c_{7} (β_{c}) = - 29.738 042 179 000 644 - 2.266 916 912 193 806 β_{c} \\ + 24.689 251 860 732 828 β_{c}^{2} + 2.271 539 156 916 986 β_{c}^{3} \\ - 1.867 706 886 189 081 3 β_{c}^{4} + 0.344 688 708 245 448 8 β_{c}^{5} \\ - 0.040 023 433 517 290 326 β_{c}^{6} \\ + 0.003 211 634 372 538 153 β_{c}^{7} \\ - 0.000 165 444 382 500 483 54 β_{c}^{8} \\ + 4.791 043 532 302 556 \times 10^{- 6} β_{c}^{9} \\ - 5.865 386 244 138 967 5 \times 10^{- 8} β_{c}^{10}, \end{matrix}

(A8)

\begin{matrix} c_{8} (β_{c}) = 4.625 148 665 477 973 - 0.505 275 164 231 022 5 β_{c} \\ - 5.397 029 748 250 34 β_{c}^{2} - 0.562 142 236 631 621 3 β_{c}^{3} \\ + 0.375 511 158 572 171 93 β_{c}^{4} - 0.060 255 387 771 924 11 β_{c}^{5} \\ + 0.006 292 183 567 809 945 β_{c}^{6} \\ - 0.000 488 834 388 113 832 5 β_{c}^{7} \\ + 0.000 025 608 788 685 571 06 β_{c}^{8} \\ - 7.651 119 404 966 455 \times 10^{- 7} β_{c}^{9} \\ + 9.652 853 992 263 516 \times 10^{- 9} β_{c}^{10} . \end{matrix}

(A9)

APPENDIX B: UVP INTEGRAL METHOD FLOW CHART

A flow chart showing the procedure for the UVP integral boundary layer method is shown in Fig. 41.

FIG. 41.

View large Download slide

Flowchart summarizing the UVP method.

References

1.

L.

Prandtl

, “

The mechanics of viscous fluids

,” in

Aerodynamic Theory

, edited by

W. F.

Durand

(

Guggenheim Fund for the Promotion of Aeronautics

,

1934

), Vol.

III

, pp.

80

–

90

.

Google Scholar

Crossref

2.

B.

Thwaites

,

Incompressible Aerodynamics

(

Oxford University Press

,

1960

), pp.

62

–

64

.

Google Scholar

3.

K.

Pohlhausen

, “

Zur naherungsweisen integration der differentialgleichungen der laminaren reilbungsschicht

,”

ZAMM

1

,

252

–

260

(

1921

).

https://doi.org/10.1002/zamm.19210010402

Google Scholar

Crossref

4.

J.

Majdalani

and

L.

Xuan

, “

On the Kármán momentum-integral approach and the Pohlhausen paradox

,”

Phys. Fluids

32

,

123605

(

2020

).

https://doi.org/10.1063/5.0036786

Google Scholar

Crossref

5.

M. R.

Head

, “

Entrainment in the turbulent boundary layer

,” Aeronautical Research Council Reports and Memoranda No. 3152,

1960

.

Google Scholar

6.

H.

Ludwieg

and

W.

Tillman

, “

Investigations of wall shearing stress in turbulent boundary layers

,” NACA Technical Memorandum No. 1285,

1950

.

Google Scholar

7.

H.

Schlichting

and

K.

Gersten

,

Boundary Layer Theory

(

Springer

,

2000

).

Google Scholar

Crossref

8.

B. J.

Cantwell

,

Introduction to Symmetry Analysis

(

Cambridge University Press

,

2002

).

Google Scholar

9.

T.

von Kármán

, “

Uber laminaire und turbulente reibung

,”

ZAMM

1

,

233

–

252

(

1921

).

https://doi.org/10.1002/zamm.19210010401

Google Scholar

Crossref

10.

B. J.

Cantwell

, “

Integral measures of the zero pressure gradient boundary layer over the Reynolds number range

0 \leq R_{τ} < \infty

,”

Phys. Fluids

33

,

085108

(

2021

).

https://doi.org/10.1063/5.0061535

Google Scholar

Crossref

11.

B. J.

Cantwell

, “

A universal velocity profile for smooth wall pipe flow

,”

J. Fluid Mech.

878

,

834

(

2019

).

https://doi.org/10.1017/jfm.2019.669

Google Scholar

Crossref

12.

T.

von Kármán

, “

Mechanical similitude and turbulence

,” NACA Technical Memorandum No. 611,

1931

.

Google Scholar

13.

L.

Prandtl

, “

The mechanics of viscous fluids

,” in

Aerodynamic Theory

, edited by

W. F.

Durand

(

Guggenheim Fund for the Promotion of Aeronautics

,

1934

), Vol III, pp.

102

–

112

.

Google Scholar

Crossref

14.

E. R.

van Driest

, “

On turbulent flow near a wall

,”

J. Aeronaut. Sci.

23

,

1007

–

1011

(

1956

).

https://doi.org/10.2514/8.3713

Google Scholar

Crossref

15.

W.

Kollmann

, “

Asymptotic properties of mixing length closures for turbulent pipe flow

,”

Phys. Fluids

32

,

115126

(

2020

).

https://doi.org/10.1063/5.0030328

Google Scholar

Crossref

16.

T.

Cebeci

and

P.

Bradshaw

,

Momentum Transfer in Boundary Layers

(

Hemisphere Publishing Corporation

,

1977

), Chap. 6.4, p.

174

.

Google Scholar

17.

M. A.

Subrahmanyam

,

B. J.

Cantwell

, and

J. J.

Alonso

, “

A universal velocity profile for turbulent wall flows including adverse pressure gradient boundary layers

,”

J. Fluid Mech.

933

,

A16

(

2022

).

https://doi.org/10.1017/jfm.2021.998

Google Scholar

Crossref

18.

D. I.

Pullin

,

M.

Inoue

, and

N.

Saito

, “

On the asymptotic state of high Reynolds number, smooth-wall turbulent flows

,”

Phys. Fluids

25

,

015116

(

2013

).

https://doi.org/10.1063/1.4774335

Google Scholar

Crossref

19.

M. A.

Subrahmanyam

,

B. J.

Cantwell

, and

J. J.

Alonso

, “

A universal velocity profile for turbulent wall flows

,” AIAA Paper No. 2021-0061,

2021

.

Google Scholar

20.

M. P.

Simens

,

J.

Jiménez

,

S.

Hoyas

, and

Y.

Mizuno

, “

A high-resolution code for turbulent boundary layers

,”

J. Comput. Phys.

228

,

4218

–

4231

(

2009

).

https://doi.org/10.1016/j.jcp.2009.02.031

Google Scholar

Crossref

21.

G.

Borrell

,

J. A.

Sillero

, and

J.

Jiménez

, “

A code for direct numerical simulation of turbulent boundary layers at high Reynolds numbers in BG/P supercomputers

,” in

23rd International Conference on Parallel Fluid Dynamics ParCFD2011, 2011

[

Comput. Fluids

80

,

37

–

43

(

2013

)].

Google Scholar

Crossref

22.

J. A.

Sillero

,

J.

Jiménez

, and

R. D.

Moser

, “

One-point statistics for turbulent wall-bounded flows at Reynolds numbers up to

δ^{+}

= 2000

,”

Phys. Fluids

25

,

105102

(

2013

).

https://doi.org/10.1063/1.4823831

Google Scholar

Crossref

23.

G.

Eitel-Amor

,

O.

Ramis

, and

P.

Schlatter

, “

Simulation and validation of a spatially evolving turbulent boundary layer up to Re = 8300

,”

Int. J. Heat Fluid Flow

47

,

57

–

69

(

2014

).

https://doi.org/10.1016/j.ijheatfluidflow.2014.02.006

Google Scholar

Crossref

24.

P. S.

Klebanoff

, “

Characteristics of turbulence in a boundary layer with zero pressure gradient

,”

NACA Report No. 1247

,

1955

.

Google Scholar

25.

D. B.

DeGraaff

and

J. K.

Eaton

, “

Reynolds-number scaling of the flat-plate turbulent boundary layer

,”

J. Fluid Mech.

422

,

319

–

346

(

2000

).

https://doi.org/10.1017/S0022112000001713

Google Scholar

Crossref

26.

R.

Baidya

,

J.

Phillip

,

N.

Hutchins

,

J. P.

Monty

, and

I.

Marusic

, “

Distance from the wall scaling of turbulent motions in wall-bounded flows

,”

Phys. Fluids

29

,

020712

(

2017

).

https://doi.org/10.1063/1.4974354

Google Scholar

Crossref

27.

A. E.

Perry

and

I.

Marusic

, “

A wall-wake model for the turbulence structure of boundary layers. Part 2. Further experimental support

,”

J. Fluid Mech.

298

,

361

–

407

(

1995

).

https://doi.org/10.1017/S0022112095003351

Google Scholar

Crossref

28.

M.

Jones

,

I.

Marusic

, and

A. E.

Perry

, “

Evolution and structure of sink-flow turbulent boundary layers

,”

J. Fluid Mech.

428

,

1

–

27

(

2001

).

https://doi.org/10.1017/S0022112000002597

Google Scholar

Crossref

29.

M. V.

Zagarola

and

A. J.

Smits

, “

Mean-flow scaling of turbulent pipe flow

,”

J. Fluid Mech.

373

,

33

–

79

(

1998

).

https://doi.org/10.1017/S0022112098002419

Google Scholar

Crossref

30.

G. B.

Schubauer

and

P. S.

Klebanoff

, “

Contribution on the mechanics of boundary layer transition

,”

NACA Technical Note TN 3489

,

1955

.

Google Scholar

31.

J.

Coupland

, “

ERCOFTAC special interest group on laminar to turbulent transition and retransition T3A and T3B test cases

,”

Technical Report No. A309514

,

1990

.

Google Scholar

32.

U.

Kaynak

,

O.

Bas

,

S. C.

Cakmakcioglu

, and

I. H.

Tuncer

, “

Transition modeling for low to high speed boundary layer flows with CFD applications

,” in

Boundary Layer Flows—Theory, Applications, and Numerical Methods

(

IntechOpen

,

2020

).

Google Scholar

Crossref

33.

T. D.

Economon

,

F.

Palacios

,

S.

Copeland

,

T.

Lukaczyk

, and

J. J.

Alonso

, “

SU2: An open-source suite for multi-physics simulation and design

,”

AIAA J.

54

,

828

–

846

(

2016

).

https://doi.org/10.2514/1.J053813

Google Scholar

Crossref

34.

J. F.

Nash

, “

A note on skin-friction laws for the incompressible turbulent boundary layer

,” Aeronautical Research Council Current Papers No. 862,

1966

.

Google Scholar

35.

H.

Nagib

,

C.

Christophorou

,

J.

Reudi

,

P.

Monkewitz

,

J.

Osterlund

, and

S.

Gravante

, “

Can we ever rely on results from wall-bounded turbulent flows without direct measurements of wall shear stress?

,” AIAA Paper No. 2004-2392,

2004

.

Google Scholar

36.

M.

Polewski

and

P. G. A.

Cizmas

, “

Several cases for the validation of turbulence models implementation

,”

J. Appl. Sci.

11

,

3377

(

2021

).

https://doi.org/10.3390/app11083377

Google Scholar

Crossref

37.

P. R.

Spalart

and

S. R.

Allmaras

, “

A one equation turbulence model for aerodynamic flows

,” AIAA Paper No. 89-0366,

1992

.

Google Scholar

38.

F. H.

Clauser

, “

Turbulent boundary layers in adverse pressure gradients

,”

J. Aeronaut. Sci.

21

,

91

–

108

(

1954

).

https://doi.org/10.2514/8.2938

Google Scholar

Crossref

39.

M.

Van Dyke

,

Perturbation Methods in Fluid Mechanics

(

Academic Press

,

1964

).

Google Scholar

40.

C. L.

Ladson

, “

Effects of independent variation of Mach and Reynolds numbers on the low-speed aerodynamic characteristics of the NACA 0012 airfoil section

,” NASA Technical Memorandum No. 4074,

1988

.

Google Scholar

41.

W. J.

McCroskey

, “

A critical assessment of wind tunnel results for the NACA 0012 airfoil

,” NASA Technical Memorandum No. 100019,

1987

.

Google Scholar

42.

D. C.

Jesperson

,

T.

Pulliam

, and

M. L.

Childs

, “

OVERFLOW turbulence modeling resource validation results

,” “

NASA Technical Report No. NAS-2016-01

,

2016

.

Google Scholar

43.

A. M. O.

Smith

and

N.

Gamberoni

, “

Transition, pressure gradient and stability theory

,”

Douglas Aircraft Report No. ES 26388

,

1956

.

Google Scholar

44.

J. L.

van Ingen

, “

A suggested semi-empirical method for the calculation of the boundary layer transition region

,”

Technische Hogeschool Vliegtuigbouwkunde Report No. VTH-74

,

1956

.

Google Scholar

45.

J.

Brasseur

and

T.

Wei

, “

Designing large-eddy simulation of the turbulent boundary layer to capture law-of-the-wall scaling

,”

Phys. Fluids

22

,

021303

(

2010

).

https://doi.org/10.1063/1.3319073

Google Scholar

Crossref

46.

X.

Yang

,

G.

Park

, and

P.

Moin

, “

Log-layer mismatch and modeling of the fluctuating wall stress in wall-modeled large-eddy simulations

,”

Phys. Rev. Fluids

2

,

104601

(

2017

).

https://doi.org/10.1103/PhysRevFluids.2.104601

Google Scholar

Crossref

PubMed

47.

S.

Kawai

and

J.

Larsson

, “

Wall-modeling in large eddy simulation: Length scales, grid resolution, and accuracy

,”

Phys. Fluids

24

,

015105

(

2012

).

https://doi.org/10.1063/1.3678331

Google Scholar

Crossref

48.

S.

Cai

and

P.

Sagaut

, “

Explicit wall models for large eddy simulation

,”

Phys. Fluids

33

,

041703

(

2021

).

https://doi.org/10.1063/5.0048563

Google Scholar

Crossref

2022

Author(s)

A new boundary layer integral method based on the universal velocity profile

I. INTRODUCTION

II. THE UNIVERSAL VELOCITY PROFILE (UVP)

A. The high Reynolds number explicit form of the UVP

B. Optimal parameters

1. The equivalent channel half height thickness

C. Evolution of the UVP from the laminar to the turbulent profile as $R_{τ}$ increases from 0 to $\infty$

D. Friction law comparisons

III. THE VON KÁRMÁN INTEGRAL EQUATION EXPRESSED IN TERMS OF $R_{τ}$

IV. THE EFFECT OF PRESSURE GRADIENT ON THE UVP

A. A modified Clauser parameter

V. LEADING EDGE TREATMENT

A. The Hiemenz solution near the leading edge

B. The UVP near the leading edge

VI. ITERATIVE PROCEDURE AND SOLUTION FOR THE J0012 AIRFOIL AT $R_{chord} = 10^{7}$ USING THE INTEGRAL FORM OF THE UVP

A. Step 1—Initial data, initial calculation of $R_{δ_{1}} (R_{τ})$ and $R_{δ_{2}} (R_{τ})$ ⁠, first integration of the Kármán equation

B. Step 2

C. Step 3—Finish

VII. HIGH REYNOLDS NUMBER

VIII. VISCOUS DRAG COEFFICIENT OF THE JOUKOWSKY AND NACA 0012 AIRFOILS

A. Drag force comparison between the UVP and SU2 NACA0012 results at $R_{chord} = 10^{7}$

IX. CONCLUDING REMARKS

ACKNOWLEDGMENTS

AUTHOR DECLARATIONS

Conflict of Interest

Author Contributions

DATA AVAILABILITY

APPENDIX A: SHAPE FUNCTION COEFFICIENTS

APPENDIX B: UVP INTEGRAL METHOD FLOW CHART

References

Citing articles via

Submit your article

Sign up for alerts

Resources

Explore

pubs.aip.org

Connect with AIP Publishing

A new boundary layer integral method based on the universal velocity profile

I. INTRODUCTION

II. THE UNIVERSAL VELOCITY PROFILE (UVP)

A. The high Reynolds number explicit form of the UVP

B. Optimal parameters

1. The equivalent channel half height thickness

C. Evolution of the UVP from the laminar to the turbulent profile as Rτ increases from 0 to ∞

D. Friction law comparisons

III. THE VON KÁRMÁN INTEGRAL EQUATION EXPRESSED IN TERMS OF Rτ

IV. THE EFFECT OF PRESSURE GRADIENT ON THE UVP

A. A modified Clauser parameter

V. LEADING EDGE TREATMENT

A. The Hiemenz solution near the leading edge

B. The UVP near the leading edge

VI. ITERATIVE PROCEDURE AND SOLUTION FOR THE J0012 AIRFOIL AT Rchord=107 USING THE INTEGRAL FORM OF THE UVP

A. Step 1—Initial data, initial calculation of Rδ1(Rτ) and Rδ2(Rτ)⁠, first integration of the Kármán equation

B. Step 2

C. Step 3—Finish

VII. HIGH REYNOLDS NUMBER

VIII. VISCOUS DRAG COEFFICIENT OF THE JOUKOWSKY AND NACA 0012 AIRFOILS

A. Drag force comparison between the UVP and SU2 NACA0012 results at Rchord=107

IX. CONCLUDING REMARKS

ACKNOWLEDGMENTS

AUTHOR DECLARATIONS

Conflict of Interest

Author Contributions

DATA AVAILABILITY

APPENDIX A: SHAPE FUNCTION COEFFICIENTS

APPENDIX B: UVP INTEGRAL METHOD FLOW CHART

References

Citing articles via

Submit your article

Sign up for alerts

Related Content

Resources

Explore

pubs.aip.org

Connect with AIP Publishing

This Feature Is Available To Subscribers Only

C. Evolution of the UVP from the laminar to the turbulent profile as $R_{τ}$ increases from 0 to $\infty$

III. THE VON KÁRMÁN INTEGRAL EQUATION EXPRESSED IN TERMS OF $R_{τ}$

VI. ITERATIVE PROCEDURE AND SOLUTION FOR THE J0012 AIRFOIL AT $R_{chord} = 10^{7}$ USING THE INTEGRAL FORM OF THE UVP

A. Step 1—Initial data, initial calculation of $R_{δ_{1}} (R_{τ})$ and $R_{δ_{2}} (R_{τ})$ ⁠, first integration of the Kármán equation

A. Drag force comparison between the UVP and SU2 NACA0012 results at $R_{chord} = 10^{7}$