The Potential of Machine Learning Methods for Separated Turbulent Flow Simulations: Classical Versus Dynamic Methods

Abstract

Feasible and reliable predictions of separated turbulent flows are a requirement to successfully address the majority of aerospace and wind energy problems. Existing computational approaches such as large eddy simulation (LES) or Reynolds-averaged Navier–Stokes (RANS) methods have suffered for decades from well-known computational cost and reliability issues in this regard. One very popular approach to dealing with these questions is the use of machine learning (ML) methods to enable improved RANS predictions. An alternative is the use of minimal error simulation methods (continuous eddy simulation (CES), which may be seen as a dynamic ML method) in the framework of partially or fully resolving simulation methods. Characteristic features of the two approaches are presented here by considering a variety of complex separated flow simulations. The conclusion is that minimal error CES methods perform clearly better than ML-RANS methods. Most importantly and in contrast to ML-RANS methods, CES is demonstrated to be well applicable to cases not involved in the model development. The reason for such superior CES performance is identified here: it is the ability of CES to properly account for causal relationships induced by the structure of separated turbulent flows.

1. Introduction

The ability to accurately predict wall-bounded turbulent flows usually involving flow separation is highly essential to most aerospace and wind energy problems. This concerns, for example, aircraft wing aerodynamics, hypersonic flow, aircraft stability and control, wind turbine blade performance, wind farm layout and interaction, and offshore wind turbines. The challenge of simulating such flows is the appearance of detached flow behind complex geometries: flow separates from the surface of an object due to adverse pressure gradients or sudden changes in geometry and a recirculation region is formed (see, e.g., the flows described in Section 2 and Section 3). The latter usually has a significant influence on the overall flow behavior.

The advantages and disadvantages of basic computational methods usually applied to address such problems are well-known. The direct solution of Navier–Stokes equations, direct numerical simulation (DNS) [1], is computationally way too expensive for simulations of most flows of practical relevance. For example, comparisons with DNS are involved in the discussions presented below only in regard to periodic hill flows up to a Reynolds number ($Re$) of 5.6 K [2,3]. Stationary Reynolds-averaged Navier–Stokes (RANS) methods are applied with the motivation to fully model such flows. This means the turbulence equations applied in simulations are supposed to properly describe the structure of turbulent flows without the simulation of instationary turbulent motions. Technically, this can be accomplished by using a relatively large characteristic turbulence length scale (which implies a relatively large turbulent viscosity that damps turbulent fluctuations). Such simulations are computationally very efficient, but proved to be incapable of appropriately representing the structure of turbulent flows involving flow separation [1,4,5]. The large eddy simulation (LES), in particular wall-resolved LES (WRLES), is applied with the motivation to almost fully resolve turbulent flows. This means the equations applied are supposed to produce turbulent fluctuations (resolved flow) as seen in reality in separated turbulent flows. Technically, this can be accomplished by using a relatively small model length scale (which implies a relatively small turbulent viscosity). WRLES is able to accurately simulate separated turbulent flows, but unfortunately, the high computational cost of LES often does not allow simulations of practical relevance [1,6,7,8]. The search for appropriate ways to deal with these challenges has continued for several decades. Several strategies were applied to address these challenges [9,10].

A first basic approach, referred to as wall-modeled LES (WMLES), is to relax the resolution requirements of LES for wall-bounded turbulent flows by involving RANS components close to the wall. The basics of this approach can be found in Refs. [6,11,12,13,14,15,16,17], more recent reviews can be found in Refs. [18,19,20,21,22,23,24,25]. This approach is known to substantially reduce the computational cost of LES [23,24,26]. However, on top of simulation performance issues [27,28], this approach suffers from functionality issues. Such WMLES results depend on the use of different (equilibrium or non-equilibrium) wall models, definitions of regions where different models and grids are applied, different mesh distributions, and setup options to manage the information exchange between such different flow regions. Conceptual questions related to the use of alternative LES, Reynolds stress-constrained LES (RSC-LES) are discussed elsewhere [9,29,30,31]. The latter methods were applied to a large variety of cases like three-dimensional turbulent channel flow and flow past a circular cylinder [32], a full commercial aircraft simulation [33], periodic hill flows [34], compressible turbulent channel flows, and compressible flow past a circular cylinder [35,36], laminar-turbulent transition in a temporally developing channel flow [37], a U-duct flow [38], different airfoil type flows [39], and turbulent flow and heat transfer in a stationary ribbed duct [40].

A second basic approach, referred to as detached eddy simulation (DES), is to expand the capability of RANS equations by switching from the RANS turbulence length scale applied close to the wall to a much smaller LES-type length scale away from the wall. The basics of this approach can be found in Refs. [41,42,43,44,45,46,47,48], and more recent reviews can be found in Refs. [49,50,51,52]. The computational cost of this approach can very much exceed the cost of RANS simulations [53]. This approach also suffers from both simulation performance and functionality issues [27,28]. For each flow considered, there is the need to consider the validity of modifications of the basic approach. This involves, e.g., comparing the validity of zonal improved delayed DES (IDDES, which applies a synthetic turbulence generator at a certain flow-dependent RANS-LES interface) and non-zonal approaches, based, for example, on a shear layer-adapted (SLA) definition of the subgrid length scale (the filter size $\Delta$) [28]. It is worth noting that corresponding issues do not only apply to aerospace and wind energy problems but to a variety of other problems, as, for example, mesoscale and microscale modeling in regard to atmospheric simulations and many technical applications [54,55,56].

Several alternative strategies were studied to overcome these issues [9,10], e.g., unified RANS-LES (UNI-LES) [57,58,59,60,61,62,63], partially averaged Navier–Stokes (PANS) [64,65,66,67,68,69,70,71,72,73,74,75], partially integrated transport modeling (PITM) [76,77,78,79,80,81,82,83,84,85,86], and scale-adaptive simulation (SAS) methods [52,87,88,89,90,91,92]. These methods aim at theory-supported modifications of two-equation RANS models to enable the generation of resolved motion. A typical feature of these methods is that the resolution imposed by the model is often poorly realized in flow simulations. The latter can significantly affect the model functioning, as, e.g., the accurate prediction of characteristic separation zone features [27,28,93].

Given this situation, what are the options to enable feasible and reliable computational predictions of separated turbulent flows? Observations reported above reveal that the implementation of RANS or LES features in different equations cannot be expected to result in methods that perform equally well (using the same model) for a variety of separated turbulent flow types. The most promising alternatives are mathematically consistent extensions of RANS equations that either enable the inclusion of resolved motion, or pure RANS modeling that includes characteristic features of flow separation, meaning the significant flow structure modifications due to nonequilibrium processes.

The first alternative mentioned in the preceding paragraph was recently realized via the introduction of a new computational approach: continuous eddy simulation (CES). The latter minimal error simulation methods (which minimize the hybridization error) [9,27,28,93,94,95,96,97,98,99,100,101] were presented as an alternative to existing hybrid RANS-LES. The methodological setup can be applied to all popular turbulence models and in several types of hybridization. The most relevant feature of these methods is their derivation based on strict mathematics. The latter has significant implications for the model mechanism. In contrast to existing computational methods, in minimal error simulation methods, the model receives information about the actual flow resolution degree, which enables the model to appropriately respond by increasing or decreasing its contribution to the simulation. The latter is the essential requirement for a functional swing between almost modeled (RANS) and almost resolved (LES) regimes. It is also the requirement for reliable predictions of very high $Re$ flow regimes, which cannot be properly studied based on existing methods.

The second alternative is this use of popular machine learning (ML) methods for the extension of RANS methods (to enable RANS to properly deal with characteristic features of flow separation) [3,102,103,104,105,106,107,108,109,110,111,112,113,114,115,116,117,118,119,120,121]. The predominant approach is the modification of RANS equations considered (usually by introducing parameters that modify the production in turbulence equations). The latter parameters are determined first depending on physical space by using field inversion (FI) techniques in conjunction with high-fidelity reference data. In the second step, these spatial distributions of introduced model parameters are determined in dependence on model parameters available through the simulation by using ML techniques.

Given the two alternatives described in the preceding two paragraphs, the use of ML-RANS seems to be the more attractive approach because of their lower computational cost (ML-RANS do not attempt to resolve flow). However, such a preference should be accompanied by convincing evidence that ML-RANS is capable of providing reliable separated turbulent flow predictions even under conditions unseen in the model training. The latter question currently represents a question mark. In particular, there is a lack of clarity on what actually causes ML-RANS to be not convincing in this regard. Based on the latest developments, comparisons of ML-RANS with CES will be used to identify the reason for corresponding ML-RANS shortcomings. Simultaneously, these comparisons will be used to demonstrate the superior predictive power of CES in comparison to ML-RANS and other usually applied computational simulation methods. Characteristic features of these two alternatives considered will be presented in the following two sections including a discussion of differences. Conclusions will be presented finally.

2. CES Methods

2.1. Minimal Error Simulation Methods

CES methods [9,27,28,93,94,95,96,97,98,99,100,101] apply a sensibilization of RANS equations such that RANS equations become aware of the flow resolution using appropriate grids. More specifically, this approach enables the model to properly respond to the amount of flow resolution: the model contribution decreases (increases) if the degree of flow resolution is high (low). The latter is the essential requirement to properly transition from RANS regimes (dominated by the model contribution) and LES regimes (dominated by resolved flow). The CES approach can be seen as a dynamic ML approach: the model learns on the fly about the flow resolution amount and adjusts its contribution to the flow simulation properly. The mathematical approach, which is based on variational analysis, ensures that the model minimizes the hybridization error; see the discussion below related to

Table 1. Minimal error $k-\omega$ models: KOS and KOK hybridizations are considered in analysis options ${\mathcal{O}}_2$. Model errors $\lambda$, first variations, and resulting mode controls ${\beta }^{*}$ are provided. Variations applied are given in brackets.

$\frac{Dk}{Dt}=P-ϵ+D_k,\frac{D\omega }{Dt}=C_{{\omega }_1}{\omega }^2\left(\frac{P}{ϵ}-{\beta }^{*}\right)+D_{\omega }$,     $D_k=\frac{\partial }{\partial x_j}{\nu }_t\frac{\partial k}{\partial x_j}$,    $D_{\omega }=\frac{\partial }{\partial x_j}\frac{{\nu }_t}{{\sigma }_{\omega }}\frac{\partial \omega }{\partial x_j}$	(KOS)
${\displaystyle \frac{Dk}{Dt}=P-{\psi }_{\beta }ϵ+D_k,\frac{D\omega }{Dt}=C_{{\omega }_1}{\omega }^2\left(\frac{P}{ϵ}-\beta \right)+D_{\omega }}$,      ${\beta }^{*}=1+\beta -{\psi }_{\beta }$	(KOK)
• Analysis option ${\mathcal{O}}_2$ [without $Dk/Dt,D\omega /Dt$]:     [$\delta D_k/D_k=3\delta k/k-\delta ϵ/ϵ$,  $\delta D_{\omega }/D_{\omega }=\delta k/k$]
${\displaystyle {\lambda }_2=C_{{\omega }_1}{\omega }^2\left(\frac{P}{ϵ}-{\beta }_2^{*}\right)+D_{\omega }}$,    ${\displaystyle \delta \left(\frac{{\lambda }_2}{k}\right)=\frac{C_{{\omega }_1}}{L^2}({\beta }_2^{*}-1)\left[\frac{\delta L^2}{L^2}-\frac{\delta {\beta }_2^{*}}{{\beta }_2^{*}-1}\right]}$,    ${\displaystyle {\int }_{\beta }^{{\beta }_2^{*}}\frac{dx}{x-1}={\int }_{L_{tot}^2}^{L^2}\frac{dy}{y}}$,     ${\displaystyle \frac{{\beta }_2^{*}-1}{\beta -1}=L_{+}^2}$

. The approach can be applied in conjunction with all popular turbulence models, including Reynolds stress equation methods and probability density function (PDF) methods [100].

To illustrate the structure of CES models, we consider incompressible flow. The incompressible continuity equation $\partial {\tilde{U}}_i/\partial x_i=0$ and momentum equation are considered,

$${\displaystyle \frac{D{\tilde{U}}_i}{Dt}}=-{\displaystyle \frac{\partial (\tilde{p}/\rho +2k/3)}{\partial x_i}}+2{\displaystyle \frac{\partial (\nu +{\nu }_t){\tilde{S}}_{ik}}{\partial x_k}}.$$

(1)

Here, $D/Dt=\partial /\partial t+{\tilde{U}}_k\partial /\partial x_k$ denotes the filtered Lagrangian time derivative and the sum convention is used throughout this paper. ${\tilde{U}}_i$ refers to the i-th component of the spatially filtered velocity. We have here the filtered pressure $\tilde{p}$, $\rho$ is the constant mass density, k is the modeled energy, $\nu$ is the constant kinematic viscosity, and ${\tilde{S}}_{ij}=(\partial {\tilde{U}}_i/\partial x_j+\partial {\tilde{U}}_j/\partial x_i)/2$ is the rate-of-strain tensor. The modeled viscosity is given by ${\nu }_t=C_{\mu }{k}^{1/2}L$. Here, $C_{\mu }$ is a model parameter with standard value $C_{\mu }=0.09$, and L is a characteristic length scale. L can be calculated in different ways using $L=k^{3/2}/ϵ=k^{1/2}\tau =k^{1/2}/\omega$, where the dissipation rate $ϵ=k/\tau$ of modeled kinetic energy, the dissipation time scale $\tau$, and the turbulence frequency $\omega =1/\tau$ are involved.

Table 1 shows CES methods in regard to the use of the $k-\omega$ model. Here, $P={\nu }_t{S}^2$ is the production of k, where $S={\left(2{\tilde{S}}_{mn}{\tilde{S}}_{nm}\right)}^{1/2}$, $C_{{\omega }_1}=0.49$, and ${\sigma }_{\omega }=1.8$. In contrast to RANS, where ${\beta }^{*}$ is a constant, ${\beta }^{*}$ is considered here an undetermined parameter that needs to be chosen to minimize the hybridization error. The abbreviations KOS and KOK refer to the consideration of the $k-\omega$ model where the hybridization is accomplished in the scale equation (involving ${\beta }^{*}$ that needs to be determined) or k-equation (involving ${\psi }_{\beta }$ that needs to be determined), respectively. Mathematically, both variants (KOS and KOK) are equivalent if the coefficient relation ${\beta }^{*}=1+\beta -{\psi }_{\beta }$ is honored. These models can be analyzed by the two options considered (option $O_1$ and option $O_2$). Option $O_1$ is an exact hybridization where total viscosities need to be applied in turbulent transport terms. Option $O_2$ is a hybridization where the usual model viscosities are involved in turbulent transport terms and an approximation is applied (given by the neglect of substantial derivatives in regard to the model coefficient calculation). Only the second option is presented in Table 1, which was found to work very well in all applications. The way to obtain ${\beta }_2^{*}$ in option $O_2$ is explained in the lower part of Table 1. The hybridization error is given by ${\lambda }_2$. We consider the variation ($\delta$) of the normalized error ${\lambda }_2/k$ under consideration of variational properties of $D_k$ and $D_{\omega }$ given in the bracket terms. This leads to the second relation in the last row of Table 1. An extremal state is given by setting the right-hand side of this relation equal to zero. The integration of the resulting $\delta L^2/L^2=\delta {\beta }_2^{*}/({\beta }_2^{*}-1)$ from a RANS state (where ${\beta }_2^{*}=\beta$) to a state with a certain degree of flow resolution implies then ${\beta }_2^{*}=1+L_{+}^2(\beta -1)$. The same analysis in regard to the determination of ${\psi }_{\beta }$ provides the relationship ${\beta }_2^{*}=1+\beta -{\psi }_{\beta }$.

The difference between CES methods and other hybrid RANS-LES is the appearance of resolution indicators like $0\le L_{+}\le 1$ in model equations, which indicates the degree of flow resolution. The definition of $L_{+}=L/L_{tot}$ follows the definition of the modeled-to-total kinetic energy ratio $k_{+}=k/k_{tot}$: L is the modeled length scale contribution and $L_{tot}$ is the total length scale contribution. $L_{+}\approx 1$ indicates an almost completely modeled (RANS) regime, and $L_{+}\approx 0$ indicates an almost completely resolved (LES) regime. The mechanism of the model response to the actual flow resolution is the following. A higher resolution ($L_{+}$ becomes smaller) decreases ${\beta }^{*}$. Thus, there is less dissipation of $\omega$, and $\omega$ increases, which decreases the model viscosity ${\nu }_t=C_{\mu }k/\omega$.

CES methods show significant differences from WMLES and DES methods in regard to their functionality and computational cost. The issues related to properly setting up WMLES and DES methods were described in the second and third paragraphs of the introduction. In contrast to that, CES methods are independent of such functionality requirements, the model can be used as is, and CES can be expected to enable reliable predictions under conditions where validation data are unavailable. The differences in CES computational cost to other usually applied methods are discussed elsewhere [27,28].

2.2. Periodic Hill Flow Simulations

One of the applications of CES methods is the simulation of periodic hill flows, as illustrated in Figure 1 [93]. This flow is a channel flow involving periodic restrictions. This flow, which is used a lot for the evaluation of turbulence models [9], involves features such as separation, recirculation, and natural reattachment [122,123]. A thorough evaluation of the performance of CES methods in regard to simulating periodic hill flows at the highest $Re$ = 37 K (based on the hill height), for which experimental data for model evaluation are still available, can be found elsewhere, including evidence for the excellent CES prediction of mean velocities and stresses [93]. Figure 1 illustrates the remarkable differences to RANS predictions (see also the corresponding streamline results obtained by ML-RANS in Section 3.1. In contrast to RANS, the CES model provides a much more accurate characterization of the recirculation zone characteristics: the reattachment point is obtained as $x/h=3.78$, which differs from the corresponding value $x/h=3.76$ obtained in experiments only by 0.5%.

More specifically, the analyses presented in Ref. [93] involved four cases: $Re=37$ K ($G_{500}$), which is almost fully resolved, $Re=500$ K ($G_{120}$), which is almost fully modeled, and the transitional cases $Re=37$ K ($G_{120}$) and $Re=500$ K ($G_{500}$). Here, $G_{500}$ and $G_{120}$ refer to the use of grids with 500 K and 120 K grid points, respectively. Corresponding plots of velocity fluctuations are presented in Figure 2. These plots, too, underline the significant difference to steady RANS simulations, where fluctuations are excluded. The almost fully resolved $Re=37$ K ($G_{500}$) case reveals a very good resolution of fine scale structures. It is remarkable to see that grid and $Re$ variations imply smooth variations in these fluctuation structures, the merging of structures due to a coarser grid or higher $Re$. Interestingly, the almost fully modeled $Re=500$ K ($G_{120}$) case still involves strong fluctuations, which are as strong as under resolving conditions. The latter clearly speaks for a stable generation mechanism of fluctuations, which is implied by the derivation of CES methods.

2.3. NASA Hump Flow Simulations

The NASA wall-mounted hump flow illustrated in Figure 3 was studied by Seifert and Pack [124], aiming at the investigation of unsteady flow separation, reattachment, and flow control at a high Reynolds number $Re=c{\rho }_{ref}{U}_{ref}/\mu \approx 936 \mathrm{K}$ based on the chord length c and freestream velocity $U_{ref}$. Here, $\mu$ is the dynamic viscosity and the abbreviation $ref$ indicates the reference freestream conditions, which are determined at the axial point $x/c=-2.14$. The model reflects the upper surface of a 20-thick Glauert–Goldschmied airfoil that was originally designed for flow-control purposes in the early twentieth century. As a benchmark for comparison, we used the experiment conducted by Greenblatt et al. [125] without flow control. This benchmark case has been extensively documented on the NASA Langley Research Center’s Turbulence Modeling Resource webpage and has been widely used for evaluating different turbulence modeling techniques, as discussed in the 2004 CFD Validation Workshop.

Advantages of the CES-KOS model are illustrated in Figure 4. It may be seen that all methods involved in this comparison show a reasonable agreement with the pressure coefficient profile obtained via experiments. The WRLES predictions match the experimental measurement profile downstream and the model is capable of properly representing the dominant features of the flow. Regarding the reattachment region, we note that the second wall pressure peak is underpredicted by WRLES compared to CES-KOS and WMLES. Figure 4 also shows the mean skin friction coefficient obtained by CES-KOS, WMLES, and WRLES simulations, demonstrating their agreement with experimental values. In the separation zone, from $0\le x/c\le 0.65$, WRLES underpredicts the skin friction coefficient, while WMLES overestimates the actual peak. In regard to the post-reattachment, the $C_f$ profiles of WRLES and CES-KOS match relatively well, despite the very different computational setup of these methods.

A demonstration of the difficulty of dealing with accurate predictions for this flow in the RANS framework is given by the following. In an attempt to design more generally applicable RANS models, Srivastava et al. [105] recently considered two augmentation strategies, one of them (the reattachment augmentation) focusing on the improvement of predictions for the point of reattachment after smooth-body separation. The modeling approach uses the Spalart–Allmaras (SA) turbulence model [129] in conjunction with an adjustable parameter $\beta$ that modifies the production term. The motivation is illustrated in terms of Figure 5. The left-hand side shows the $\beta$ distribution obtained from an FIML method. Three regions can be seen: region 1 with $\beta >1$, region 2 with $\beta <1$, and region 3 with $\beta >1$. The focus is on region 3, leaving the modeling of regions 1 and 2 for future work. A very essential fact is that there is no use of ML in this approach, but an analytical model is provided for $\beta$. The latter includes a parameter ${\eta }_r$ which roughly correlates to region 3. The authors then apply a smoothed version of the expression $\beta =1$ if ${\eta }_r<c_{{\eta }_r}$ and $\beta =c_{\beta }^{max}$ otherwise. Here, the constants $c_{{\eta }_r}=0.9$ and $c_{\beta }^{max}=5$ are applied. A characteristic skin-friction result result is also provided in Figure 5 (SST + fix refers to the enhanced RANS model). It may be seen that this approach leads to improvements but significant model deficiencies remain; see the comparison with Figure 4.

2.4. Bachalo and Johnson Axisymmetric Transonic Bump Flow Simulations

The experimental configuration and the computational domain of the axisymmetric transonic bump considered by Bachalo and Johnson [130,131], along with the applied boundary conditions is illustrated in Figure 6. This case pertains to shock-triggered boundary layer separation induced by an axially symmetric bump mounted on a slim spherical cylinder, which extends 61cm upstream. The case reflects the upper surface of a transonic wing. It is characterized by a Mach number ($M_{\infty }$) of $0.875$ and a Reynolds number ($Re$) of $2.763$ M relative to the airfoil’s chord length c.

Streamwise velocity profiles obtained by CES-KOS, WMLES [132], and WRLES [133] are shown in the first row of Figure 7. The CES-KOS model predicts the streamwise velocity more accurately than WMLES and WRLES. In regard to turbulent shear stress profiles shown in the second row, we see a reasonable agreement of WMLES, WRLES, and CES-KOS with experimental data. Interestingly, CES-KOS and WRLES provide very similar results. In attached flow regions, WMLES overpredicts the turbulent shear stress. Due to its delayed reattachment point, WMLES predicts a faster separated shear layer growth and a higher maximum Reynolds stresses compared to CES-KOS. These observations are fully in line with past numerical studies using alternative models [134].

Figure 7 also shows pressure coefficient distributions obtained by CES-KOS, WMLES, and WRLES. Based on their sufficient flow resolution ability, CES-KOS and WRLES accurately predict pressure coefficient profiles. In contrast to that, WMLES fails to accurately capture the separation zone characteristics. Both CES-KOS and WRLES reasonably predict the shock location and post-shock pressure recovery. We note that WRLES results agree slightly better with the experimental data downstream of the bump ($1.1\le x/c\le 1.3$) compared to the CES-KOS model. Skin-friction coefficient distributions obtained by CES-KOS, WMLES, and WRLES are also shown in this figure. We observe that WMLES significantly underestimates the skin-friction coefficient in the separation region and fails to accurately represent the post-separation flow characteristics. The predictions of CES-KOS and WRLES are very similar, although CES-KOS agrees better with the experimental data in the $C_f$ plateau region upstream of separation. Hence, the most accurate flow predictions are provided by CES-KOS.

3. Ml Methods for Separated Turbulent Flows

Characteristic features of ML-RANS methods will be described next by focusing on hill-type simulations followed by an analysis of differences to CES methods. The usual methodological setup of such ML methods was described at the end of the introduction: with the exception of Ref. [106], a correction parameter $\beta$ modifies the production in turbulence equations. Then, the spatial distribution of this parameter is found via FI techniques followed by an ML representation of $\beta$.

3.1. 2D Periodic Hill Flow Simulations

By focusing on 2D periodic hill flow simulations, three ML-based extensions of RANS equations [3,106,107] are presented in this subsection. Essential characteristics of these studies are presented in

Table 2. ML-RANS simulations of hill flows: references and characteristics. Here, SA refers to the Spalart–Allmaras (SA) turbulence model [129]. The parameter $\alpha$ refers to different hill geometries.

Reference	Characteristics (Geometry, $\mathit{Re}$, Validation Data, Training and Validation, Output Variables)
Köhler et al. [107]	• Periodic hills, fixed geometry.
	• DNS/LES [2] [$Re=(0.7,1.4,2.8,5.6,10.6)$ K], Exp. [122] [$Re=37$ K].
	• Training: $Re=(0.7,1.4,2.8,10.6)$ K, validation: $Re=(5.6,37)$ K.
	• $k-\omega$ model, ML production [streamwise velocity profiles, streamlines, $C_f$, k, $\langle u^{\prime }{v}^{\prime }\rangle$].
Yan et al. [3]	• Periodic hills, variable geometry $\alpha =(0.5,0.8,1.0,1.2,1.5)$.
	• DNS [135] [$Re=5.6$ K, $\alpha =(0.5–1.5)$], LES [136] [$Re=10.6$ K, $\alpha =0.875$].
	• Training: $Re=5.6$ K, $\alpha =(0.5,0.8)$, validation: $Re=10.6$ K, $\alpha =0.875$.
	• SA, ML production [streamwise velocity profiles, velocity contours, $C_f$].
Volpiani et al. [106]	• Periodic hills, variable geometry $\alpha =(0.8,1.0,1.2)$.
	• DNS [135] [$Re=5.6$ K, $\alpha =(0.8–1.2)$], LES [136] [$Re=(2.8,10.6,19)$ K, $\alpha =1$].
	• Training: $\alpha =1$, $Re=(2.8,10.6,19)$ K (setup I), $Re=(2.8,5.6,10.6)$ K (setup II).
	• Validation: [$Re=19$ K, $\alpha =1$, setup II], [$Re=5.6$ K, $\alpha =(0.8,1.2)$, setup I and II].
	• SA, ML momentum eq. [streamwise and vertical velocity profiles, velocity contours, $C_f$].
Yan et al. [108]	• Fixed geometry: 3D symmetric bump and 3D FAITH hill.
	• Exp. of sym. bump [137,138] [$Re=130$ K], Exp. of FAITH hill [139] [$Re=500$ K].
	• Training: 3D sym. bump: full (ANN1), truncated (ANN2), sampled (ANN3) data.
	• Validation: 3D FAITH hill.
	• SA, ML production [velocity profiles, velocity contours, $C_p$].

.

One of these periodic hill flow simulations was presented by Köhler et al. [107]. The FIML model is referred to as an augmented $k-\omega$ model, the term inverse $k-\omega$ model refers to the model obtained by the FI method. The number of grid points was not reported. The $Re=5600$ case simulations show that there is basically no difference between the predictions of the augmented $k-\omega$ model and the inverse $k-\omega$ model. Streamwise velocity profiles show a clear improvement over RANS results, but the comparison with LES results reveals deficiencies. Corresponding observations can be made in regard to modeled kinetic energy (k) and Reynolds shear stress ($\langle u^{\prime }{v}^{\prime }\rangle$) predictions; the agreement with LES is even less satisfactory. The skin-friction coefficient ($C_f$) predictions show significant improvements compared to RANS results. The most relevant conclusions may be seen by considering the model application to the $Re=37$ K case. The $Re=37$ K case simulations show substantial differences between the inverse and augmented $k-\omega$ model: the inverse $k-\omega$ model is performing better.

• In particular, Figure 8 shows separation zone characteristics for the $Re=37$ K case unseen in the training. Although the augmented $k-\omega$ model is performing much better than the baseline $k-\omega$ model, there are clear disagreements with corresponding LES results: the separation zone characteristics are improperly described.

Figure 8. Periodic hill flow streamlines at $Re=37$ K obtained by Köhler et al. [107] for the baseline $k-\omega$ model (left) and the augmented $k-\omega$ model (right) with the reattachment point at $x/h=6.37$ and $x/h=5.02$, respectively (black dot), compared to the reattachment point from the experiment at $x/h=3.76$ (red dot). Reprinted by permission of the American Institute of Aeronautics and Astronautics, Inc.

Figure 8. Periodic hill flow streamlines at $Re=37$ K obtained by Köhler et al. [107] for the baseline $k-\omega$ model (left) and the augmented $k-\omega$ model (right) with the reattachment point at $x/h=6.37$ and $x/h=5.02$, respectively (black dot), compared to the reattachment point from the experiment at $x/h=3.76$ (red dot). Reprinted by permission of the American Institute of Aeronautics and Astronautics, Inc.

Another ML-based analysis of periodic hill flows was presented by Yan et al. [3], which covers geometry variations (see the different hill sizes depending on $\alpha$ in Figure 9). The ML-supported model was referred to as the FIML method. It involves the adjustable parameter $\beta$, which modifies the production term in the Spalart–Allmaras (SA) model [129]. The results depend on finding an appropriate regularization (a weight used in the objective function to limit variations of $\beta$ from the corresponding RANS value), which is important for guiding the FI to a physically realizable field [3]. One out of three options considered was found to be optimal. The grid applied used 81 points in the streamwise direction and 77 points in the normal direction. The $Re=5600$ case simulations show a relatively good model performance for the $\alpha =(0.8,1.2)$ cases involved in the training. There is an obvious improvement in RANS results, although the FIML model does not provide highly accurate predictions, in particular, close to the lower wall.

• The effect of geometry variations is shown in Figure 9. The FIML model performance is good for $\alpha =0.5$ (involved in model training), reasonable for $\alpha =1$ (slightly outside of the training range), and unsatisfactory for $\alpha =1.5$ (clearly outside of the model training). In particular, there is hardly any improvement compared to the performance of the RANS model.

The effect of an increased $Re$ is shown in Figure 10: the model is applied to the $Re=\mathrm{10,595}$ case in conjunction with $\alpha =0.875$ instead of the $Re=5600$ case considered so far. The FIML model obviously improves RANS predictions. Its performance is acceptable, although the FIML predictions are not highly accurate. The acceptable model performance may be related to the fact that the $\alpha =0.875$ case considered is not too far from the range considered in the model training. The authors state that the results show a certain extent of generalization ability both on different geometry configurations and Reynolds numbers.

Another comprehensive analysis of periodic hill flows was presented by Volpiani et al. [106]. Using a neural network (NN) and high-fidelity numerical data, a vectorial source correction term in the momentum equations (output quantity) applied in conjunction with the SA turbulence model is determined. The corresponding model is referred to as NN-RANS. Two scenarios are considered: setup I applies model training based on $Re$ = (2800, 10,595, 19,000) simulations, whereas setup II applies model training based on $Re$ = (2800, 5600, 10,595) simulations. The number of grid points was not reported.

• Using the same $\alpha =1$ geometry as in the training, the model used at a higher $Re$ =19,000 than that considered in training leads to significant improvements compared to RANS, although NN-RANS velocity discrepancies can be found in the upper channel region, where the extrapolation overestimated the volume forces.
• Similar observations are found in regard to applying the model for the $Re=5600$ case for different geometries than considered in the model training: there are significant improvements compared to RANS. It is of interest to note that setup I (trained using different geometries and Reynolds numbers) performs slightly better in the $\alpha =0.8$ case, but slightly underperforms setup II in the $\alpha =1.2$ case.
• We note that the range of modifications of unseen cases is smaller than in the approaches considered in the preceding paragraphs: the $Re=37$ K and $\alpha =(0.5,1.5)$ cases were not involved.

With respect to modifications of NN input parameter variations, the authors note that the choice of input features to improve NN-based RANS models is still a subject that needs further investigation.

3.2. 3D Hill-Type Flow Simulations

Next, let us consider the performance of 3D hill-type simulations based on ML extensions of RANS models, as presented by Yan et al. [108] and Ho and West [109]. Essential characteristics of these simulations are shown in Table 2 and

Table 3. Ho and West analysis [109]: RF refers to random forest, OTC refers to one-time correction, ITER refers to the iterative approach, 2D refers to training on 2D cases, 3D refers to training on all 2D and 3D cases. An NN approach also considered by Ho and West [109] is excluded here because of its poor performance. $\omega$ equation with ML production.

Model and Training	Validation
RF-2D-OTC
• 2D Channel [140] [DNS, $R{e}_{\tau }=2$ K]	• 2D periodic hills [141] [LES, $Re=10.6$ K]
• 2D NASA hump [125,142] [Exp. and LES, $Re=936$ K]	• 3D FAITH hill [139] [Exp., $Re=500$ K]
• 2D Curved backstep [143] [LES, $Re=13.7$ K]	• 3D non-axisymmetric bump [144]
• 2D Cylinder [145,146] [Exp. and LES, $Re=3.9$ K]	[Exp., $Re=16$ K]
RF-2D-ITER
• Same as above	• 2D periodic hills [141] [LES, $Re=10.6$ K]
RF-3D-OTC
• Same as above	• 3D non-axisymmetric bump [144]
• 3D FAITH hill [139] [Exp., $Re=500$ K]	[Exp., $Re=16$ K]

, respectively. The grid applied for the FAITH hill simulation by Yan et al. was $165\times 69\times 93=1.06$ M, whereas Ho and West [109] applied a grid of 0.4 M cells for this simulation (and 374,400 cells for the 3D non-axisymmetric bump simulation).

Yan et al. [108] report that the effect of the three ANN considered on the pressure coefficient distributions is very similar, which enables us to focus on the ANN1 results. Reference to the term artificial NN (ANN) is made here to provide a direct link to the cited work, bearing in mind that the terms ANN and NN have the same meaning. Characteristic results obtained by Yan et al. [108] for the unseen geometry of the 3D FAITH hill are shown in Figure 11 and Figure 12. Figure 11 shows substantial differences between the flow field structures seen in experiments and using the FIML model, although the ANN-enhanced SA model shortens the separation region compared to the original SA model.

• Figure 12 shows that there are only very minor differences between the prediction of the original SA model and the ANN-enhanced SA model for the unseen geometry of the 3D FAITH hill, both the original and enhanced SA models underestimate the length of the separation region.

Ho and West [109] applied the ML model called on individual cells to predict the correction factor during the simulation. Two options were considered: the correction factor calculation during every iteration in the iterative approach or just once in the one-time correction approach. The one-time correction (OTC) approach (having the lowest computational cost) was found to be most appropriate. The authors note that the iterative approach led to non-unique solutions based on the initialization of the flow field while the NN model produced large errors. In regard to 2D periodic hill flow simulations, the RF-2D-OTC model performed better than the RANS model as seen in mean velocities and k profiles, but the results are not very accurate in comparison to LES. The most interesting validations concern 3D FAITH hill and 3D non-axisymmetric bump applications.

• In regard to the 3D FAITH hill simulations, the results are very similar to the results of Yan et al. [108], shown in Figure 11 and Figure 12. A minor difference is that the velocity profiles obtained by Ho and West using their RF-2D-OTC model agree a little better with the experimental results.
• The 3D non-axisymmetric bump flow simulations and streamwise velocity contour plots are shown in Figure 13. The conclusion is that there are improvements compared to the RANS result, but significant differences to the experimental results remain. It may be seen that the inclusion of the 3D FAITH hill training data via the RF-3D-OTC model has a positive effect on the 3D non-axisymmetric bump results.

An interesting attempt to improve the ML approach was recently presented by Ho et al. [147]. This work is based on a probabilistic ML approach which reduces unwanted corrections, leading to more robust predictions. The latter is accomplished through only accepting predictions from the model above a certain confidence [147].

3.3. CES vs. ML-RANS Methods

The general idea of developing ML-RANS methods is highly attractive: such well-functioning methods can be computationally highly efficient because no attempt is made to resolve flows. The core of this idea is the development of equations that represent well the physics of the flows considered. More specifically, such equations are supposed to represent causal relationships in close agreement with first principles. An illustration of such causal relationships is given in Figure 14: this figure illustrates regions of high and low turbulence production (compared to the underlying standard RANS model) implied by different hill shapes. The task to improve RANS predictions means to explain how such areas of high and low production are related to the specifics of the hill considered.

In regard to the use of ML-RANS models, we observe the following:

• The intellectual breeding ground of ML methods (the data applied in the model development) are reflections of causal relationships, i.e., not the causal relationships themselves. For example, by taking reference to Figure 14, ML methods are capable of explaining how the high turbulence production relates to flow variables (like velocity gradients and velocity differences) and other model variables. Such reflections of causal relationships have to be expected to vary (significantly) depending on the geometry and $Re$ considered: see Figure 14 in regard to the geometry effect.
• As a consequence, ML methods can be expected to improve RANS results, but they cannot be expected to reliably predict flows involving geometries and $Re$ unseen in model testing. The latter expectation is fully confirmed by the results presented above. As described in Section 3.1 and Section 3.2, there may be hardly any improvements of RANS results.
• From a general viewpoint, the development of equations for separated turbulent flows in the steady RANS framework involving all relevant causal relationships seems to be out of reach. At the end, this would require, e.g., to explain how the geometry considered implies the production increase in Figure 14. There is no mathematical basis for doing this.

This situation is different in the frame of CES methods (an illustration of differences between ML-RANS methods and CES methods can be found in Figure 15):

4.

The geometry of separated turbulent flows implies the appearance of instantaneous detached flow (behind hills), which essentially determines the flow structure. In simulations, such flow corresponds to resolved flow, which is the instantaneous component of flow simulations. This leads to the task of explaining the causal relationship between the relative amount of resolved flow and the variation in production or dissipation in turbulence equations.

5.

The latter is not doable in the RANS framework, which excludes resolved flow, but it can be conducted in the CES framework. This is the core idea of the CES approach, which uses exact mathematics to determine the variation in the dissipation in turbulence equations as a consequence of resolved flow implied, e.g., by obstacles in the flow. In particular, the generality of such causality relationships is confirmed via the general mathematical derivation.

6.

From a general viewpoint, the use of CES methods may be seen as an experiment: causalities are involved to test whether the missing inclusion of causalities is indeed the origin of problems seen in regard to the use of ML-RANS methods for cases unseen in model testing. The results reported in Section 2 speak a clear language: the inclusion of causalities in CES methods enables much better flow predictions than given by ML-RANS models.

Figure 15. An illustration of conceptual differences between CES methods and ML-RANS methods in regard to hill-type separated turbulent flows. The illustration of required model modifications compared to standard RANS methods references Figure 14.

Figure 15. An illustration of conceptual differences between CES methods and ML-RANS methods in regard to hill-type separated turbulent flows. The illustration of required model modifications compared to standard RANS methods references Figure 14.

4. Summary

The issues of existing computational approaches related to the feasible and reliable simulation of high $Re$-separated turbulent flows led to the search for mathematically better-founded simulation methods. In a nonresolving framework, one way to address these questions is the very popular use of ML methods to enable improved RANS simulations. Another way is the use of minimal error CES simulation methods in the framework of methods being applicable under conditions ranging from resolving to almost modeled simulations.

The facts reported in Section 3, in regard to the performance of ML-RANS methods for the flows considered, speak a clear language. The setup of such ML methods usually involves a variety of adjustable model settings (see above) that need to be chosen to enhance the model performance. Such methods can provide improvements in RANS predictions (ML-RANS results are usually equivalent or better than RANS predictions), but in general, there is currently no indication that such methods are capable of providing reliable predictions for flow geometries or $Re$ not involved in model training (it is worth noting that the results reported here were presented in the very last years). One may argue that this situation can be improved by involving larger amounts of data in model training, hoping that the range of model applicability can be enhanced in this way. As discussed in regard to Figure 13, the latter can lead to improvements, which agrees with observations made by Fang et al. [148]. But there is no indication so far that this approach can lead to systematic improvements (a gradual increase in the predictive power of ML-RANS).

Table 4. ML-RANS results presented by Rumsey et al. [103]. Upper part: cases used in training. ZPG refers to zero pressure gradient, APG refers to adverse pressure gradient, $Ma$ refers to Mach number. Lower part: ML-RANS performance for cases unseen in training. All the lower part cases were trained by involving all upper part cases.

Case	$\mathit{Ma}$	$\mathit{Re}$	Unseen Case Predictions
ZPG-APG plate	0.1	$0.18\times {10}^6$
ZPG plate	0.2	$5.0\times {10}^6$
NACA 0012, $\alpha ={10}^{\circ }$	0.15	$6.0\times {10}^6$
NASA hump	0.1	$0.936\times {10}^6$
ZPG-APG plate, lower $Re$	0.1	$0.08\times {10}^6$	worse than RANS
Channel at very high $Re$	0.2	$80\times {10}^6$	worse than RANS
NACA 0012, other $\alpha$	0.15	$6.0\times {10}^6$	worse than RANS
Periodic hill	0.2	10,595	same as RANS
Axisymmetric jet	0.01	5601	much worse than RANS
Curved backstep	0.1	13,700	slightly better than RANS

shows facts of interest in this regard [103]. On average, this strategy simply fails. In line with Ref. [148] results, these results indicate that competing types of ML corrections are required for different cases.

The CES results reported in Section 2 also speak a clear language in regard to the three complex flow simulations considered so far. Such CES methods perform very well, better than WMLES or specifically designed DES methods, and at least as good or better than WRLES (at a little fraction of WRLES computational cost). These applications include evidence that these methods can work as resolving or partially resolving simulation methods. A remarkable feature of these CES methods is their well-balanced performance in contrast to other methods which usually provide appropriate predictions of some characteristics at the cost of other simulation characteristics. It is worth mentioning that such CES methods are independent of a variety of model setup options: they work on the basis of the RANS model considered, taking the mathematically determined hybridization into account.

The reason for the different features of ML-RANS and CES methods was identified in Section 3.3: it is given by the different way of accounting for causal relationships. The latter are a characteristic feature of separated turbulent flows, the flow geometry implies characteristic separation region features. The setup of ML-based methods enables the identification of correlations between model and flow variables, but the latter can vary a lot from flow to flow depending on the specific separation features. In contrast, CES methods can account for such causality via their inclusion of resolved motion implied by the type of flow separation. The generality of this inclusion of causality is confirmed via the general mathematical derivation and applications performed.

The conclusion of these facts is that the use of CES methods for separated turbulent flow simulations provides a much better alternative compared to the use of ML-RANS methods. Given that, what is the way to take maximum advantage of ML methods? First, in the frame of nonresolving methods, a more detailed investigation of the predictive power of ML methods in regard to relatively minor case variations would be of significant interest. This concerns, e.g., the prediction of $Re$ effects for the same geometry under conditions where different $Re$ cases are involved in the model training. The latter matters with respect to questions about asymptotic flow structures, and the known difficulty of resolving computational methods (DNS and LES) with extreme $Re$ predictions [149,150]. Second, in the frame of partially resolving methods, a merging of ML methods as considered here with CES concepts would be highly beneficial to take advantage of both approaches. The current development of CES methods is based on two-equation turbulence equations, which include an ($\omega$ or $ϵ$) scale equation, or an equation for the turbulent viscosity [96]. However, there are several usually applied codes that do not involve dissipation transport equations, as is currently the case with the widely used Weather Research and Forecasting Model (WRF) dealing with atmospheric flow simulations. By following the CES concept, ML can help to identify corresponding causal relationships between parameters or source terms in simplified turbulence equations (without scale equation, or even without an k transport equation). The latter can significantly contribute to the availability of corresponding computational methods via available codes for a variety of equation structures.