Nomenclature
- a
-
sound speed
- d
-
distance to nearest solid wall
- l
-
characteristic length
- M
-
Mach number
- N
-
number of grid nodes
- p
-
static pressure
-
$Re$
-
Reynolds number based on the characteristic length,
$\rho Ul/\mu $
-
$\tilde R{e_{\theta t}}$
-
transported quantity of the transition onset momentum-thickness Reynolds number
-
$R{e_{\theta c}}$
-
critical momentum-thickness Reynolds number,
$\rho U{\theta _c}/\mu $
-
$R{e_{\theta t}}$
-
transition onset momentum-thickness Reynolds number,
$\rho U{\theta _t}/\mu $
-
$R{e_S}$
-
strain-rate magnitude Reynolds number,
$\rho {d^2}S/\mu $
-
${R_T}$
-
eddy viscosity ratio,
${\mu _t}/\mu $
- s
-
streamwise coordinate
- S
-
strain-rate magnitude,
${(2{S_{ij}}{S_{ij}})^{1/2}}$
-
${S_{ij}}$
-
strain-rate tensor,
$0.5(\partial {u_i}/\partial {x_j} + \partial {u_j}/\partial {x_i})$
-
$\tilde S$
-
modified strain-rate magnitude
- Tu
-
turbulence intensity
-
${u_\tau }$
-
friction velocity,
${({\tau _w}/\rho )^{1/2}}$
-
${u_i}$
-
Cartesian velocity component
- U
-
local velocity magnitude,
${({u^2} + {v^2} + {w^2})^{1/2}}$
-
$y + $
-
non-dimensional wall distance,
$\rho {u_\tau }d/\mu $
Greek symbol
-
$\Omega $
-
vorticity magnitude,
${(2{\Omega _{ij}}{\Omega _{ij}})^{1/2}}$
-
${\Omega _{ij}}$
-
vorticity tensor,
$0.5(\partial {u_i}/\partial {x_j} - \partial {u_j}/\partial {x_i})$
-
$\gamma $
-
transported quantity of intermittency
-
$\rho $
-
density
-
$\mu $
-
molecular viscosity
-
${\mu _t}$
-
eddy viscosity
-
$\tilde \nu $
-
modified eddy viscosity
-
${\lambda _\theta }$
-
pressure gradient parameter,
$(\rho {\theta ^2}/\mu )(dU/ds)$
-
$\theta $
-
momentum thickness
-
${\theta _c}$
-
critical momentum thickness
-
${\theta _t}$
-
transition onset momentum thickness
-
$\phi $
-
smooth maximum/minimum function
-
$\psi $
-
compressibility correction for Tollmien-Schlichting instabilities
-
${\psi _{{\rm{scf}}}}$
-
compressibility correction for stationary crossflow instabilities
-
$\kappa $
-
heat capacity ratio
-
${\tau _w}$
-
wall shear stress,
$\mu {(\partial U/\partial d)_{d = 0}}$
Subscripts
- e
-
boundary-layer edge
- comp
-
compressible
-
$\infty $
-
far field
- scf
-
stationary crossflow
1.0 Introduction
The location of boundary-layer transition from laminar to turbulent flow can have a significant impact on the aerodynamic performance of a transonic wing. An upstream laminar boundary layer directly influences shock location and strength, while laminar separation bubbles can lead to severe drag penalties, stability and control issues, and adverse stall characteristics. Recent studies suggest the application of laminar flow control to large commercial aircraft can reduce aerodynamic drag by approximately 10% [Reference Malik, Crouch, Saric, Lin and Whalen1]. In order to design wings optimised for low drag for commercial transport aircraft, which typically fly at transonic Mach numbers and high Reynolds numbers, it is necessary to be able to accurately and efficiently predict boundary-layer transition in transonic flow regimes.
The primary transition mechanisms for a swept transport-aircraft wing result from Tollmien-Schlichting wave growth leading to natural transition, crossflow instabilities resulting from highly swept wings, concave curvature producing Görtler instabilities, and attachment-line instabilities as a result of large leading-edge radius and sweep at the root of the wing [Reference Reed and Saric2]. While the last two mechanisms can be prevented by appropriate profile design, balancing natural and crossflow-induced transition can be difficult, as favourable pressure gradients used to stabilise streamwise instabilities destabilise crossflow instabilities.
The effects of boundary-layer transition can be included in RANS-based simulations using various transition prediction and modelling techniques, consisting of either local or non-local operations [Reference Arnal, Casalis and Houdeville3–Reference Pasquale, Rona and Garrett5]. Local correlation-based transition models, initially developed by Langtry and Menter [Reference Langtry and Menter6], provide a framework for predicting boundary-layer transition that can be easily integrated in modern highly parallel solution algorithms. However, the empirical correlations of Langtry and Menter (LM2009) [Reference Langtry and Menter6] and Langtry et al. (LM2015) [Reference Langtry, Sengupta, Yeh and Dorgan7] for two-dimensional mechanisms (natural, bypass and separation-induced transition) and stationary crossflow instabilities, respectively, were developed based on results from subsonic experiments and stability analysis. Recent studies have demonstrated that models based on these correlations severely under-predict the extent of laminar flow when applied to transonic test cases [Reference Venkatachari, Paredes, Derlaga, Buning, Choudhari, Li and Chang8, Reference Fehrs9], such as the NASA CRM-NLF [Reference Campbell and Lynde10–Reference Lynde, Campbell and Viken12]. This behaviour was foreseen by Arnal who stated, “simple transition criteria developed for low-speed flows cannot be used with confidence in configurations where compressibility effects become significant” [Reference Arnal13].
The stabilising effect of Mach number was recently observed experimentally by Risuis et al. [Reference Risuis, Costantini, Koch, Hein and Klein14] for high subsonic Mach numbers. However, this phenomenon can be more efficiently studied using stability theory due to the difficulties of obtaining detailed transition measurements at transonic and high Reynolds number flow conditions [Reference Van Driest15]. One of the first investigations into the effects of flow compressibility on the stability of laminar boundary layers was performed by Lees and Lin [Reference Lees and Lin16] who developed results for the generalised inflection point and for the effects of wall heating and cooling in subsonic and transonic flow. More recently, the stabilising effect of flow compressibility on laminar boundary layers was investigated first by Mack [Reference Mack17] and then by Arnal [Reference Arnal13, Reference Arnal18]. By applying the
${e^N}$
method with linear stability theory to adiabatic flat plates, Mack [Reference Mack17] and Arnal [Reference Arnal18] observed that compressibility has a strong stabilising effect in transonic flow, before becoming destabilising from Mach 2 to 3.5, and stabilising again for hypersonic flow. However, applying stability analysis to three-dimensional configurations, such as in the design of commercial transport aircraft, can be challenging: the framework for applying stability analysis is complex and can require significant modifications to the flow solver, while simulations often require intervention from experienced and knowledgeable users.
The non-local AHD transition criterion [Reference Arnal, Habiballah and Coutols19] was developed as a simpler alternative to the
${e^N}$
method [Reference Smith and Gamberoni20, Reference Van Ingen21], based on stability analysis of Falkner-Skan attached self-similar boundary-layer velocity profiles. The criterion was extended to include the stabilising effect of compressibility on Tollmien-Schlichting instabilities based on stability analysis of compressible local similarity solutions first up to Mach 1.6 [Reference Arnal, Houdeville, Seraudie and Vermeersch22] and, more recently, up to Mach 4 [Reference Perraud and Durant23]. However, the AHD transition criterion requires non-local flow quantities, which can complicate its integration in highly parallel and unstructured solution algorithms. Two transport-equation-based transition models have recently been developed that create a local framework for evaluating the AHD criterion [Reference Pascal, Delattre, Deniau and Cliquet24–Reference Ströer, Krimmelbein, Krumbein and Grabe26]. However, both models consist of four transport equations coupled to the two-equation Shear Stress Tensor (SST) turbulence model [Reference Menter27], resulting in six-equation coupled systems. Two one-equation AHD-based transition models, each coupled to the SST turbulence model resulting in three-equation coupled systems, have recently been developed; however, they have not yet been extended to include crossflow instabilities [Reference François, Krumbein, Krimmelbein and Grabe28, Reference Ströer, Krimmelbein, Krumbein and Grabe29].
Crossflow instabilities are dominated by the properties of the inflection point, which are not affected by flow compressibility to the same degree as viscous instabilities, such as Tollmien-Schlichting waves [Reference Dagenhart30, Reference Arnal31]. Malik et al. [Reference Malik, Balakumar and Chang32] developed a compressibility correction for a crossflow Reynolds number criterion. This correction was used by Kroo and Sturdza [Reference Kroo and Sturdza33] to develop a compressible crossflow Reynolds number criterion for the design of laminar supersonic swept wings, which was subsequently used by Lee and Jameson [Reference Lee and Jameson34] to perform aerodynamic shape optimisation of transonic natural-laminar-flow (NLF) wings. More recently, the C1 criterion [Reference Arnal, Habiballah and Coutols19], a crossflow criterion developed based on stability analysis of the solutions of the Falkner-Skan and Cooke equations, was extended to take into account compressibility effects using this correction [Reference Arnal, Houdeville, Seraudie and Vermeersch22]. A similar correction was also used by Xu et al. [Reference Xu, Qiao and Bai35] to extend a local model for stationary crossflow instabilities to transonic flows.
The SA-sLM2015 smooth transition model [Reference Piotrowski and Zingg36] is a local two-equation model which couples the LM2009 and LM2015 empirical correlations with the Spalart-Allmaras (SA) one-equation turbulence model [Reference Spalart and Allmaras37], resulting in a three-equation coupled system. Smooth approximations were introduced to replace discontinuous and stiff functions in the transition model source terms, which improves the iterative convergence of the model and provides continuous gradients that facilitate integration in gradient-based optimisation algorithms [Reference Piotrowski and Zingg36, Reference Piotrowski and Zingg38]. The model provides a less expensive and simpler alternative to non-local and higher-fidelity approaches described above.
The goal of this work is to develop a transition model suitable for use in the design of commercial aircraft. This is achieved by extending the LM2009 [Reference Langtry and Menter6] and LM2015 [Reference Langtry, Sengupta, Yeh and Dorgan7] empirical correlations in the SA-sLM2015 transition model to transonic flight regimes typical of a commercial transport aircraft. A compressibility correction for Tollmien-Schlichting instabilities is developed based on results from the compressible extension to the AHD criterion [Reference Arnal, Houdeville, Seraudie and Vermeersch22, Reference Perraud, Deniau and Casalis39] and the stationary crossflow instability compressibility correction developed by Malik et al. [Reference Malik, Balakumar and Chang32] is implemented through a new crossflow source term function. The transition model is investigated with and without the compressibility corrections through comparison with two- and three-dimensional transonic test cases for which transition onset has been measured.
The work is presented as follows: the development of the Tollmien-Schlichting compressibility correction and the modifications made to the transition model to integrate the stationary crossflow instability compressibility correction are presented in Section 2. The comparisons of the model with and without the compressibility corrections to experimental data for two- and three-dimensional transonic transition test cases are presented in Section 3, along with a discussion of the results. Conclusions are presented in Section 4.
2.0 Compressibility corrections
This section presents the compressibility correction developed for the Tollmien-Schlichting instability empirical correlation,
$\psi $
, and the framework created for applying the crossflow compressibility correction developed by Malik et al. [Reference Malik, Balakumar and Chang32],
${\psi _{{\rm{scf}}}}$
. The SA-sLM2015 transition model equations are presented in Appendix A, while a full description of the model, including validation with two- and three-dimensional subsonic transition test cases, can be found in Piotrowski and Zingg [Reference Piotrowski and Zingg36] and Piotrowski [Reference Piotrowski40]. All equations are presented in non-dimensional form to be consistent with previous work [Reference Piotrowski and Zingg36]. The non-dimensionalisation procedure follows that of the Navier-Stokes equations, which is described by Pulliam and Zingg [Reference Pulliam and Zingg41]. The corrections extend the domain of applicability for these commonly used empirical correlations and do not affect the predictive capabilities of the model in the incompressible flow regime, as both compressibility corrections,
$\psi $
and
${\psi _{{\rm{scf}}}}$
, approach unity for low Mach numbers.
2.1 Tollmien-Schlichting instabilities
In an effort to compare the behaviour of the LM2009 empirical correlation [Reference Langtry and Menter6] to stability analysis, Perraud et al. [Reference Perraud, Deniau and Casalis39] developed a simplified model based on the database of stability analysis results that were used to develop the compressible extension to the AHD criterion [Reference Arnal, Houdeville, Seraudie and Vermeersch22]. The model was developed for Mach numbers up to 1.1 and low turbulence intensities (
$Tu \lt 1\% $
). Their model, which approximates the sensitivity of the linear stability equations to varying turbulence intensities, Tu, pressure gradients,
${\lambda _{\theta ,{\rm{e}}}}$
, and flow compressibility,
${M_{\rm{e}}}$
, is defined as,
\begin{equation}R{e_{\theta t}} = f(Tu,{\lambda _\theta },{M_{\rm{e}}}),\quad {\lambda _{\theta ,{\rm{e}}}} = \frac{{\rho {\theta ^2}}}{\mu }\frac{{d{U_{\rm{e}}}}}{{ds}}R{e_\infty },\quad Tu = T{u_\infty },\end{equation}
\begin{equation}R{e_{\theta t}} = - \left(177M_{\rm{e}}^2 - 22{M_{\rm{e}}} + 210\right)\ln \!\left((7{M_{\rm{e}}} + 4.8)Tu/100\right)\exp \!\left((5{M_{\rm{e}}} + 27){\lambda _{\theta ,{\rm{e}}}}\right),\end{equation}
where
$R{e_\infty }$
is introduced as part of the non-dimensionalisation procedure.
A comparison of the sensitivities of the LM2009 empirical correlation and the stability-based model to pressure gradient and Mach number is illustrated in Fig. 1. The LM2009 empirical correlation [Reference Langtry and Menter6] is based on the empirical correlations by Abu-Ghannam et al. [Reference Abu-Ghannam and Shaw42], which were developed based on experimental data with a focus on higher turbulence intensity conditions. In addition, the AHD criterion is calculated using an averaged pressure gradient parameter that is computed by integrating from the critical point to the location being evaluated, while the LM2009 empirical correlation is calculated using a fully local formulation. These factors can help to explain the different sensitivities to pressure gradient in Fig. 1(a). However, as Fig. 1(b) demonstrates, the LM2009 empirical correlation under-predicts the transition onset momentum-thickness Reynolds number at higher Mach numbers by not including the stabilising effect of flow compressibility.
Figure 1. Sensitivity of the LM2009 empirical correlation [Reference Langtry and Menter6], the stability-based model (Equations (1) and (2)), and stability analysis [Reference Perraud, Deniau and Casalis39] to pressure gradient and Mach number. A higher transition onset momentum-thickness Reynolds number delays boundary-layer transition.
To better understand the stabilising effect of flow compressibility as predicted by linear stability theory, the stability-based model is investigated over a range of boundary-layer edge Mach numbers,
${M_{\rm{e}}}$
, at various turbulence intensities, Fig. 2(a), and pressure gradients, Fig. 2(b). To isolate the effects of flow compressibility, the values of transition onset momentum-thickness Reynolds number predicted by both the stability-based model and stability analysis [Reference Perraud, Deniau and Casalis39] are normalised by the value at a boundary-layer edge Mach number of zero.
Figure 2. Sensitivity of the stability-based model (Equations (1) and (2)) and stability analysis [Reference Perraud, Deniau and Casalis39] to Mach number with varying turbulence intensity and pressure gradient. The results are normalised by the values at a Mach number of zero to isolate the effects of flow compressibility.
The results demonstrate that an increasing Mach number produces a more significant stabilising effect at low turbulence intensities and in regions with less adverse and more favourable pressure gradients (here corresponding to a more positive
${\lambda _{\theta ,{\rm{e}}}}$
). This behaviour is consistent with the analysis by Ströer et al. [Reference Ströer, Krimmelbein, Krumbein and Grabe26], who investigated the compressible AHD criterion on a flat plate with varying pressure gradients. As the Reynolds number increases, transition often occurs in regions of favourable pressure gradient. Therefore, the stabilising effect of flow compressibility is most relevant at the cruise conditions typical of commercial transport aircraft, which is characterised by low turbulence intensity and high Reynolds number flow.
To include the stabilising effect of flow compressibility on the LM2009 empirical correlation for Tollmien-Schlichting instabilities, an initial compressibility correction,
${\psi _{{\rm{init}}}}$
, is introduced. The correction was designed to reproduce the trends produced by the normalised stability analysis-based results (Fig. 2). Specifically, it was developed to approximate the stabilising effect of flow compressibility on Tollmien-Schlichting instabilities, and its relationship with the flow disturbance environment, including the effects of turbulence intensity (Fig. 2(a)) and pressure gradient (Fig. 2(b)). The correction is evaluated using the Mach number at the edge of the boundary layer,
${M_{\rm{e}}}$
, the freestream turbulence intensity,
$T{u_\infty }$
, and a local approximation to the boundary-layer edge pressure gradient,
${\lambda _{\theta ,{\rm{e}}}}$
, as follows:
\begin{equation}{\psi _{{\rm{init}}}} = \left({a_1}M_{\rm{e}}^2 + {a_2}{M_{\rm{e}}} + {a_3}\right)\exp \!\left[({b_1}{\lambda _{\theta ,{\rm{e}}}}){M_{\rm{e}}}\right]\exp \!\left[\left({c_1} + {c_2}\sqrt {T{u_\infty }} \right){M_{\rm{e}}}\right],\end{equation}
\begin{equation}{a_1} = 0.44,\;{a_2} = - 0.38,\;{a_3} = 1.00;\quad {b_1} = 5.00;\quad {c_1} = 0.41,\;{c_2} = - 0.27.\end{equation}
It is compared with the normalised stability analysis-based results in Fig. 3, with the results presented in grey. The results demonstrate that the correction does a reasonable job of reproducing the trends predicted by stability analysis.
Figure 3. Sensitivity of the stability-based model (blue) (Equations (1) and (2)), stability analysis [Reference Perraud, Deniau and Casalis39], and the initial compressibility correction,
${\psi _{{\rm{init}}}}$
(grey) (Equations (3) and (4)), to Mach number with varying turbulence intensity and pressure gradient. The stability analysis-based results are normalised by the values at a Mach number of zero to isolate the effects of flow compressibility.
The correction was initially applied to the local transition onset momentum-thickness Reynolds number,
$R{e_{\theta t}}$
, which is convected and diffused into the boundary layer by the
$\tilde R{e_{\theta t}}$
transport equation. However, simulations using the test cases presented in Section 3 demonstrated that this approach significantly over-predicted the laminar extent of the boundary layer for the CAST10-2 and VA-2 aerofoils. The rapid acceleration at the upper surface leading edge produced large upstream values of the compressibility correction, and therefore
$\tilde R{e_{\theta t}}$
. These large values convect downstream in the boundary layer and delay transition. This behaviour is not reproduced by methods based on the AHD criterion where the pressure gradient parameter is integrated from the critical point.
Applying the correction to the critical momentum-thickness Reynolds number,
$R{e_{\theta c}}$
(Equation (A.33)), produced transition locations that agree better with the CAST10-2 and VA-2 experimental results by avoiding this behaviour, while still agreeing well with the CRM-NLF test case. This approach is similar to scaling the local value of the transported momentum-thickness Reynolds number,
$\tilde R{e_{\theta t}}$
, as
$R{e_{\theta c}}$
is a near-linear function of
$\tilde R{e_{\theta t}}$
. A consequence of this approach is that the
${F_{{\rm{length}}}}$
empirical correlation (Equation (A.34)) is not scaled by the compressibility correction. However, the stability-based model and the AHD criterion do not predict the effects of compressibility on the transition region length. For example, the AHD-based transition models developed by Pascal et al. [Reference Pascal, Delattre, Deniau and Cliquet24] and Ströer et al. [Reference Ströer, Krimmelbein, Krumbein and Grabe25, Reference Ströer, Krimmelbein, Krumbein and Grabe26] rely on user-specified inputs for the length of the transition region. Further work is required to develop a correlation for the transition length in compressible flows.
Using this approach, simulations revealed that the smooth transition model with the initial compressibility correction,
${\psi _{{\rm{init}}}}$
, still over-predicted the extent of the laminar boundary layer for simulations of the VA-2 aerofoil at lower angles of attack, and appeared to over-estimate the effect of favourable and adverse pressure gradients on the stabilising effect of compressibility. This could again be related to differences between the local model and the AHD criterion, such as the use of a fully local pressure gradient parameter versus an integrated value. Therefore, the correction was modified to improve agreement with the experimental results. Specifically, the
${a_1}$
and
${b_1}$
constants were reduced relative to the initial correction in order to reduce the sensitivity to Mach number and pressure gradient, respectively. Care was taken to ensure that these modifications did not have a significant adverse effect on the simulated CAST10-2 and CRM-NLF transition fronts. The modified compressibility correction,
$\psi $
, produces transition locations that agree well with experimental results for the test cases and flow conditions investigated (
$0.71 \le M \le 0.856$
,
$2 \times {10^6} \le Re \le 15 \times {10^6}$
), and is given by,
\begin{equation}R{e_{\theta c,{\rm{comp}}}} = \psi R{e_{\theta c}},\end{equation}
\begin{equation}\psi = \left({a_1}M_{\rm{e}}^2 + {a_2}{M_{\rm{e}}} + {a_3}\right)\exp \!\left[\left({b_1}{\lambda _{\theta ,{\rm{e}}}}\right){M_{\rm{e}}}\right]\exp \!\left[\left({c_1} + {c_2}\sqrt {T{u_\infty }} \right){M_{\rm{e}}}\right],\end{equation}
\begin{equation}{a_1} = 0.34,\;{a_2} = - 0.38,\;{a_3} = 1.00;\quad {b_1} = 3.00;\quad {c_1} = 0.41,\;{c_2} = - 0.27.\end{equation}
The correction is compared with the initial correction,
${\psi _{{\rm{init}}}}$
(Equations (3) and (4)), in Fig. 4.
The compressibility corrections are evaluated using the Mach number at the edge of the boundary layer,
${M_{\rm{e}}}$
, which can be approximated locally using isentropic relations and the compressible Bernoulli equation as follows [Reference Coder and Maughmer43]:
\begin{equation}{M_{\rm{e}}} = \frac{{{U_{\rm{e}}}}}{{{a_{\rm{e}}}}},\quad {a_{\rm{e}}} = \sqrt {\frac{{\kappa {p_\infty }}}{{{\rho _\infty }}}{{\left(\frac{p}{{{p_\infty }}}\right)}^{1 - \frac{1}{\kappa }}}} ,\end{equation}
\begin{equation}{U_{\rm{e}}} = \sqrt {U_\infty ^2 + \frac{{2\kappa }}{{\kappa - 1}}\frac{{{p_\infty }}}{{{\rho _\infty }}}\left[1 - {{\left(\frac{p}{{{p_\infty }}}\right)}^{1 - \frac{1}{\kappa }}}\right]} ,\end{equation}
where
$\kappa $
is the heat capacity ratio.
The pressure gradient parameter,
${\lambda _\theta }$
(Equations (A.13)–(A.17)), used in the LM2009 empirical correlation is based on the local streamwise velocity gradient formed using Cartesian velocity gradients, which is not valid in the boundary layer. To avoid problems with this formulation, the LM2009 transition model introduces the
${F_{\theta t}}$
function (Equation (A.4)) in the
${P_{\theta t}}$
source term (Equation (A.2)) to disable the Tollmien-Schlichting empirical correlation inside the boundary layer. However, the compressibility corrections are most active in the middle of the boundary layer where the transition onset process begins. To address this, the pressure gradient parameter used to evaluate the compressibility correction is calculated using a local approximation to the boundary-layer edge pressure gradient,
${\lambda _{\theta ,{\rm{e}}}}$
, adopted from Grabe et al. [Reference Grabe, Shengyang and Krumbein44, Reference Grabe, Shengyang and Krumbein45], which is approximated using Cartesian pressure gradients and isentropic relations as follows:
\begin{equation}{\lambda _{\theta ,{\rm{e}}}} = \frac{{\rho {\theta ^2}}}{\mu }\frac{{d{U_{\rm{e}}}}}{{ds}}R{e_\infty },\end{equation}
\begin{equation}\frac{{d{U_{\rm{e}}}}}{{ds}} = \left[\left(\frac{u}{U}\right)\frac{{d{U_{\rm{e}}}}}{{dx}} + \left(\frac{v}{U}\right)\frac{{d{U_{\rm{e}}}}}{{dy}} + \left(\frac{w}{U}\right)\frac{{d{U_{\rm{e}}}}}{{dz}}\right],\end{equation}
\begin{equation}\frac{{d{U_{\rm{e}}}}}{{dx}} = - \frac{1}{{{\rho _\infty }{U_{\rm{e}}}}}{\left(\frac{p}{{{p_\infty }}}\right)^{ - \frac{1}{\kappa }}}\frac{{dp}}{{dx}},\end{equation}
\begin{equation}\frac{{d{U_{\rm{e}}}}}{{dy}} = - \frac{1}{{{\rho _\infty }{U_{\rm{e}}}}}{\left(\frac{p}{{{p_\infty }}}\right)^{ - \frac{1}{\kappa }}}\frac{{dp}}{{dy}},\end{equation}
\begin{equation}\frac{{d{U_{\rm{e}}}}}{{dz}} = - \frac{1}{{{\rho _\infty }{U_{\rm{e}}}}}{\left(\frac{p}{{{p_\infty }}}\right)^{ - \frac{1}{\kappa }}}\frac{{dp}}{{dz}}.\end{equation}
Similar to the original pressure gradient parameter,
${\lambda _\theta }$
(Equations (A.13)–(A.17)),
${\lambda _{\theta ,{\rm{e}}}}$
is not Galilean invariant.
It is important to note that the boundary-layer edge pressure gradient parameter,
${\lambda _{\theta ,{\rm{e}}}}$
, defined above is not used to calculate the LM2009 empirical correlation (Equations (A.8)–(A.12)). Although the two pressure gradient formulations produce similar results for favourable and adverse pressure gradients, the original, velocity-gradient-based pressure gradient parameter,
${\lambda _\theta }$
, is used to be consistent with the Langtry-Menter model [Reference Langtry and Menter6], as it provides better agreement with experimental data for zero-pressure-gradient flat plate cases. This is likely due to the calibration of the model. For numerical robustness,
${\lambda _{\theta ,{\rm{e}}}}$
is limited using the bounds introduced by Langtry and Menter [Reference Langtry and Menter6] for
${\lambda _\theta }$
(
$ - 0.1 \le {\lambda _{\theta ,{\rm{e}}}} \le 0.1$
).
2.2 Stationary crossflow instabilities
The stationary crossflow transition onset momentum-thickness Reynolds number correlation,
$R{e_{\theta t,{\rm{scf}}}}$
, produced by the LM2015 empirical correlation (Equations (A.18)–(A.23)) was originally implemented using an additional source term in the
$\tilde R{e_{\theta t}}$
transport equation,
${D_{{\rm{scf}}}}$
[Reference Langtry, Sengupta, Yeh and Dorgan7]. The LM2009 empirical correlation influences
$\tilde R{e_{\theta t}}$
in the freestream, which is then convected and diffused into the boundary layer where the stationary crossflow source term is active. The combined effect of these correlations is realised through the
${F_{{\rm{onset}}}}$
(Equation (A.31)) function in the intermittency source terms through the critical momentum-thickness Reynolds number,
$R{e_{\theta c}}$
.
In order to include the stabilising effects of compressibility on each set of correlations, and subsequent Reynolds numbers,
$R{e_{\theta t}}$
and
$R{e_{\theta t,{\rm{scf}}}}$
, separately,
${F_{{\rm{onset}}}}$
is separated into two functions, one for two-dimensional mechanisms and one for stationary crossflow instabilities. A similar approach to that developed by Grabe et al. [Reference Grabe, Shengyang and Krumbein45] and Carnes and Coder [Reference Carnes and Coder46] is adopted. A new onset function,
${F_{{\rm{onset,scf}}}}$
, is introduced for the intermittency transport equation, replacing the
${D_{{\rm{scf}}}}$
source term in the
$\tilde R{e_{\theta t}}$
transport equation. The
${F_{{\rm{onset,scf}}}}$
function was developed using the same approach as the original
${F_{{\rm{onset}}}}$
function, and therefore the smooth approximation introduced in previous work [Reference Piotrowski and Zingg36] (Equation (A.31)), leading to the following:
\begin{equation}{F_{{\rm{onset,scf}}}} = \frac{{\tanh (6({F_{{\rm{onset,scf,1}}}} - 1.35)) + 1}}{2},\end{equation}
\begin{equation}{F_{{\rm{onset,scf,1}}}} = \sqrt {{{\left(\frac{{R{e_S}}}{{2.60R{e_{\theta c,{\rm{scf}}}}}}\right)}^2} + {{({R_T})}^2}} ,\quad R{e_{\theta c,{\rm{scf}}}} = 0.623R{e_{\theta t,{\rm{scf}}}},\end{equation}
where
$R{e_S}$
and
${R_T}$
represent the strain-rate magnitude Reynolds number and eddy viscosity ratio, respectively. The transition process is triggered when
${F_{{\rm{onset,scf,1}}}}$
exceeds unity and is sustained as the eddy viscosity ratio,
${R_T}$
, increases. The stationary crossflow critical momentum-thickness Reynolds number,
$R{e_{\theta c,{\rm{scf}}}}$
, correlation in the new
${F_{{\rm{onset,scf}}}}$
function was calibrated using the NLF2-0415 Infinite Swept Wing [Reference Dagenhart and Saric47, Reference Radeztsky, Reibert and Saric48] transition test case. The value of 0.623 applied to
$R{e_{\theta t,{\rm{scf}}}}$
in order to initiate the transition process upstream of the predicted transition location closely resembles a linear fit of the smooth
$R{e_{\theta c}}$
function (Equation (A.33)) and is consistent with the value of 0.62 used by Medida [Reference Medida49].
The
${F_{{\rm{onset,scf}}}}$
function is incorporated in the intermittency source terms (Equations (A.29) and (A.30)) by taking a smooth maximum of the two onset functions:
\begin{equation}{F_{{\rm{onset}}}} \leftarrow {\phi _{300}}({F_{{\rm{onset}}}},{F_{{\rm{onset,scf}}}}),\end{equation}
where
${\phi _{ \pm 300}}$
represent smooth approximations of minimum/maximum operators introduced in previous work [Reference Piotrowski and Zingg36]. A positive value of the smoothing parameter reproduces a maximum, while a negative value approximates a minimum, and a larger magnitude value more closely reproduces the minimum/maximum operator. A value of
$ \pm 300$
is used for the results presented. The vorticity limiting procedure in the intermittency source terms introduced in previous work to improve the iterative convergence of the model [Reference Piotrowski and Zingg36] is re-written in Appendix A to clarify the implementation in dimensional solvers. This formulation is equivalent to that presented by Piotrowski and Zingg [Reference Piotrowski and Zingg36] which was non-dimensionalised following the procedure for the Navier-Stokes equations(41), where velocity and distance are normalised by the sound speed and reference length, l, respectively. The new crossflow source term formulation was validated using the TU Braunschweig Sickle Wing transition test case [Reference Petzold and Radespiel50], with the results presented by Piotrowski [Reference Piotrowski40].
Simulations of transonic swept wings revealed that both crossflow source-term approaches,
${D_{{\rm{scf}}}}$
and
${F_{{\rm{onset,scf}}}}$
, produce a transition front that is dependent on the initialisation of the flow field. A solution initialised with far-field conditions converges to a transition front upstream of the converged transition front produced when initialised with a converged solution obtained with the crossflow correlation inactive. It is important to note that this behaviour does not appear in subsonic transition test cases with crossflow, such as the NASA NLF2-0415 Infinite Swept Wing [Reference Dagenhart and Saric47, Reference Radeztsky, Reibert and Saric48] and the TU Braunschweig Sickle Wing [Reference Petzold and Radespiel50]. To prevent these inaccurate upstream transition fronts, for transonic cases the crossflow source term is activated after the total residual drops several orders of magnitude without the crossflow source term,
${F_{{\rm{onset,scf}}}}$
, active. A relative residual drop tolerance of five orders of magnitude for the total residual was found to be sufficient and is used for the results presented in the current work. This strategy produces the same results as initialising with a converged simulation performed without the crossflow correlations active.
The crossflow Reynolds number compressibility correction developed by Malik et al. [Reference Malik, Balakumar and Chang32],
${\psi _{{\rm{scf}}}}$
, is applied by scaling the incompressible stationary crossflow transition onset momentum-thickness Reynolds number,
$R{e_{\theta t,{\rm{scf}}}}$
, in the
${F_{{\rm{onset,scf}}}}$
function as follows:
\begin{equation}R{e_{\theta t,{\rm{scf}},{\rm{comp}}}} = {\psi _{{\rm{scf}}}}R{e_{\theta t,{\rm{scf}}}},\end{equation}
\begin{equation}{\psi _{{\rm{scf}}}} = 1 + \frac{{\kappa - 1}}{2}M_{\rm{e}}^2.\end{equation}
The compressibility corrections for Tollmien-Schlichting and stationary crossflow instabilities,
$\psi $
and
${\psi _{{\rm{scf}}}}$
, respectively, are presented in Fig. 5 for
$T{u_\infty } = 0.05\% $
and
${\lambda _{\theta ,{\rm{e}}}} = 0.04$
. As expected, the compressibility correction for Tollmien-Schlichting instabilities produces a stronger stabilising effect than the crossflow instability compressibility correction, especially in a favourable flow environment.
Figure 5. Effects of Mach number on the Tollmien-Schlichting and stationary crossflow instability compressibility corrections,
$\psi $
and
${\psi _{{\rm{scf}}}}$
, respectively.
3.0 Results
The transition model with and without the compressibility corrections, SA-sLM2015cc and SA-sLM2015, respectively, is applied to two- and three-dimensional transonic transition test cases. These cases consist of the CAST10-2 and VA-2 aerofoil test cases and the NASA CRM-NLF wing-body geometry, with the Reynolds number increasing for each successive test case. As the Reynolds number increases and transition occurs in regions with less adverse and more favourable pressure gradients, the stabilising effect of flow compressibility on Tollmien-Schlichting instabilities is expected to become more significant (see Fig. 3(b)).
A second-order discretisation is used for the mean-flow equations using a matrix-based dissipation model [Reference Swanson and Turkel51], with a first-order upwind scheme applied to the turbulence and transition model convective terms. Care is taken for each test case to ensure that all results are sufficiently iteratively and grid converged to produce the low levels of numerical error required to investigate the modelling error and validity of the transition model variants. Machine-zero residual convergence for the total residual is achieved for the cases presented unless stated otherwise. More details on the solution strategy and the iterative and grid convergence of the smooth transition model with the compressibility corrections can be found in previous work [Reference Piotrowski40, Reference Piotrowski and Zingg52].
3.1 CAST10-2 aerofoil
The CAST10-2 aerofoil [Reference Stanewsky and Zimmer53] was investigated by Hebler et al. [Reference Hebler, Schojda and Mai54] in the Transonic Wind Tunnel Göttingen (DNW-TWG) at a Mach number of 0.75, Reynolds number of
$2 \times {10^6}$
, and angles of attack ranging from
$ - 0.79$
to 1.41 degrees. The aerofoil has also been investigated numerically by several researchers [Reference Arnal, Casalis and Houdeville3, Reference Ströer, Krimmelbein, Krumbein and Grabe26, Reference Fehrs, van Rooij and Nitzsche55, Reference Fehrs, Helm and Kaiser56]. The turbulence intensity of the DNW-TWG wind tunnel was not reported in the experiment. However, Fehrs et al. [Reference Fehrs, van Rooij and Nitzsche55] suggest a range of turbulence intensities for the wind tunnel of between 0.25–0.40%, while Ströer et al. [Reference Ströer, Krimmelbein, Krumbein and Grabe26] assumed a value of 0.22%. Hebler et al. simulated the aerofoil with the 2D coupled Euler-boundary-layer solver MSES [Reference Drela and Giles57], which is capable of predicting transition using a database
${e^N}$
method. The results demonstrate that an N-factor of 6 produced a transition front that agrees well with the experiment, which, using Mack’s relation [Reference Mack58], corresponds to a turbulence intensity of approximately 0.24% [Reference Hebler, Schojda and Mai54].
The blunt trailing-edge aerofoil coordinates are provided by Dress et al. [Reference Dress, Johnson, McGuire, Stanewsky and Ray59]. As discussed by Fehrs et al. [Reference Fehrs, van Rooij and Nitzsche55], the trailing edge of the aerofoil is located below the reference line, with the aerofoil coordinates exhibiting a
${0.88^ \circ }$
angle-of-attack. All angles of attack are provided with respect to this reference system. A Mach number shift of
$ - 0.01$
and angle-of-attack shift of
$ - 0.30$
degrees were identified by Hebler et al. [Reference Hebler, Schojda and Mai54] in order to match the experimental data with free-air simulations. These shifts are applied in the current work, along with an assumed turbulence intensity of 0.25% based on the lower bound estimated by Fehrs et al. [Reference Fehrs, van Rooij and Nitzsche55]. However, the results are plotted and referenced using the uncorrected angles of attack in order to be consistent with the results presented by Hebler et al. [Reference Hebler, Schojda and Mai54].
Four grid levels are investigated using structured multi-block grids with the characteristics presented in Table 1. A grid-refinement study is presented using the transition model both with and without the compressibility corrections at the
$ - 0.39$
and 0.82 degree angles of attack, with the lift and drag grid convergence illustrated in Fig. 6. The largest difference in the
${C_l}$
and
${C_d}$
between the second-finest and finest grid levels is less than
$0.30\% $
for both transition model variants at each angle-of-attack.
Table 1. CAST10-2 structured O-grid dimensions
Figure 6. Grid-convergence results for the CAST10-2 aerofoil simulations at
$ - 0.39$
and
$ 0.82 $
degrees angle-of-attack and
$M = 0.74$
,
$Re = 2 \times {10^6}$
, and
$Tu = 0.25\% $
.
An angle-of-attack sweep is presented for the three finest grid levels both fully turbulent using the SA turbulence model and with free transition using the transition model with and without the compressibility corrections. The pressure and upper-surface skin friction coefficient profiles at angles of attack of
$ - 0.39$
, 0.05, and 0.82 degrees are overlaid with the pressure profiles from the experiment in Fig. 7.
Figure 7. Pressure and upper-surface skin friction coefficient profiles for the CAST10-2 aerofoil produced at
$M = 0.74$
,
$Re = 2 \times {10^6}$
, and
$Tu = 0.25\% $
at three angles of attack overlaid with the pressure profiles from the experiment [Reference Hebler, Schojda and Mai54].
The results demonstrate that the pressure and skin friction profiles are sufficiently grid converged. The free-transition pressure coefficient profiles demonstrate good agreement with the experiment at
$ - 0.39$
and 0.82 degrees angle-of-attack. However, for the 0.05 degree angle-of-attack case, in the nonlinear lift region (see Fig. 8), the fully turbulent pressure profiles appear to provide better agreement with the experiment. In general, the wind-tunnel corrections suggested by Hebler et al. [Reference Hebler, Schojda and Mai54] appear to reproduce the wind-tunnel environment well.
Figure 8. Upper-surface transition locations (dashed) and lift curves (solid) for the CAST10-2 aerofoil obtained on the L3 grid at
$M = 0.74$
,
$Re = 2 \times {10^6}$
, and
$Tu = 0.25\% $
over a range of angles of attack compared with the results from the experiment [Reference Hebler, Schojda and Mai54].
As expected, the Tollmien-Schlichting instability compressibility correction does not have a significant effect on the transition onset locations for this low Reynolds number case, as demonstrated by the upper-surface skin friction coefficient profiles, as transition primarily occurs either due to a laminar separation bubble or a strong adverse pressure gradient. Figure 3(b) demonstrates that flow compressibility is not predicted to have a strong stabilising effect in these flow environments. This behaviour is also demonstrated in the work by Ströer et al. [Reference Ströer, Krimmelbein, Krumbein and Grabe26], who simulated the CAST10-2 aerofoil at similar flight conditions using a transport-equation-based formulation of the compressible AHD criterion.
The predicted transition locations on the upper surface of the CAST10-2 aerofoil and lift curves, both obtained on the L3 grid, are presented in Fig. 8 and are compared with the experimental results from Hebler et al. [Reference Hebler, Schojda and Mai54]. The results demonstrate that the free-transition simulations do a reasonable job of predicting the nonlinear transition and lift curves produced by the experiment. However, at angles of attack between
$ - 0.17$
and 0.64 degrees, the free-transition simulations under-predict the drop in the lift curve by over-predicting the extent of the laminar boundary layer relative to experimental data. This could be a result of the uncertainty in the turbulence intensity for the wind tunnel at these conditions.
As illustrated in Fig. 8, the CAST10-2 aerofoil produces highly nonlinear transition and lift curves, with transition moving forward significantly on the upper surface of the aerofoil over a small range of angles of attack. For the 0.05 degree angle-of-attack case, steady-state simulations using the transition model with the compressibility corrections produced an oscillatory force history, with the upper-surface transition front moving upstream and downstream in a periodic manner. The red areas in Fig. 7 represent the regions between the pressure and skin friction coefficient profiles for the most upstream and downstream transition locations on each grid level, while the error bars in Fig. 8 represent the range in the transition location and lift coefficient over one cycle of this oscillatory behaviour.
The results in Figs 7 and 8 demonstrate the importance of including the relationship between pressure gradient and flow compressibility in the Tollmien-Schlichting instability compressibility correction. Due to the low Reynolds number, transition occurs in regions of adverse pressure gradient where stability analysis does not predict the effects of flow compressibility to be significant. A compressibility correction that does not take into account the effects of pressure gradient, such as that investigated by Venkatachari et al. [Reference Venkatachari, Paredes, Choudhari, Li and Chang60], will over-predict the stabilising effect of flow compressibility at these flow conditions, resulting in predicted transition and lift curves that do not capture the nonlinear trends produced by Hebler et al. [Reference Hebler, Schojda and Mai54]. As the Reynolds number increases, the stabilising effect of flow compressibility becomes more significant. This is demonstrated in the following test cases.
3.2 VA-2 supercritical aerofoil
The VA-2 supercritical aerofoil was recently investigated in the DNW-TWG wind tunnel by Costantini et al. [Reference Costantini, Lee, Nonomura, Asai and Klein61]. The aerofoil was examined at a Mach number of 0.72, Reynolds number of
$10 \times {10^6}$
, and angles of attack from
$ - 0.4$
to 2.0 degrees. The upper-surface skin friction coefficient distributions and transition locations were determined using a global luminescent oil-film skin friction field estimation method (GLOFSFE), which measures the development of the thickness of an oil film, the distribution of which can be used to calculate the skin friction, based on its luminescent intensity. The pressure distributions used for comparison in the current work were obtained before applying the oil film.
Two transition locations were identified in the experiment: the transition onset location, identified where the skin friction coefficient increased beyond a value of
$4 \times {10^{ - 4}}$
, and the end of the transition region, defined as the location where the skin friction coefficient reached a value of
$3 \times {10^{ - 3}}$
. While the authors state that the end of the transition region can only be considered as qualitative due to uncertainties in the turbulent skin friction estimation originating from thick oil-film distributions in these regions, quantitative results were obtained for transition onset locations with lower uncertainty. Transition locations were also determined using an automatic temperature-sensitive paint (TSP) transition detection method, where the recorded transition locations are expected to be between the transition onset and end locations [Reference Costantini, Lee, Nonomura, Asai and Klein61, Reference Costantini, Henne, Risius and Klein62].
The Mach number and angle-of-attack shifts identified by Hebler et al. for the DNW-TWG wind tunnel [Reference Hebler, Schojda and Mai54],
$ - 0.01$
Mach and
$ - 0.30$
degrees, respectively, are applied for the VA-2 simulations. Similar to the previous case, the results are presented and referenced using the uncorrected angles of attack. The turbulence intensity for the VA-2 experimental investigations in the DNW-TWG was not provided. Therefore, the turbulence intensity is assumed to be consistent with the CAST10-2 DNW-TWG conditions, with a value of 0.25% used in the current work.
Four grid levels are investigated using structured multi-block grids with the characteristics presented in Table 2. Grid-refinement studies are presented both with and without the compressibility corrections at angles of attack of
$ - 0.40$
and 1.80 degrees in Fig. 9. The difference in the drag coefficient between the finest grid level and the grid-converged value was evaluated using Richardson extrapolation for each transition model variant at both angles of attack. For the
$ - 0.40$
degree case, the differences are 0.18% and 0.80% for the transition model with and without the compressibility corrections, respectively, while for the 1.80 degree case the differences are 0.24% and 0.37%. The iterative and grid convergence for this test case is presented in more detail in previous work [Reference Piotrowski and Zingg52].
Table 2. VA-2 structured O-grid dimensions
Figure 9. Grid-convergence results for the VA-2 aerofoil simulations at
$ - 0.40$
and 1.80 degrees angle-of-attack and
$M = 0.71$
,
$Re = 10 \times {10^6}$
, and
$Tu = 0.25\% $
.
An angle-of-attack sweep is presented for the three finest grid levels both with and without the compressibility corrections. The pressure and upper-surface skin friction coefficient profiles are overlaid with the profiles from the experiment [Reference Costantini, Lee, Nonomura, Asai and Klein61] in Fig. 10. For the lower angles of attack,
$ - 0.40$
degrees to 1.20 degrees, there is good agreement between the simulated pressure coefficient profiles and the profiles from the experiment. However, Costantini et al. observed that the algorithm used to control the DNW-TWG adaptive wind tunnel walls, which are used to reduce the effect of the walls on the pressure distributions, failed to converge at angles of attack above 1.20 degrees [Reference Costantini, Lee, Nonomura, Asai and Klein61], which helps to explain the difference in the simulated and experimental pressure profiles observed at angles of attack of 1.50 and 2.00 degrees. The results confirm that the Mach number and angle-of-attack corrections applied in the simulations accurately reproduce the wind-tunnel environment.
Figure 10. Pressure and upper-surface skin friction coefficient profiles for the VA-2 aerofoil produced at
$M = 0.71$
,
$Re = 10 \times {10^6}$
, and
$Tu = 0.25\% $
over a range of angles of attack overlaid with the results from the experiment [Reference Costantini, Lee, Nonomura, Asai and Klein61].
The simulated upper-surface skin friction coefficient profiles are plotted and compared with the profiles obtained using the GLOFSFE method [Reference Costantini, Lee, Nonomura, Asai and Klein61]. The experimental skin friction coefficient profiles are plotted both with the raw data, which is clipped at a value of
$6 \times {10^{ - 3}}$
, and with the data scaled by a factor of 17. Costantini et al. note that the oil-film thickness in the turbulent boundary-layer regions was larger than the height of the viscous sublayer, which produced a hydraulically rough surface that increased the skin friction relative to the clean configuration [Reference Costantini, Lee, Nonomura, Asai and Klein61]. The clipped experimental data provides a quantitative comparison of the laminar skin friction coefficient profiles and the transition onset locations, which were not as affected by the oil film due to the thinner oil-film thickness in these regions, while the scaled data allows for a qualitative comparison of the skin friction coefficient profiles in the transition and turbulent regions. The clipped experimental skin friction coefficient profiles agree well with the simulated profiles in the laminar boundary layer.
At the lower angles of attack, below 0.80 degrees, there is a weak favourable pressure gradient on the upper surface of the aerofoil that develops into a weak adverse pressure gradient as the angle-of-attack increases. In these regions, where Fig. 3(b) demonstrates that flow compressibility is expected to be stabilising, the compressibility correction pushes transition aft to better agree with the scaled GLOFSFE-estimated skin friction coefficient profiles [Reference Costantini, Lee, Nonomura, Asai and Klein61]. As the angle-of-attack increases and the flow decelerates on the upper surface, the stabilising effect of flow compressibility is reduced, and the transition front moves upstream to match both the transition model without the compressibility corrections and the experimental data. Above an angle-of-attack of 1.20 degrees, the upper-surface pressure coefficient profiles plateau and the strong adverse pressure gradient moves aft in the form of a shock wave. The stabilising effect of compressibility increases with increasing angle-of-attack as the adverse pressure gradient weakens, pushing transition downstream. Although wind-tunnel effects at the 1.50 and 2.00 degree angles of attack appear to be significant, good agreement is achieved for both the pressure and skin friction coefficient profiles at 1.80 degrees.
The predicted transition locations on the L3 grid are compared with values from the experiments performed using the GLOFSFE and TSP methods [Reference Costantini, Lee, Nonomura, Asai and Klein61] in Fig. 11. The results demonstrate that the transition model without the compressibility corrections produces a flatter transition curve that is upstream of the experimental GLOFSFE-estimated transition onset locations for all angles of attack except
$ - 0.40$
and 1.20 degrees. However, at
$ - 0.40$
degrees the transition location produced without the compressibility corrections is upstream of the location produced by the TSP method. The TSP transition locations should lie between the GLOFSFE-estimated transition onset and end locations [Reference Costantini, Lee, Nonomura, Asai and Klein61], such as for the
$ - 0.40$
, 0.80, and 1.80 degrees angle-of-attack. There is increased uncertainty where this is not the case, such as at the 0.00 and 0.40 degree angles of attack. In addition, at 1.20 and 1.50 degrees the experimental upper-surface shocks are upstream of the simulated results, which explains the upstream transition locations. The Tollmien-Schlichting instability compressibility correction delays the predicted transition onset locations to agree best with the GLOFSFE-estimated end of transition region locations, except for 1.80 and 2.00 degrees angle-of-attack, where they are closer to the GLOFSE-estimated transition onset locations.
Figure 11. Upper-surface transition locations for the VA-2 aerofoil obtained on the L3 grid at
$M = 0.71$
,
$Re = 10 \times {10^6}$
, and
$Tu = 0.25\% $
over a range of angles of attack compared with results from the experiments [Reference Costantini, Lee, Nonomura, Asai and Klein61].
3.3 NASA CRM-NLF wing-body geometry
The NASA CRM-NLF configuration was recently investigated as part of the First AIAA Transition Modelling and Prediction WorkshopFootnote
*
. Transition visualisations, pressure coefficient profiles, and integrated forces were provided from the experiment over a range of angles of attack [Reference Lynde, Campbell, Rivers, Viken, Chan, Watkins and Goodliff11, Reference Lynde, Campbell and Viken12]. The test conditions for the data provided are presented in Table 3. The turbulence intensity in the NASA Langley National Transonic Facility (NTF) wind tunnel was determined through comparison with results obtained using stability analysis [Reference Crouch, Sutanto, Witkowski, Watkins, Rivers and Campbell63]. The Tollmien-Schlichting N-factor was found to vary from 4 to 8 [Reference Crouch, Sutanto, Witkowski, Watkins, Rivers and Campbell63], with a critical N-factor of 6 assumed for the workshop. Using Mack’s relation [Reference Mack58], this corresponds to a turbulence intensity of approximately
$0.24\% $
, which is the value used in the current work. Surface roughness has a significant effect on crossflow instabilities [Reference Dagenhart and Saric47, Reference Radeztsky, Reibert and Saric48]. To prevent roughness elements from destabilising the boundary layer, the surface of the CRM-NLF wind tunnel model was frequently sanded and polished, with the average size of the roughness elements after testing measured to vary from 0.83 to
$1.10\mu {\rm{in}}$
[Reference Lynde, Campbell, Rivers, Viken, Chan, Watkins and Goodliff11]. A value of
$1.00\mu {\rm{in}}$
is assumed in the current work.
Table 3. CRM-NLF wind tunnel test conditions*
Four grid levels are investigated using structured multi-block grids following the gridding guidelines provided by the workshop committeeFootnote
†
, with the grid characteristics listed in Table 4. A grid-refinement study is presented at the 2524 (
$ \approx {2.0^ \circ }$
) test conditions (see Table 3) using three modelling strategies: fully turbulent using the SA turbulence model [Reference Spalart and Allmaras37, Reference Allmaras and Johnson64] with the QCR2000 correction [Reference Spalart65], and free transition using the smooth transition model both with, SA-QCR2000-sLM2015cc, and without, SA-QCR2000-sLM2015, the Tollmien-Schlichting and stationary crossflow instability compressibility corrections. The original
${D_{{\rm{scf}}}}$
crossflow instability source term (Equation (A.3)) is used for the cases without the compressibility corrections, with the new
${F_{{\rm{onset,scf}}}}$
source term presented in Section 2.2 (Equations (15) and (16)) used for the cases with the compressibility corrections.
Table 4. CRM-NLF structured grid characteristics†
The grid-refinement study results at the 2524 (
$ \approx {2.0^ \circ }$
) test conditions obtained with free transition using both transition model variants are presented in Fig. 12. The difference in the drag coefficient on the finest grid level relative to the grid-converged value was calculated for the simulations with and without the compressibility corrections to be 0.65% and 0.66%, respectively. The iterative and grid convergence for this test case is also presented in more detail in previous work [Reference Piotrowski and Zingg52].
Figure 12. Grid-convergence results for the CRM-NLF simulations at the 2524 test conditions (
$\alpha \approx {2.0^ \circ }$
).
The residual convergence histories, illustrated by the norm of the total residual normalised by the value at the first iteration, for the free-transition simulations are presented in Fig. 13. The residuals are plotted against equivalent residual evaluations, which is the total wallclock time normalised by the average cost for computing a residual evaluation. It is important to note that the same solver parameters are used for both transition model variants on all four grids. The results demonstrate that the compressibility corrections and new
${F_{{\rm{onset,scf}}}}$
source term do not significantly affect the iterative and grid convergence of the model.
Figure 13. Grid-refinement study residual convergence histories for the CRM-NLF simulations at the 2524 test conditions (
$\alpha \approx {2.0^ \circ }$
).
The pressure and skin friction coefficient profiles for the three finest grids are extracted at nine spanwise stations across the wing and compared with the experimental profiles in Fig. 14. The upper-surface shocks help to differentiate the upper and lower surface skin friction coefficient profiles, as the shock causes a local reduction in the upper-surface skin friction coefficient due to the rapid flow deceleration. The compressibility corrections delay transition on the lower surface, which is dominated by crossflow instabilities, and on the upper surface in regions with less adverse and more favourable pressure gradients. The upper-surface transition locations appear sufficiently grid converged. However, at the
$\eta = 0.910$
station the lower-surface transition front produced by the transition model with the compressibility corrections moves downstream from the L1 to the L2 grid level.
Figure 14. Pressure and skin friction coefficient profiles for the CRM-NLF grid-refinement study at the 2524 test conditions (
$\alpha \approx {2.0^ \circ }$
) compared with the pressure profiles from the experiment [Reference Lynde, Campbell, Rivers, Viken, Chan, Watkins and Goodliff11] at varying spanwise stations
$\eta $
.
The pressure coefficient profiles produced by the fully turbulent simulations provide better agreement with the experiment at the inboard and midspan sections (
$\eta \lt 0.640$
) relative to the simulations with free transition. This counter-intuitive behaviour can be explained by the TSP images by Lynde et al. [Reference Lynde, Campbell, Rivers, Viken, Chan, Watkins and Goodliff11, Reference Lynde, Campbell and Viken12], which demonstrate a significant amount of bypass transition produced by surface imperfections, including the leading-edge pressure ports. Lynde et al. [Reference Lynde, Campbell, Rivers, Viken, Chan, Watkins and Goodliff11] attempted to estimate the natural transition front behind these bypass transition-induced turbulent wedges. As the simulated transition front moves downstream to better match this estimated natural transition front, a worse agreement is seen with the experimental pressure coefficient profiles. This behaviour was investigated in detail by Helm et al. [Reference Helm, Fehrs, Krimmelbein and Krumbein66]. At the further outboard stations, where there is less bypass transition observed in the experiment due to the reduced chord Reynolds number, the increased laminar extent of the boundary layer produced by the transition model with the compressibility corrections moves the shock locations aft to better agree with experimental data.
The free-transition skin friction coefficient profiles on the upper surface of the wing obtained on the three finest grid levels are overlaid with the estimated natural transition front from Lynde et al. [Reference Lynde, Campbell, Rivers, Viken, Chan, Watkins and Goodliff11] in Fig. 15. The compressibility corrections successfully push the transition front aft to better agree with the experimental results, particularly in the outboard regions. The results again demonstrate that the upper-surface transition front is sufficiently grid converged. The wing lower surface and fuselage nose boundary layers were tripped in the experiment. Therefore, transition fronts on the fuselage and lower surface of the wing are not available for comparison.
Figure 15. Upper-surface skin friction coefficient profiles for the CRM-NLF grid-refinement study at the 2524 test conditions (
$\alpha \approx {2.0^ \circ }$
) overlaid with the estimated natural transition front from the experiment [Reference Lynde, Campbell, Rivers, Viken, Chan, Watkins and Goodliff11].
Skin friction coefficient profiles produced by free-transition simulations at the 2523 (
$ \approx {1.5^ \circ }$
), 2525 (
$ \approx {2.5^ \circ }$
), and 2526 (
$ \approx {3.0^ \circ }$
) test conditions obtained on the L1 grid are illustrated in Fig. 16. Again, the compressibility corrections delay transition, producing transition fronts that better agree with experimental data. Regions of the transition fronts produced with the compressibility corrections remain upstream of the experiment. This is particularly evident for the 2523 (
$ \approx {1.5^ \circ }$
) test case, where transition is upstream of experimental data inboard near the Yehudi break, and at approximately mid-span at the 2525 (
$ \approx {2.5^ \circ }$
) and 2526 (
$ \approx {3.0^ \circ }$
) conditions. However, this behaviour is also present to a lesser extent in the work of Paredes et al. [Reference Paredes, Venkatachari, Choudhari, Li, Hildebrand and Chang67], Venkatachari et al. [Reference Venkatachari, Paredes, Choudhari, Li and Chang60], and Krimmelbein and Krumbein [Reference Krimmelbein and Krumbein68], who each investigate the CRM-NLF configuration using a variety of stability analysis methods. These differences could be due to uncertainties in the turbulence intensity, flow conditions, or aeroelastic deflection of the wind-tunnel model. In general, there is better agreement between the transition model with the compressibility corrections and the experimental data over the outboard portion of the wing.
Figure 16. Upper-surface skin friction coefficient profiles for the CRM-NLF obtained on the L1 grid overlaid with the estimated natural transition fronts from the experiment [Reference Lynde, Campbell, Rivers, Viken, Chan, Watkins and Goodliff11].
The residual convergence histories for the free-transition simulations are presented in Fig. 17. The same solver parameters are used for both transition model variants at all three angles of attack, which also match the solver settings used in the grid-refinement study at the 2524 test conditions shown in Fig. 13. The results again demonstrate similar iterative convergence between the two transition model variants.
Figure 17. Residual convergence histories for the CRM-NLF simulations on the L1 grid.
The pressure and skin friction coefficient profiles produced by the fully turbulent simulations and the free-transition simulations both with and without the compressibility corrections at the 2523 (
$ \approx {1.5^ \circ }$
), 2525 (
$ \approx {2.5^ \circ }$
), and 2526 (
$ \approx {3.0^ \circ }$
) test conditions are extracted at nine spanwise locations on the L1 grid and illustrated in Figs 18, 19, and 20, respectively. Similar to the 2524 results presented in Fig. 14, there is a consistent trend with the compressibility corrections delaying crossflow transition on the lower surface, and significantly delaying Tollmien-Schlichting transition on the upper surface in regions of favourable pressure gradient. There is again significant bypass transition observed in the experiment at the 2523, 2525, and 2526 test conditions, especially at the leading-edge pressure ports [Reference Lynde, Campbell, Rivers, Viken, Chan, Watkins and Goodliff11, Reference Lynde, Campbell and Viken12], which makes a detailed comparison of the pressure coefficient profiles difficult.
Figure 18. Pressure and skin friction coefficient profiles for the CRM-NLF at the 2523 test conditions (
$\alpha \approx {1.5^ \circ }$
) obtained on the L1 grid compared with the pressure profiles from the experiment [Reference Lynde, Campbell, Rivers, Viken, Chan, Watkins and Goodliff11] at varying spanwise stations
$\eta $
.
Figure 19. Pressure and skin friction coefficient profiles for the CRM-NLF at the 2525 test conditions (
$\alpha \approx {2.5^ \circ }$
) obtained on the L1 grid compared with the pressure profiles from the experiment [Reference Lynde, Campbell, Rivers, Viken, Chan, Watkins and Goodliff11] at varying spanwise stations
$\eta $
.
Figure 20. Pressure and skin friction coefficient profiles for the CRM-NLF at the 2526 test conditions (
$\alpha \approx {3.0^ \circ }$
) obtained on the L1 grid compared with the pressure profiles from the experiment [Reference Lynde, Campbell, Rivers, Viken, Chan, Watkins and Goodliff11] at varying spanwise stations
$\eta $
.
4.0 Conclusions
Compressibility corrections have been developed and applied to extend the SA-sLM2015 smooth local correlation-based transition model to transonic flow regimes. A compressibility correction for Tollmien-Schlichting instabilities has been developed to reproduce trends predicted by stability analysis and a crossflow source term function has been developed to apply a compressibility correction for stationary crossflow instabilities. These corrections and modifications do not impact the predictive capability of the model in the incompressible flow regime and do not have a significant impact on its iterative and grid convergence behaviour.
The smooth transition model both with and without the compressibility corrections is applied to a range of transonic external aerodynamic transition test cases. The results demonstrate that at higher Reynolds numbers the compressibility corrections successfully delay the predicted transition onset locations to better agree with experimental data. Specifically, the Tollmien-Schlichting compressibility correction delays transition in regions of favourable pressure gradient where stability analysis predicts a strong stabilising effect. The correction is less stabilising in adverse pressure gradients, enabling accurate prediction of transition locations at lower Reynolds number conditions. The crossflow compressibility correction prevents an inaccurate, upstream transition front from forming on the upper surface of the CRM-NLF wing where crossflow instabilities are not expected to be dominant.
The results demonstrate that the compressibility corrections successfully extend the empirical correlations in the smooth local correlation-based transition model to transonic flow regimes. Specifically, the corrections and the underlying transition model are investigated over a range of Mach numbers from 0.71 to 0.856, and for a range of Reynolds numbers from
$2 \times {10^6}$
to
$15 \times {10^6}$
. However, further work is required to validate RANS-based transition models over a wider range of flow conditions, specifically at the higher Reynolds numbers typical of commercial transport aircraft — the increased availability of high-quality data from experiments at these flight conditions is a necessary prerequisite for this development.
Acknowledgments
This work was partially funded by the NASA Transformational Tools and Technologies (TTT) and Advanced Air Transport Technology (AATT) projects, the Natural Sciences and Engineering Research Council (NSERC), and the University of Toronto. Computations were performed on the Niagara supercomputer at the SciNet HPC Consortium, a part of Compute Canada. The authors gratefully acknowledge the input and support provided by Dr. Thomas Reist, Dr. Holger Mai, and Dr. Cetin Kiris. In addition, the authors sincerely appreciate the generous support provided by Dr. Marco Costantini in providing high-quality data for the VA-2 aerofoil experiment.
Appendix A: SA-sLM2015 transition model equations
A.1 Transition onset momentum-thickness Reynolds number transport equation
\begin{equation}\frac{{\partial \tilde R{e_{\theta t}}}}{{\partial t}} + {u_j}\frac{{\partial \tilde R{e_{\theta t}}}}{{\partial {x_j}}} = {P_{\theta t}} + {D_{{\rm{scf}}}} + \frac{1}{{R{e_\infty }}}\frac{\partial }{{\partial {x_j}}}\left[{\sigma _{\theta t}}(\nu + {\nu _t})\frac{{\partial \tilde R{e_{\theta t}}}}{{\partial {x_j}}}\right]\end{equation}
\begin{equation}{P_{\theta t}} = \frac{{{c_{\theta t}}}}{{\rm{t}}}(R{e_{\theta t}} - \tilde R{e_{\theta t}})(1 - {F_{\theta t}})\end{equation}
\begin{equation}{D_{{\rm{scf}}}} = \frac{{{c_{\theta t}}}}{{\rm{t}}}{c_{{\rm{crossflow}}}}{\phi _{300}}(R{e_{\theta t,{\rm{scf}}}} - \tilde R{e_{\theta t}},0)({F_{\theta t}})\end{equation}
\begin{equation}{F_{\theta t}} = {F_{{\rm{wake}}}}{e^{ - {{\left(\frac{d}{\delta }\right)}^4}}},\quad {F_{{\rm{wake}}}} = {e^{ - R{e_S} \times {{10}^{ - 6}}}},\quad R{e_S} = \frac{{\rho {d^2}S}}{\mu }R{e_\infty }\end{equation}
\begin{equation}\delta = \frac{{50d\Omega }}{U}{\delta _{BL}},\quad {\delta _{BL}} = \frac{{15}}{2}{\theta _{BL}},\quad {\theta _{BL}} = \frac{{\tilde R{e_{\theta t}}\mu }}{{\rho U}}\frac{1}{{R{e_\infty }}}\end{equation}
\begin{equation}{c_{\theta t}} = 0.03,\quad {\sigma _{\theta t}} = 2.0,\quad {c_{{\rm{crossflow}}}} = 0.6\end{equation}
\begin{equation}{\rm{t}} = \frac{{500\mu }}{{\rho {U^2}}}\frac{1}{{R{e_\infty }}}\end{equation}
\begin{equation}R{e_\theta }_t = \frac{{\rho U\theta }}{\mu }R{e_\infty } = \left\{ {\begin{array}{*{20}{c}}{\left[1173.51 - 589.428T{u_\infty } + \dfrac{{0.2196}}{{Tu_\infty ^2}}\right]F({\lambda _\theta })} & {T{u_\infty } \le 1.3}\\[15pt]{331.50[T{u_\infty } - {{0.5658]}^{ - 0.671}}F({\lambda _\theta })} & {T{u_\infty } \gt 1.3}\end{array}} \right.\end{equation}
\begin{equation}F{({\lambda _\theta })_1} = 1 + 0.275\left[1 - {e^{[ - 35{\lambda _\theta }]}}\right]{e^{-\left[\frac{{T{u_\infty }}}{{0.5}}\right]}}\end{equation}
\begin{equation}F{({\lambda _\theta })_2} = {\phi _{300}}(F{({\lambda _\theta })_1},1)\end{equation}
\begin{equation}F{({\lambda _\theta })_3} = 1 - \left[ - 12.986{\lambda _\theta } - 123.66\lambda _\theta ^2 - 405.689\lambda _\theta ^3\right]{e^{ - {{\left[\frac{{T{u_\infty }}}{{1.5}}\right]}^{1.5}}}}\end{equation}
\begin{equation}F({\lambda _\theta }) = {\phi _{ - 300}}(F{({\lambda _\theta })_2},F{({\lambda _\theta })_3})\end{equation}
\begin{equation}{\lambda _\theta } = \frac{{\rho {\theta ^2}}}{\mu }\frac{{dU}}{{ds}}R{e_\infty },\quad U = ({u^2} + {v^2} + {w^2}{)^{\frac{1}{2}}}\end{equation}
\begin{equation}\frac{{dU}}{{ds}} = \left[\left(\frac{u}{U}\right)\frac{{dU}}{{dx}} + \left(\frac{v}{U}\right)\frac{{dU}}{{dy}} + \left(\frac{w}{U}\right)\frac{{dU}}{{dz}}\right]\end{equation}
\begin{equation}\frac{{dU}}{{dx}} = ({u^2} + {v^2} + {w^2}{)^{ - \frac{1}{2}}} \cdot \left[u\frac{{du}}{{dx}} + v\frac{{dv}}{{dx}} + w\frac{{dw}}{{dx}}\right]\end{equation}
\begin{equation}\frac{{dU}}{{dy}} = ({u^2} + {v^2} + {w^2}{)^{ - \frac{1}{2}}} \cdot \left[u\frac{{du}}{{dy}} + v\frac{{dv}}{{dy}} + w\frac{{dw}}{{dy}}\right]\end{equation}
\begin{equation}\frac{{dU}}{{dz}} = ({u^2} + {v^2} + {w^2}{)^{ - \frac{1}{2}}} \cdot \left[u\frac{{du}}{{dz}} + v\frac{{dv}}{{dz}} + w\frac{{dw}}{{dz}}\right]\end{equation}
\begin{equation}R{e_{\theta t,{\rm{scf}}}} = \frac{{\rho \!\left(\frac{U}{{0.82}}\right){\theta _t}}}{\mu }R{e_\infty } = - 35.088\ln \!\left(\frac{h}{{{\theta _t}}}\right) + 319.51 + f\!\left(\Delta H_{{\rm{crossflow}}}^ + \right) - f\!\left(\Delta H_{{\rm{crossflow}}}^ - \right)\end{equation}
\begin{equation}\Delta {H_{{\rm{crossflow}}}} = {H_{{\rm{crossflow}}}}(1.0 + {\phi _{ - 300}}({R_T},0.4)),\quad {R_T} = \frac{{{\mu _t}}}{\mu }\end{equation}
\begin{equation}\Delta H_{{\rm{crossflow}}}^ + = {\phi _{300}}(0.1066 - \Delta {H_{{\rm{crossflow}}}},0)\end{equation}
\begin{equation}f\!\left(\Delta H_{{\rm{crossflow}}}^ + \right) = 6200\!\left(\Delta H_{{\rm{crossflow}}}^ + \right) + 50000\!\left(\Delta H_{{\rm{crossflow}}}^ + \right)^2\end{equation}
\begin{equation}\Delta H_{{\rm{crossflow}}}^ - = {\phi _{300}}\!\left({-} \left(0.1066 - \Delta {H_{{\rm{crossflow}}}}\right),0\right)\end{equation}
\begin{equation}f\!\left(\Delta H_{{\rm{crossflow}}}^ - \right) = 75{\rm{tanh}}\left(\frac{{\Delta H_{{\rm{crossflow}}}^ - }}{{0.0125}}\right)\end{equation}
\begin{equation}\vec U = \left(\frac{u}{{\sqrt {{u^2} + {v^2} + {w^2}} }},\frac{v}{{\sqrt {{u^2} + {v^2} + {w^2}} }},\frac{w}{{\sqrt {{u^2} + {v^2} + {w^2}} }}\right)\end{equation}
\begin{equation}\vec \Omega = \left(\frac{{\partial w}}{{\partial y}} - \frac{{\partial v}}{{\partial z}},\frac{{\partial u}}{{\partial z}} - \frac{{\partial w}}{{\partial x}}\frac{{\partial v}}{{\partial x}} - \frac{{\partial u}}{{\partial y}}\right)\end{equation}
\begin{equation}\quad {\Omega _{{\rm{streamwise}}}} = |\vec U \cdot \vec \Omega |\end{equation}
\begin{equation}\quad \;\;{H_{{\rm{crossflow}}}} = \frac{{d{\Omega _{{\rm{streamwise}}}}}}{U}\end{equation}
A.2 Intermittency transport equation
\begin{equation}\frac{{\partial \gamma }}{{\partial t}} + {u_j}\frac{{\partial \gamma }}{{\partial {x_j}}} = {P_\gamma } - {E_\gamma } + \frac{1}{{R{e_\infty }}}\frac{\partial }{{\partial {x_j}}}\left[\left(\nu + \frac{{{\nu _t}}}{{{\sigma _f}}}\right)\frac{{\partial \gamma }}{{\partial {x_j}}}\right]\end{equation}
\begin{equation}{P_\gamma } = {c_{a1}}{F_{{\rm{length}}}}{F_{{\rm{onset}}}}\left[{\phi _{ - 300}}\left(\Omega ,\frac{{{U_\infty }}}{l}\frac{{\sqrt {{M_\infty }R{e_\infty }} }}{{20}}\right)\right]\sqrt \gamma (1 - {c_{e1}}\gamma )\end{equation}
\begin{equation}{E_\gamma } = {c_{a2}}{F_{{\rm{turb}}}}\left[{\phi _{ - 300}}\left(\Omega ,\frac{{{U_\infty }}}{l}\frac{{\sqrt {{M_\infty }R{e_\infty }} }}{{20}}\right)\right]\gamma ({c_{e2}}\gamma - 1)\end{equation}
\begin{equation}{F_{{\rm{onset}}}} = \frac{{\tanh \left(6\left({F_{{\rm{onset,1}}}} - 1.35\right)\right) + 1}}{2}\end{equation}
\begin{equation}{F_{{\rm{onset,1}}}} = \sqrt {{{\left(\frac{{R{e_S}}}{{2.60R{e_{\theta c}}}}\right)}^2} + {{({R_T})}^2}} \end{equation}
\begin{equation}R{e_{\theta c}} = 0.67\tilde R{e_{\theta t}} + 24\sin\! \left(\frac{{\tilde R{e_{\theta t}}}}{{240}} + 0.5\right) + 14\end{equation}
\begin{equation}{F_{{\rm{length}}}} = 44 - \frac{{44 - (0.50 - 3 \cdot {{10}^{ - 4}}(\tilde R{e_{\theta t}} - 596))}}{{{{(1 + {F_{{\rm{length,1}}}})}^{\frac{1}{6}}}}}\end{equation}
\begin{equation}{F_{{\rm{length,1}}}} = \exp \!\left({-} 3 \cdot {10^{ - 2}}(\tilde R{e_{\theta t}} - 460)\right)\end{equation}
\begin{equation}{F_{{\rm{turb}}}} = \left(1 - {F_{{\rm{onset}}}}\right)\exp \!\left({-} {R_T}\right)\end{equation}
\begin{equation}{c_{e1}} = 1.0,\quad {c_{a1}} = 2.0,\quad {\sigma _f} = 1.0,\quad {c_{e2}} = 50,\quad {c_{a2}} = 0.06\end{equation}
A.3 Spalart-Allmaras turbulence model (SA-neg-noft2)
\begin{equation}\frac{{\partial \tilde \nu }}{{\partial t}} + {u_j}\frac{{\partial \tilde \nu }}{{\partial {x_j}}} = {P_{\tilde \nu }} - {D_{\tilde \nu }} + \frac{1}{{\sigma R{e_\infty }}}\left[\frac{\partial }{{\partial {x_j}}}\left((\nu + \tilde \nu )\frac{{\partial \tilde \nu }}{{\partial {x_j}}}\right) + {c_{b2}}\frac{{\partial \tilde \nu }}{{\partial {x_i}}}\frac{{\partial \tilde \nu }}{{\partial {x_i}}}\right]\end{equation}
\begin{equation}{P_{\tilde \nu }} = \gamma \frac{{{c_{b1}}}}{{R{e_\infty }}}\tilde S\tilde \nu ,\quad {D_{\tilde \nu }} = \frac{{{c_{w1}}{f_w}}}{{R{e_\infty }}}{\left(\frac{{\tilde \nu }}{d}\right)^2}\end{equation}
Open access
