2.3.1. Tensor basis neural network

The TBNN, originally introduced by Ling et al. (Reference Ling, Kurzawski and Templeton2016), involves computing the invariants and the tensor basis assuming that the subgrid stress tensor is a function of the resolved rate of strain tensor, $\boldsymbol{\bar {S}_{\textit{ij}}}$ , and the resolved rotation rate tensor, $\boldsymbol{\bar {R}_{\textit{ij}}}$ . A detailed derivation of the tensor basis and invariants for the subgrid stress tensor can be found in works such as Stallcup et al. (Reference Stallcup, Kshitij and Dahm2022) and Wu & Lele (Reference Wu and Lele2025). Below, the invariants and tensor basis are listed for reference.

If one assumes that

(2.6) \begin{align} \tau _{\textit{ij}} = F(\boldsymbol{\bar {S}_{\textit{ij}}},\boldsymbol{\bar {R}_{\textit{ij}}}) , \\[-28pt] \nonumber \end{align}

(2.7) \begin{align} \boldsymbol{\bar {S}} \equiv \boldsymbol{\bar {S}}_{\textit{ij}} = \frac {1}{2}\left(\frac {\partial \bar {u}_i}{\partial x_{\kern-1pt j}}+\frac {\partial \bar {u}_{\kern-1pt j}}{\partial x_i}\right)\!, \\[-28pt] \nonumber \end{align}

(2.8) \begin{align} \boldsymbol{\bar {R}} \equiv \boldsymbol{\bar {R}}_{\textit{ij}} = \frac {1}{2}\left(\frac {\partial \bar {u}_i}{\partial x_{\kern-1pt j}}-\frac {\partial \bar {u}_{\kern-1pt j}}{\partial x_i}\right) \!, \\[0pt] \nonumber \end{align}

then the invariants can be written as

(2.9) \begin{align} \big\{tr \big(\boldsymbol{\bar {S}} \big), tr \big(\boldsymbol{\bar {S}^2} \big), tr \big(\boldsymbol{\bar {S}^3} \big), tr\big(\boldsymbol{\bar {R}^2}\big), tr\big(\boldsymbol{\bar {S}\bar {R}^2}\big), tr \big(\boldsymbol{\bar {S}^2\bar {R}^2} \big), tr \big(\boldsymbol{\bar {S}^2\bar {R}^2\bar {S}\bar {R}} \big) \big\} ,\end{align}

and the tensor basis can be written as

(2.10) \begin{align} \big\{\boldsymbol{I}, \boldsymbol{\bar {S}}, \boldsymbol{\bar {S}^2}, \boldsymbol{\bar {R}^2}, \boldsymbol{\bar {S}\bar {R} - \bar {R}\bar {S}}, \boldsymbol{\bar {R}\bar {S}\bar {R}}, \boldsymbol{\bar {S}^2\bar {R} - \bar {R}\bar {S}^2}, \boldsymbol{\bar {R}\bar {S}\bar {R}^2 - \bar {R}^2\bar {S}\bar {R}} \big\}. \end{align}

As seen, the invariants involve extremely high powers of the resolved strain rate and rotation rate tensors, which are not ideal due to aliasing errors, which increase with increasing power (Roberts & Mullis Reference Roberts and Mullis1987). Thus, instead of using the invariants, which are mathematically complete, one can simply use the velocity gradient tensor in place of the invariants, which is an input commonly used for neural network subgrid stress models (Park & Choi Reference Park and Choi2021; Kang et al. Reference Kang, Jeon and You2023; Cho et al. Reference Cho, Park and Choi2024). It is noted that there is no performance degradation using the velocity gradient tensor instead of the invariants as inputs to the neural network a priori. Even if the input to the neural network is changed, the backend of the neural network remains the same. The neural network will predict coefficients of the tensor basis expansion, listed here as $c_1-c_8$ (eight scalar values are predicted given the filtered velocity gradient tensor as an input)

(2.11) \begin{align} \tau _{\textit{ij}} & = c_1 \boldsymbol{I} + c_2 \boldsymbol{\bar {S}} + c_3 \boldsymbol{\bar {S}^2} + c_4 \boldsymbol{\bar {R}^2} + c_5 \big(\boldsymbol{\bar {S}\bar {R} - \bar {R}\bar {S}}\big) \nonumber \\[6pt] & \quad + c_6 \boldsymbol{\bar {R}\bar {S}\bar {R}} + c_7 \big(\boldsymbol{\bar {S}^2\bar {R} - \bar {R}\bar {S}^2}\big) + c_8 \big(\boldsymbol{\bar {R}\bar {S}\bar {R}^2 - \bar {R}^2\bar {S}\bar {R}}\big) . \end{align}

2.3.2. Overall architecture

Now that the backend TBNN is introduced, the overall neural network architecture can be seen in figure 2.

Figure 2.

The overall neural network architecture. The whole three-dimensional (3-D) domain is inputted into the neural network, and the output is the corresponding subgrid stress at each point. For each S4ND layer, a state space dimension of eight was chosen to keep the trainable parameters comparable to existing models in the literature. The average pooling operation will downsample all spatial dimensions by a factor of two. As seen, after the final S4ND layer, the neural network architecture splits into two ‘prediction heads’ to predict the ‘magnitude’ and ‘structure’ of the subgrid stress tensor, and both prediction heads use the tensor basis (each prediction head predicts its own set of coefficients to multiply the tensor basis components). See Wu (Reference Wu2025) for downloadable version of the model and an example input/output.

The inputs to the neural network are the 3-D, spatially varying velocity gradient tensor, as well as the strain rate and rotation rate tensor magnitudes. These inputs are passed into a U-Net S4ND architecture that allows for multiscale fusion of different spatial scales of turbulence. After the U-Net architecture, the neural network splits into two streams. Wu & Lele (Reference Wu and Lele2025) found that using two neural networks, one to predict the magnitude of $\tau _{\textit{ij}}$ while the other predicts the structure (normalized) $\tau _{\textit{ij}}$ , results in better performance. This neural network architecture retains the idea of predicting a ‘structure’ component and a ‘magnitude’ component, but is able to do it in one single neural network instead of splitting the task of predicting a ‘structure’ or ‘magnitude’ component to separate neural networks. As seen, the neural network splits into two paths after the U-Net architecture, with both paths using the TBNN (predicting eight coefficients that multiply the tensor basis). However, instead of using fully connected layers, which are unable to handle different domain sizes, 1 × 1 × 1 convolutions are used, which can emulate fully connected layers while still being able to generalize to various domain sizes.

Non-dimensionalizing the model and scaling the inputs and outputs are crucial to training a neural network model. Here, the inputs (velocity gradient tensor, strain rate and rotation rate tensor magnitudes, so a total of 11 quantities at each point in the domain) are normalized with a max–min scaling,

(2.12) \begin{align} x_{\textit{norm}} = \frac {x-x_{\textit{min}}}{x_{\textit{max}}-x_{\textit{min}}}2-1 \end{align}

where $x$ corresponds to every component of the velocity gradient tensor (or the strain rate/rotation rate tensor magnitudes). The maximum and minimum correspond to the 3-D domain’s maximum and minimum. This input normalization ensures that the neural network model inputs are always scaled between [−1, 1]. The tensor basis are also scaled, although with the Frobenius norm instead,

(2.13) \begin{align} W^{\textit{norm}}_{\textit{ij}} = \frac {W_{\textit{ij}}}{|W_{\textit{ij}}|} , \end{align}

where $W_{\textit{ij}}$ is a tensor corresponding to a tensor in the tensor basis expansion. Note that in (2.13), repeated indices do not imply summation. Furthermore, the output is also normalized accordingly to ensure that the model is non-dimensional. Following observations made by Kang et al. (Reference Kang, Jeon and You2023), the Clark’s gradient model term is taken as a normalization for $\tau _{\textit{ij}}$ , as seen in figure 2. However, one small modification is made. Work done by Scotti et al. (Reference Scotti, Meneveau and Lilly1993) introduced the anisotropic mesh length factor, $f$ , which multiplies the effective grid width $\varDelta _e$ as grid anisotropy may play a role in determining $\tau _{\textit{ij}}$ . Therefore, the normalization used for $\tau _{\textit{ij}}$ becomes

(2.14) \begin{align} \frac {\tau _{\textit{ij}}}{\left|\dfrac {1}{12}\varDelta _f^2\dfrac {\partial \bar {u}_{l}}{\partial x_k}\dfrac {\partial \bar {u}_{m}}{\partial x_k} \right|} \end{align}

where $\varDelta _f = f \varDelta _e$ , $f = 1 + ({2}/{27})[\ln(a_1)^2-\ln(a_1)\ln(a_2)+\ln(a_2)^2]$ and $\varDelta _e = (\varDelta _1 \varDelta _2 \varDelta _3)^{{1}/{3}}$ . Here $\varDelta _i$ denotes the grid spacing in the ith dimension, while $a_1$ and $a_2$ denote the grid spacing ratios (both less than 1).

2.4.1. Training data

Both filtered HIT and channel flow data are queried from the Johns Hopkins Turbulence Database (JHTDB), a service that provides high-quality DNS data for research purposes (Perlman et al. Reference Perlman, Burns, Li and Meneveau2007; Li et al. Reference Li, Perlman, Wan, Yang, Meneveau, Burns, Chen, Szalay and Eyink2008; Graham et al. Reference Graham2016). From Park & Choi (Reference Park and Choi2021), training on two different filter widths provides the neural network with more generalizability to different grid spacings. Thus, since this paper is focused on applications to aggressive, coarse grid spacing settings, after querying the DNS data, the DNS data is filtered to a filter width corresponding to $\varDelta _g = 16 \varDelta _{\textit{DNS}}$ and $\varDelta _g = 32 \varDelta _{\textit{DNS}}$ , where $\varDelta _{\textit{DNS}}$ denotes the DNS grid spacing. The filtering operation involves using a box filter in physical space. Following Wu & Lele (Reference Wu and Lele2025), two different domain sizes are also queried in order to help the neural network generalize to various domain sizes. The querying procedure used here is identical to that used in Wu & Lele (Reference Wu and Lele2025), and is also listed here for completeness.

The training data is split into a $80:10:10$ split, where 80 % of the total training data is used for training the neural network, 10 % of the data is used as a validation set during training for early stopping and 10 % of the data is used as the hold-out test set in the a priori testing section. To produce the training data for forced HIT, five snapshots are queried from the JHTDB database. To assemble the channel flow training data, four snapshots are queried. Each snapshot is queried at equal time intervals that are statically decorrelated from each other (time intervals different for HIT and channel flow). The DNS HIT resolution queried is $1024^3$ , with a Taylor Reynolds number ( $\textit{Re}_{\lambda }$ ) of 433 and a domain size of $[2\pi , 2\pi , 2\pi ]$ . The channel flow resolution queried is $2048 \times 512 \times 1536$ , with a friction Reynolds number ( $\textit{Re}_{\tau }$ ) of 1000 and a domain size of $[8\pi , 2, 3\pi ]$ when non-dimensionalized by the channel half-height. After the training data is queried, it is filtered. For example, filtering at a $16 \varDelta _{\textit{DNS}}$ filter width would result in the filtered HIT resolution being $64^3$ and the filtered channel flow resolution being $128\times 32\times 96$ . Then, the full domains are split into ‘subdomains’ of sizes $16^3$ and $32^3$ . This results in the total training data content in table 1.

Table 1.

Training data.

Table 1. Long description

The table presents the parameters for training data used for forced HIT and channel flow. It includes the number of realizations, filter width, domain size, and flow type. The table has 8 rows and 4 columns. The columns are labeled as Number of realizations, Filter width, Domain size, and Flow type. The rows are as follows: Row 1: Number of realizations, 320; Filter width, 16 delta DNS; Domain size, 16 cubed; Flow type, HIT. Row 2: Number of realizations, 320; Filter width, 32 delta DNS; Domain size, 8 cubed; Flow type, HIT. Row 3: Number of realizations, 40; Filter width, 16 delta DNS; Domain size, 32 cubed; Flow type, HIT. Row 4: Number of realizations, 40; Filter width, 32 delta DNS; Domain size, 16 cubed; Flow type, HIT. Row 5: Number of realizations, 320; Filter width, 16 delta DNS; Domain size, 16 cubed; Flow type, Channel. Row 6: Number of realizations, 320; Filter width, 32 delta DNS; Domain size, 8 cubed; Flow type, Channel. Row 7: Number of realizations, 40; Filter width, 16 delta DNS; Domain size, 32 cubed; Flow type, Channel. Row 8: Number of realizations, 40; Filter width, 32 delta DNS; Domain size, 16 cubed; Flow type, Channel.

As stated before, the training data is divided into 80 per cent for training, 10 per cent for validation and 10 per cent as a held-out test set. Each of the categories of training data is split according to the percentages above and then recombined into the total training, total validation, and total test set. In addition, filter width $\varDelta _g = 24 \varDelta _{\textit{DNS}}$ and filter width $\varDelta _g = 64 \varDelta _{\textit{DNS}}$ data are also queried, but are not included in the training set and are instead held out as part of the test set to test the neural network generalization capabilities.

2.4.2. Loss function and metrics

The loss function used for training follows from the loss function introduced in Wu & Lele (Reference Wu and Lele2025), involving a two-term loss function that enforces both the root mean squared error (RMSE) between the predicted and actual subgrid stress tensor, as well as the RMSE between the predicted and actual interscale energy transfer field, $\epsilon = - \tau _{\textit{ij}} \bar {S}_{\textit{ij}}$ . As such, the loss function involves

(2.15) \begin{align} L = \boldsymbol{{RMSE}}(\tau _{ij,\textit{pred}}-\tau _{\textit{ij}})/\boldsymbol{{RMS}}(\tau _{\textit{ij}})+ \lambda \boldsymbol{{RMSE}}(\epsilon _{\textit{pred}}-\epsilon )/\boldsymbol{{RMS}}(\epsilon ) , \end{align}

where $\lambda$ is a constant that can bias the neural network towards emphasizing the energy dissipation or the predicted values of $\tau _{\textit{ij}}$ . In the present work, following Wu & Lele (Reference Wu and Lele2025), $\lambda = 1$ is used, which indicates an equal bias between both loss function terms. The first term of this loss function is termed ‘RMSE’ (unbolded) and the second term of this loss function is termed ‘DRMSE’ (unbolded). Meanwhile, the bolded terms signify the operation of taking the RMSE or RMS of the particular quantity in the parenthesis. During training, the neural networks are trained with a batch size of eight, and the denominators of the loss function are computed per batch. As seen, this loss function directly enforces that the predicted and actual subgrid stress tensors are close, while simultaneously ensuring that the predicted and actual interscale energy transfer fields are also close. This ensures that the overall energy dissipation is accurate, meaning that the subgrid stress model is dissipating the correct energy from the resolved field.

For a posteriori tests, both qualitative and quantitative evaluations are important. However, quantitative metrics are often case dependent. Therefore, a few metrics are introduced that can be used to compare the various subgrid stress models with the filtered DNS solution that work for both forced HIT and channel flow. Note that these are necessary but not sufficient metrics to evaluating a subgrid stress model’s quality. The first metric involves computing a sum-absolute error (SAE) between the filtered DNS spectra and the LES spectra,

(2.16) \begin{align} \textit{SAE}_{E} = \int _k \big|\big(\log(E_{\textit{DNS}}(k))-\log(E_{\textit{LES}}(k))\big)\big| {\rm d}k \approx \sum _{k} \big|\log\!\big(E_{\textit{DNS}_k}\big)-\log\!\big(E_{\textit{LES}_k}\big)\big|, \end{align}

where the energy spectra is summed over a specific wavenumber band, with $k$ denoting the wavenumber. The logarithm ensures that differences in the spectra at large wavenumbers can still meaningfully increase this metric. Another metric involves computing the SAE between the filtered DNS and LES correlations (spatial or temporal),

(2.17) \begin{align} \textit{SAE}_{\textit{corr}} = \int _{\delta } \big|C_{\textit{DNS}}(\delta )-C_{\textit{LES}}(\delta )\big| {\rm d}\delta \approx \sum _{\delta } \big|C_{\textit{DNS}_{\delta }}-C_{\textit{LES}_{\delta }}\big|, \end{align}

where $C$ denotes the temporal or spatial correlation and $\delta$ represents the spacing (temporal or spatial) used to calculate the value of the correlation function. For channel flow, additional metrics involving the Reynolds stress profiles can also be developed. One involves computing the SAE over the various Reynolds stress profiles summed in the wall-normal direction,

(2.18) \begin{align} \textit{SAE}_{R_{xx}} = \int _y \big|R_{xx_{\textit{DNS}}}(y)-R_{xx_{\textit{LES}}}(y)\big| {\rm d}y \approx \sum _{y} \big|R_{xx_{\textit{DN}S_{y}}}-R_{xx_{\textit{LES}_{y}}}\big| ,\end{align}

where $R_{xx}$ represents the Reynolds stress over the ‘xx’ direction (i.e. 11, 22, 33), while $y$ denotes the wall-normal direction. Another quantitative metric for channel flow involves computing the SAE for the mean velocity profiles, summed over the wall-normal direction,

(2.19) \begin{align} \textit{SAE}_{u} = \int _{y} \big|\overline {u(y)}_{\textit{DNS}}-\overline {u(y)}_{\textit{LES}} \big| {\rm d}y \approx \sum _{y} \big|\overline {u(y)}_{\textit{DNS}}-\overline {u(y)}_{\textit{LES}} \big|, \end{align}

where $y$ once again is the wall-normal direction and $\overline {u(y)}$ is the LES or filtered DNS mean velocity profile (averaged over the homogenous directions, in this case x and z). The above four quantitative metrics do not represent a complete characterization of turbulence, but serve as ‘first-order’ metrics in order to compare various subgrid stress models quantitatively in an a posteriori setting.

4.1.1. Forced HIT at $16 \varDelta _{\textit{DNS}}$

The forced HIT simulation is run with a RK4–5 time stepping scheme (Bogacki & Shampine Reference Bogacki and Shampine1996) on a $64^3$ grid, which corresponds to a $16\varDelta _{\textit{DNS}}$ filter width. The spatial discretization scheme is Fourier collocation, possible since the boundary conditions are triply periodic, with a domain length of $2\pi$ in each direction. The forcing scheme forces the energy in the wavenumber band below 2 to be constant, and the non-dimensional dissipation, $\epsilon$ is set to 0.103 for all simulations. Due to Fourier collocation, $ {2}/{3}$ dealiasing is applied. For this simulation, all cases have a Taylor Reynolds number ( $\textit{Re}_{\lambda }$ ) of approximately 433. As seen in figure 8, all the LES simulations with various models have reasonable spectra. The eddy viscosity models, as expected, are more dissipative as compared with filtered DNS, while the neural network models better approximate the overall shape of the spectra but are slightly under-dissipative. From the ‘bump’ in the spectra, the TBNN subgrid stress model is the most under-dissipative, followed by the NN model, then the S4ND and NN-V model. The $ \textit{SAE}_{E}$ table validates this, as this metric suggests that the S4ND model and the NN-V model are the closest to the filtered DNS spectra, as lower $ \textit{SAE}_{E}$ values are closer to filtered DNS statistical quantities. Overall, the NN-V and S4ND models have spectra that are extremely similar, suggesting that S4ND, even with bandlimiting, is not significantly reducing the performance of the neural network in forced HIT.

Figure 8.

Forced HIT spectra for various SGS models a posteriori on a grid size corresponding to a filter width of 16 $\varDelta _{\textit{DNS}}$ .

Figure 9.

Forced HIT spatial, temporal, and space–time correlations for various SGS models a posteriori on a grid size corresponding to a filter width of 16 $\varDelta _{\textit{DNS}}$ .

Figure 9. Long description

The image contains multiple graphs depicting spatial, temporal, and spacetime correlations for various subgrid-scale models. The graphs include line plots and contour plots. Panel A shows temporal autocorrelation versus temporal spacing, with different models represented by distinct lines. The models include S4ND, NN, TBNN, Sigma, Static, and NN-V. Panel B displays contour plots of temporal spacing versus spatial spacing for different models, with contour lines indicating correlation values. Panel C presents another set of contour plots similar to Panel B but with a different range of spatial spacing. Panel D shows spatial autocorrelation versus spatial spacing, with the same models as in Panel A. The legends and labels are consistent across the panels, and the color scheme represents different correlation values. The graphs compare the performance of various models in capturing correlations at different grid spacings.

From figure 9, all neural network models are fairly good at predicting the temporal decorrelation, as they all follow the overall expected decay shape with the S4ND, NN-V and NN neural network models being slightly more accurate. When evaluating the spatial correlations, all models are fairly accurate for small spacings, but the S4ND and NN-V models are less accurate towards the larger spacings compared with both the sigma and NN models. While the overall spatial correlation shapes are preserved for all models, as seen from the $ \textit{SAE}_{\textit{corr}}$ for the spatial decorrelation, the S4ND and NN-V models have the largest error compared with the filtered DNS decorrelation. This suggests that the performance degradation for S4ND is not necessarily due to the state space model, but is likely a product of the normalization and inputs. From the space–time correlations, one can conclude that all the neural network models are far more accurate than either of the traditional models, especially for small spatial and temporal spacings. It is hard to determine which neural network model is more superior when it comes to space–time correlations. Overall, the benefits of neural network complexity is not clear in this set of simulations, as sometimes the S4ND model is better, while other times (such as when evaluating spatial correlations) it is worse. However, this simulation is run with an in-training set filter width of $16\varDelta _{\textit{DNS}}$ , where it has already been established that one does not expect much difference between the S4ND and the CNN-based models. This is reflected in this set of simulations as the S4ND model is sometimes better, and sometimes worse than the NN model (although it is quite consistent with the NN-V model due to the same input/output normalizations). The S4ND model’s strength, as seen a priori, is in extrapolative scenarios.

4.1.2. Channel flow at $16 \varDelta _{\textit{DNS}}$

Figure 10.

Channel flow streamwise spectra at $z/H = 0.46875$ for various SGS models compared with filtered DNS. Both the static Smagorinsky model and TBNN fail at predicting the overall spectra shape.

Similarly to the HIT cases, the channel flow simulation is also discretized with Fourier collocation in the homogeneous directions, while a sixth-order compact scheme is used in the wall-normal direction. The wall model used is not data-driven, and instead follows from the wall model used by Bou-Zeid et al. (Reference Bou-Zeid, Meneveau and Parlange2005). The simulation grid size is 128 × 96 × 32, while the domain size is [ $8\pi$ , $3\pi$ , 2], with the streamwise velocity pointing in the ‘first’ direction and the wall-parallel direction being the ‘second’ direction. All simulations run in this section correspond to a filter width of $16 \varDelta _{\textit{DNS}}$ with a friction Reynolds number of $\textit{Re}_{\tau } = 1000$ , with $\textit{Re}_{\tau } = {u_{\tau }H}/{\nu }$ , where $u_{\tau }$ is the friction velocity, $H$ is the channel half-height and $\nu$ is the kinematic viscosity. As seen in figure 10, in a more complex flow (as compared with HIT), the utility of more complicated neural network models is starting to show. The TBNN is unable to predict even the spectra properly, and is even worse than static Smagorinsky. Meanwhile the overall shapes associated with the NN and S4ND model is significantly closer to that of the filtered DNS profile. While the S4ND model is slightly over-dissipative at small scales, it follows the inertial range of the spectra better than the NN model when comparing these two profiles with the filtered DNS profile. However, as seen with NN-V, the S4ND model being over-dissipative is not due to the state space model, as the NN-V model also follows a similar spectra, suggesting that the differences in the spectra are more due to normalization and different inputs into the neural network. Overall, from the $ \textit{SAE}_{E}$ table, one can see that the S4ND model, NN-V and the sigma model are more similar to filtered DNS as compared with the other models.

Figure 11.

Channel flow macroscale quantities analysed for various SGS models. The mean velocity profiles and spatial correlations for the S4ND and NN model look the closest compared with filtered DNS. As seen, the TBNN completely fails in this setting, predicting jagged spatial correlations. (a) Mean streamwise velocity profiles. (b) Spatial correlations at $z/H = 0.46875$ .

Figure 11. Long description

The image contains two line graphs and two tables. Panel A shows a line graph of mean streamwise velocity profiles, with the x-axis labeled z/H and the y-axis labeled 〈U〉. The graph includes multiple data series represented by different colors and line styles, such as S4ND, NN, TBNN, FDNS, Static Smag., Sigma, and NN-V. Panel B shows a line graph of spatial correlations at z/H = 0.46875, with the x-axis labeled Spatial spacing and the y-axis labeled Spatial autocorrelation. This graph also includes multiple data series with different colors and line styles. Below the graphs, two tables present data for different models. The first table lists models and their corresponding SAE_u values, while the second table lists models and their corresponding SAE_corr values. The tables include models such as S4ND, NN, TBNN, Sigma, Static, and NN-V, with their respective values. The graphs and tables together analyze the performance of various SGS models in predicting mean velocity profiles and spatial correlations, comparing them with filtered Direct Numerical Simulation (DNS) data.

From figure 11(a) to 11(b), one can see a similar story when comparing the TBNN model with its more complex counterparts. The mean profiles associated with the TBNN model are inaccurate as compared with filtered DNS, not even predicting the shape properly, while the other three neural network models perform far better. In fact, the spatial correlations associated with the TBNN model at aggressive grid spacings highlights that it is completely inaccurate, as the spatial correlations become jagged. Meanwhile, the mean profiles predicted by all the other models that are not the TBNN are reasonable, with the NN model being more accurate as compared with the S4ND model. The NN-V and S4ND models have very similar mean profiles, which highlights that the S4ND architecture and the CNN architectures have relatively similar performance when the training set is still in filter-width. For the spatial correlations, both the S4ND and the NN model are far better than the traditional models, with the NN model being more accurate. When comparing NN and NN-V, one can isolate the effects of the normalization scheme and inputs. The improvements associated with the NN model is due to the inputs and normalization scheme used, not the architecture as NN-V has the same exact architecture as the NN model. From figure 11(b), it is seen that the NN-V model has worse spatial correlations as compared with the S4ND model, which highlights that even for an in-training set filter width, when comparing across the same normalizations, the convolutional models can start to lose performance. This is because channel flow operates on anisotropic grids (different ${\rm d}x, {\rm d}y, {\rm d}z$ spacings, but constant spacing in each direction), while convolutional-based models implicitly assume isotropic grids (HIT grids). This means that for anisotropic grids, the convolutional models will not perform as well when all things are considered equal (such as same normalization and same inputs).

Figure 12.

Channel flow Reynolds stress profiles for various SGS models at a filter width corresponding to 16 $\varDelta _{\textit{DNS}}$ . (a) Here $R_{11}$ varying as a function of wall normal distance. (b) Here $R_{13}$ varying as a function of wall normal distance. (c) Here $R_{33}$ varying as a function of wall normal distance.

Reynolds stress profiles are also plotted, which are seen in figure 12(a) to 12(c). Overall, both the S4ND model and NN model are more accurate compared with the other models, with the S4ND model and NN model trading off in accuracy. When comparing the NN model with the S4ND model, it would seem like for channel flow at this resolution, the NN model is slightly better than the S4ND model. This is seen in the $R_{13}$ profile, as the overall shape of the Reynolds stress profile is better captured by this model (although S4ND is better at capturing the $R_{11}$ profile). A different trend emerges when comparing the S4ND and NN-V model. In all cases, the NN-V model is either slightly worse ( $R_{11}$ ) or worse than the S4ND model ( $R_{13}, R_{33}$ ). This suggests that when accounting only for architectural effects, the S4ND model provides more accuracy on anisotropic grids than the NN-V model, which uses discrete convolutions. Since discrete convolutions assume isotropic grids, they will suffer a decrease in accuracy in flows with anisotropic grids such as channel flow when compared with the S4ND model.

While the NN-based model (with different inputs and normalization) is indeed slightly better than the S4ND model at this resolution, LES grid spacing is ultimately determined by the user, and therefore the models must also be evaluated on an out-of-training set filter width. As seen in Nguyen et al. (Reference Nguyen, Goel, Gu, Downs, Shah, Dao, Baccus and Ré2022), the S4ND model, through the use of bandlimiting, explicitly learns smoother kernels that would generalize better. However, this is an implicit form of regularization that will slightly decrease the performance of the S4ND model on filter widths that it was trained on. In the in-training set filter width a posteriori results for HIT and channel flow, the bandlimiting used actually does not decrease the performance of the S4ND model. For HIT, the S4ND model is on par with NN-V and for channel flow, since it accepts the different grid spacings in each dimension, the S4ND model is actually better than NN-V. This suggests that even though the S4ND model is biasing more towards extrapolation by regularizing the continuous convolution kernel to be smoother, the larger scales of the flow that are captured with the regularized convolution kernel are still enough to infer the subgrid stress tensor a posteriori.

4.2.1. Forced HIT at $64 \varDelta _{\textit{DNS}}$

Now, fixing all the previous simulation parameters to be identical to the $16 \varDelta _{\textit{DNS}}$ HIT simulation except for grid spacing, all models are run on a significantly coarser LES where the number of points is $16^3$ , corresponding to a filter width of $64 \varDelta _{\textit{DNS}}$ . The spectra is shown in figure 13.

Figure 13.

Spectra for various SGS models on a grid corresponding to an extrapolative filter width of $64 \varDelta _{\textit{DNS}}$ .

Figure 13. Long description

The image contains a table and a line graph. Panel A: The table lists different subgrid-scale (SGS) models and their corresponding values. The models listed are SAE_E, S4ND, NN, Sigma, Static, and NN-V, with values 0.67, 1.71, 2.71, 3.31, and 0.72 respectively. Panel B: The line graph displays the energy spectrum E(k) against the wavenumber k. The x-axis represents the wavenumber k ranging from 2 to 12, and the y-axis represents the energy spectrum E(k) on a logarithmic scale from 10^−4 to 10^0. The graph includes multiple lines representing different SGS models: S4ND (blue), NN (red), NN-V (dashed blue), Filtered DNS (black), Static Smag. (magenta), and Sigma (green). Each line shows the energy spectrum for the respective model, with trends indicating how each model performs in terms of energy dissipation and accuracy.

From the spectra, one can conclude that the sigma and static Smagorinsky models, consistent with the in-training set filter widths, are still over-dissipative. In addition, the S4ND model produces more accurate spectra as compared with the NN model. Furthermore, the S4ND and NN-V models produce similar spectra. Overall, the spectra suggest that the data-driven models are better at predicting the filtered DNS spectra, even when extrapolating, when compared with the traditional models. Furthermore, additional quantities such as correlations are shown in figure 14.

Figure 14.

Forced HIT spatial, temporal and space–time correlations for various SGS models on a grid corresponding to an extrapolative filter width of $64 \varDelta _{\textit{DNS}}$ .

Figure 14. Long description

The image contains multiple line graphs and contour plots depicting spatial, temporal, and spacetime correlations for various subgrid-scale models on a grid corresponding to an extrapolative filter width. The graphs are organized into four panels, each showing different aspects of the data. Panel A: This line graph shows temporal autocorrelation versus temporal spacing. The x-axis is labeled Temporal spacing, and the y-axis is labeled Temporal autocorrelation. The graph includes multiple lines representing different models: S4ND, NN, NN-V, Filtered DNS, Static Smag., and Sigma. Each line shows a different trend, with the Filtered DNS line generally showing higher autocorrelation values. Panel B: This contour plot shows temporal spacing versus spatial spacing with contour lines representing different levels of correlation. The x-axis is labeled Spatial spacing, and the y-axis is labeled Temporal spacing. The contour lines are color-coded, with a color bar on the right indicating the correlation values from 0 to 1. Different models are represented by different line styles and colors. Panel C: This contour plot is similar to Panel B but focuses on different data aspects. The x-axis is labeled Spatial spacing, and the y-axis is labeled Temporal spacing. The contour lines and color coding are similar to Panel B, representing different levels of correlation for various models. Panel D: This line graph shows spatial autocorrelation versus spatial spacing. The x-axis is labeled Spatial spacing, and the y-axis is labeled Spatial autocorrelation. The graph includes multiple lines representing different models: S4ND, NN, NN-V, Filtered DNS, Static Smag., and Sigma. Each line shows a different trend, with the Filtered DNS line generally showing higher autocorrelation values. The tables on the left and right sides of the image list the models and their corresponding SAE_corr values. The left table includes models S4ND, NN, Sigma, Static, and NN-V with their respective SAE_corr values. The right table includes models S4ND, NN, Sigma, Static, and NN-V with different SAE_corr values. The caption and citation provide context about the data sources and the methodology used for filtering and querying the DNS data.

At this more aggressive grid spacing, the temporal autocorrelations suggest that the S4ND model is far better at extrapolating than the NN-based model, and outperforms both traditional models. First, the traditional models are unable to predict the initial temporal decay properly, significantly overpredicting the temporal autocorrelation. Meanwhile, the NN model decorrelates far too quickly, and then decorrelates far too slowly as the temporal spacing increases. When comparing NN-V with S4ND, the temporal correlations for NN-V are surprisingly better than S4ND. This suggests that with this output normalization, some sense of scale-similarity is ‘baked-in’. In addition, since the HIT grid is isotropic, and the models are trained on two grid resolutions, the discrete convolutions in NN-V are able to learn a sense of scale similarity. This also highlights the importance of output normalizations and inputs, as the NN model, with different inputs and output normalization, was unable to extrapolate while the NN-V model was able to, despite both using CNNs. In general, the top three lowest $ \textit{SAE}_{\textit{corr}}$ include the S4ND, NN-V and static Smagorinsky models. However, both quantitative and qualitative analysis must be used, as $ \textit{SAE}_{\textit{corr}}$ for the static Smagorinsky model is lower than the S4ND model, but this can be attributed to the static Smagorinsky model crossing the filtered DNS line and thus the integrated distance between the two curves is smaller due to this crossover. As seen from the temporal decorrelation graph, the overall shape of the static Smagorinsky correlation is quite different from filtered DNS.

For the spatial correlations, when comparing the S4ND model with the NN model and traditional models, one can see that the S4ND model produces the most similar decorrelation shape to the filtered DNS profile. Meanwhile, the $ \textit{SAE}_{\textit{corr}}$ for the NN model is artificially lowered due to the line crossing the filtered DNS line. When comparing the S4ND model and the NN-V model, one can see that the S4ND model is more accurate than the NN-V model. For the spatial correlations, the S4ND model is the most accurate out of all the models. From the space–time correlations, one can conclude that the traditional models completely fail in predicting the correlations, as the black lines denote the predicted contour for filtered DNS, and the coloured lines denote the contours associated with a particular SGS model. For example, the blue contours denote the 0.9 contours for all the SGS models. Meanwhile, the traditional model ‘thin blue contours’ in the plots are extremely off, showing a significant overprediction of the correlation as compared with both the neural network models. While both the NN and S4ND model do well, the NN-V model performs a little worse, especially for small decorrelation scales. This is seen as it is slower to decorrelate than the S4ND model, as the blue contour is farther away from the 0.9 filtered DNS contour. Thus, one can see that the S4ND model, with the continuous kernel and bandlimiting, overall produces more accurate space–time correlations and spatial decorrelations as compared with all the other SGS models evaluated. Meanwhile, while the NN-V model produces good temporal decorrelations, it still suffers in accuracy when comparing the small-scale space–time decorrelations and the overall spatial autocorrelations, suggesting that S4ND’s architectural effects are directly responsible for this increase in accuracy. By comparing NN and NN-V, one can also conclude that correct inputs and output normalization is essential for scaling to coarser grid spacings.

4.2.2. Channel flow at $64 \varDelta _{\textit{DNS}}$

Channel flow is also evaluated at a more aggressive, extrapolative grid spacing. Keeping all other simulation parameters fixed, the grid spacing is changed to $32 \times 24 \times 8$ , corresponding to a filter width of $64 \varDelta _{\textit{DNS}}$ . Otherwise, the simulation parameters are kept the same as the finer grid simulation, corresponding to a $\textit{Re}_{\tau }$ of 1000. From the mean profiles, shown in figure 15, the S4ND model is able to preserve the slope of the mean streamwise velocity profile as compared with filtered DNS, while the NN model fails. Meanwhile, the two traditional models are more inaccurate as they are farther from the filtered DNS profile. When comparing the S4ND model with the NN-V model, one can see that while the NN-V model is able to improve upon the NN model form by using different inputs and output normalization, the overall slope is still incorrect when compared with the S4ND model. This is also seen in the $ \textit{SAE}_u$ tables, with the S4ND being the lowest (and the NN model being not far behind only because it crosses the filtered DNS profile).

Figure 15.

Channel flow mean profiles for various SGS models on a grid corresponding to an extrapolative filter width of $64 \varDelta _{\textit{DNS}}$ .

The Reynolds stress profiles are also plotted, as seen in figure 16. From these profiles, the accuracy of the S4ND model in an extrapolative setting is highlighted. Comparing the S4ND model with the NN model, the S4ND model is far more accurate as the NN model has significant trouble extrapolating. The NN model profile shapes are extremely off and the magnitudes are incorrect even when compared with the traditional models. As seen from every single profile, the S4ND is far better at preserving the overall shape and magnitude as compared with both the NN model and all the traditional models. However, this could also potentially be because the NN model’s output scaling, which was different from the S4ND’s output scaling, was just not as optimized for extrapolation to coarse grids. This hypothesis can be ruled out by looking at the NN-V model, which has the same exact inputs and output scaling as the S4ND model. Here, the NN-V model also is significancy worse than the S4ND model. While it does not overpredict as much as the NN model, suggesting that the output scaling used for the NN-V model and S4ND model are better at scaling, the overall profile shapes are oftentimes incorrect, especially $R_{33}$ , where it overpredicts the Reynolds stress peak such that it is over five times that of filtered DNS. This suggests that it is not the output scaling that is causing the performance improvement of S4ND for channel flow extrapolation, but solely due to the S4ND architecture. This highlights the extrapolatory capability of the S4ND model as it outperforms the NN and NN-V model while also being significantly better than the traditional models. In all situations for coarse-grid channel flow, the S4ND-based model, when integrated into a LES simulation, gives results closest to the filtered DNS profiles.

Overall, combined with the HIT results from the previous section, the S4ND model is a more flexible model since it is more accurate in extrapolating than any of the other neural network models tested (especially for channel flow), and is more accurate than the traditional SGS models even when extrapolating in grid resolution. Furthermore, by comparing with NN-V, one can highlight that while the output scaling does play a role in helping to improve extrapolation, the S4ND architecture is still key in extrapolating, especially in channel flow. This also suggests that the in-training set filter width results, where oftentimes the NN model was better, highlight that the normalization used for the NN model helped contribute to the NN model ‘overfitting’ to the in-training set filter width cases, producing better results a posteriori on these cases and then completely failing to extrapolate to HIT or channel flow on coarser grids. This can be contrasted with NN-V and S4ND, where the inputs and output normalization allowed for less a posteriori performance for the in-training set filter width cases, but then allowed for better extrapolation to coarser grids. In the simpler HIT case with isotropic grids, the normalization allowed the NN-V model, which uses CNNs, to extrapolate well to HIT, but it still failed in the more complex channel flow case with anisotropic grids. The S4ND meanwhile, was able to extrapolate to both coarser grid cases reasonably well, as a consequence of its architecture. This validates what was brought up in the previous section, where the use of bandlimiting is effectively a regularizer that smooths out the convolution kernel, allowing the S4ND model to be more accurate in out of training set filter widths, while the models with CNNs are more likely to ‘overfit’ to the specific grid resolutions that it has been trained on, particularly for channel flow, where due to anisotropic grids, the primary assumption of CNNs is already violated (isotropic grids). Physically, the S4ND coarse grid extrapolation cases highlight that the continuous convolution kernels with bandlimiting allow the model to respect the scale-dependent structure of turbulence across filter widths. When the filter width increases (grid gets coarser), more energy is filtered out, and the subgrid stress model must account for a larger fraction of the unresolved energy. Particularly for channel flow, since the S4ND model produces more accurate Reynolds stress profiles and mean profiles on coarser grids when compared with the NN-V (which has identical inputs, outputs and normalization), this suggests that the continuous formulation adjusts to the balance between resolved and modelled scales more naturally than discrete convolutions. The improved performance of the S4ND model on coarser grid resolutions is therefore attributed to the continuous formulation with bandlimiting, which ensures that the learned closure adapts to changing filter width, maintaining filter-consistent spectral content while capturing the interscale energy transfer (which can be seen directly in the a priori results).

Figure 16.

Channel flow Reynolds stress profiles for various SGS models on a grid corresponding to an extrapolative filter width of $64 \varDelta _{\textit{DNS}}$ . (a) Here $R_{11}$ varying as a function of wall normal distance. (b) Here $R_{13}$ varying as a function of wall normal distance. (c) Here $R_{33}$ varying as a function of wall normal distance.

5.1.1. The HIT

The Reynolds number scaling results for forced HIT are shown in figures 17(a) and 17(b). In these figures, the S4ND model is compared with both the static Smagorinsky and sigma models, as well as filtered DNS to evaluate the extrapolatory capability of the S4ND model.

Figure 17.

Macroscale quantities analysed for forced HIT at a higher Reynolds number not seen in the training set. (a) Forced HIT spectra. (b) Forced HIT spatial correlations.

Figure 17. Long description

The image contains two line graphs and two tables. Panel A shows a line graph titled Forced HIT spectra. The x-axis is labeled k and ranges from 10^0 to 10^1. The y-axis is labeled E(k) and ranges from 10^−4 to 10^−1. The graph includes four lines representing different models: S4ND, FDNS, Sigma, and Static Smag. Panel B shows a line graph titled Forced HIT spatial correlations. The x-axis is labeled Spatial spacing and ranges from 0 to 0.30. The y-axis is labeled Spatial autocorrelation and ranges from 0.75 to 1.00. This graph also includes four lines representing the same models: S4ND, FDNS, Sigma, and Static Smag. Below the graphs, there are two tables. The first table lists models and their corresponding SAE_E values: S4ND with 4.37, Sigma with 6.22, and Static with 12.8. The second table lists models and their corresponding SAE_corr values: S4ND with 0.12, Sigma with 0.16, and Static with 0.20.

Figure 18.

Compensated spectra for forced HIT at a higher Reynolds number not seen in the training set.

As seen, when the Reynolds number is larger than the training set Reynolds number, the S4ND model still produces reasonable results for both spectra and spatial correlations. Overall, the same trends hold when comparing the S4ND simulations that were run with an in-training set Reynolds number versus an extrapolatory Reynolds number. For example, the energy spectra for S4ND is still slightly underdissipative, while the spatial correlations are more accurate for the S4ND model compared with the spatial correlations produced by the traditional subgrid stress models. In this case, the spatial correlations are truncated since there is only one ensemble at this Reynolds number for the filtered DNS profile (JHTDB provides one snapshot), and thus the spatial correlation is truncated at a point where it is still reasonably converged. When taking a look at the SAE for both the energy spectra and the correlations, the S4ND model has the lowest error out of all the models compared. The compensated spectra, which premultiplies the energy spectra with a $k^{ {5}/{3}}$ factor, are shown in figure 18. As seen, the filtered DNS reference has an inertial range, but the inertial range varies from being slightly less than the $k^{{5}/{3}}$ scaling to slightly more than the $k^{{5}/{3}}$ scaling. This is also seen in the SGS models, with the S4ND model being consistent with previous results (slightly under-dissipative) while the traditional SGS models are also consistent with previous results (slight over-dissipative). Overall, the S4ND model is able to successfully extrapolate to a higher Reynolds number for forced HIT.

5.1.2. Channel flow

Figure 19.

Macroscale quantities analysed for channel flow at a higher Reynolds number not seen in the training set. (a) Streamwise spectra at z/H = 0.46875. (b) Spatial correlations at $z/H = 0.46875$ .

Figure 19. Long description

The image contains two line graphs comparing streamwise spectra and spatial correlations for different models. Panel A shows streamwise spectra at z/H = 0.46875. The x-axis represents the wavenumber (k) on a logarithmic scale, and the y-axis represents the energy spectrum (E(k)) on a logarithmic scale. The graph includes four lines representing different models: S4ND, Static Smag., Sigma, and FDNS. The legend indicates the corresponding colors for each model. The S4ND model has a value of 2.32, Sigma has 1.41, and Static has 5.79. Panel B shows spatial correlations at z/H = 0.46875. The x-axis represents spatial spacing, and the y-axis represents spatial autocorrelation. The graph includes four lines representing different models: S4ND, Filtered DNS, Sigma, and Static Smag. The legend indicates the corresponding colors for each model. The S4ND model has a value of 1.59, Sigma has 2.16, and Static has 5.79. The graphs illustrate the performance of different models in capturing streamwise spectra and spatial correlations.

For channel flow Reynolds number scaling, first the spectra at a specific wall normal location and the spatial correlations are shown in figures 19(a) and 19(b). As seen in the plots, the S4ND model produces spectra and spatial correlations that match filtered DNS profiles when the Reynolds number is out of the training set. There is no obvious performance degradation when compared with the in-training set filter width, and the overall profiles match the filtered DNS profiles. Furthermore, the S4ND Reynolds stress profiles are seen in figures 20(a), 20(b) and 20(c). Overall, the S4ND model is able to capture the overall shape and magnitude of the Reynolds stress profiles associated with channel flow. From all these plots, the S4ND model is in general, consistent with the in-training set Reynolds number, with no obvious performance degradation when the Reynolds number is out of the training set. This is corroborated by looking at the SAE, where the S4ND model produces reasonable SAE for all quantities. Combining both the channel flow and forced HIT results, one can conclude that the S4ND model is able to extrapolate to higher Reynolds number flow cases and still match filtered DNS reasonably well.

Figure 20.

Channel flow Reynolds stress profiles for various SGS models at a higher Reynolds number not seen in the training set. (a) Here $R_{11}$ varying as a function of wall normal distance. (b) Here $R_{13}$ varying as a function of wall normal distance. (c) Here $R_{33}$ varying as a function of wall normal distance.

5.2.1. HIT

The HIT simulation is evaluated on a $16^3$ grid, corresponding to a grid spacing in each direction that is twice as large as the training set. Furthermore, five Reynolds numbers based on the Taylor microscale are tested: $\textit{Re}_{1} = 4\,800$ , $\textit{Re}_{2} = 19\,500$ , $\textit{Re}_{3} = 56\,900$ , $\textit{Re}_{4} = 515\,000$ , $\textit{Re}_{5} = 508\,900\,000$ . The S4ND model is only trained on $\textit{Re} = 433$ . As seen in figures 21(a) and 21(b), the neural network is able to ensure that the LES simulation is stable, even when extrapolating in both grid spacing and $\textit{Re}$ . Overall, the spectra and spatial correlations look reasonable for all $\textit{Re}$ , showing that the S4ND model produces reasonable results even if the Reynolds number is over 1 000 000 times that of the training set. Given the coarseness of the grid, the results should not change drastically with Reynolds number, which is demonstrated in these plots. In addition, up to $\textit{Re}_3$ , the trends are monotonic, suggesting that the S4ND model can capture trends up to around 70 times the Reynolds number in the training set. From $\textit{Re}_1$ to $\textit{Re}_5$ , at no point is a ‘tail-up’ in the spectra seen, which is expected for ANN-based subgrid stress models that are trained on only one Reynolds number, as seen in Cho et al. (Reference Cho, Park and Choi2024).

Figure 21.

Extreme Reynolds number scaling for HIT. (a) Energy spectra at various Reynolds numbers. (b) Spatial correlations at various Reynolds numbers.

Figure 22.

Various mean profiles for channel flow extrapolation to extreme Reynolds numbers. (a) Channel flow mean profiles for various Reynolds numbers. (b) Channel flow mean profiles for various Reynolds numbers with sigma model and fits.

5.2.2. Channel flow

The same neural network is also integrated into a channel flow simulation with varying Reynolds numbers to test the extrapolatory capabilities of the neural network, on a grid of size $32 \times 42 \times 8$ on a domain of size $8 \pi \times 3 \pi \times 2$ (grid spacing is coarser than the training set). The Reynolds number based on friction velocity tested involve $\textit{Re}_1 = 5200$ , $\textit{Re}_2 = 20\,000$ , $\textit{Re}_3 = 50\,000$ , $\textit{Re}_4 = 500\,000$ and $\textit{Re}_5 = 500\,000\,000$ . The neural network is trained only on Reynolds number based on friction velocity of $\textit{Re} = 1000$ .

From figure 22(a), when scaled with friction velocity, the channel mean streamwise velocity profiles collapse nicely, suggesting that these profiles are reasonable for all Reynolds numbers. To fully verify this, the log-law profiles are plotted, along with the sigma model mean profile results in figure 22(b). Due to the large offset from the S4ND model, the sigma model was also run as verification to isolate if this offset was associated with the subgrid stress model or the wall model. As seen, using the classic offset of $B=5.0$ is insufficient, likely due to the error in the wall model. This error can be attributed to the wall model since two different SGS models (sigma model and data-driven S4ND model) produce extremely similar results at the first off-wall grid point, which is dominated by the wall model error. In this extremely coarse and extremely aggressive Reynolds number extrapolation case, the S4ND model is able to produce a profile that is closer to the log-law profile than the sigma model. From all this, the neural network is seen to provide reasonable results even when the Reynolds number is extrapolated by over 500 000 times the Reynolds number of the training set when considering mean profiles. This is further corroborated in figures 23(a), 23(b) and 23(c), where the Reynolds stresses also collapse reasonably well when scaled with friction velocity. There is no significant performance degradation visible at first glance as the Reynolds number increases, suggesting that the S4ND model remains stable even at Reynolds numbers significantly larger than the training set Reynolds number. While exact comparisons cannot be made as there is no DNS at these higher Reynolds number cases ( $\textit{Re}_2-\textit{Re}_5$ ), if the neural network is unable to extrapolate, as seen in Cho et al. (Reference Cho, Park and Choi2024), the change in shape of the profiles of turbulent statistical quantities is quite obvious to see.

Figure 23.

Reynolds stress profiles for various extreme Reynolds numbers for the S4ND model. (a) Here $R_{11}$ varying as a function of wall normal distance. (b) Here $R_{13}$ varying as a function of wall normal distance. (c) Here $R_{33}$ varying as a function of wall normal distance.

Figure 23. Long description

Three line graphs depict Reynolds stress profiles for various extreme Reynolds numbers for the S4ND model. Panel A: A line graph shows R11 varying as a function of wall normal distance. The x-axis is labeled z/H and ranges from 0 to 1. The y-axis is labeled R11/(uτ)^2 and ranges from 0.5 to 2.5. The graph includes five lines representing different Reynolds numbers (Re1, Re2, Re3, Re4, Re5) with a legend indicating each line’s color. The lines show a decreasing trend as the wall normal distance increases. Panel B: A line graph shows R13 varying as a function of wall normal distance. The x-axis is labeled z/H and ranges from 0 to 1. The y-axis is labeled R13/(uτ)^2 and ranges from −0.4 to 0. The graph includes five lines representing different Reynolds numbers (Re1, Re2, Re3, Re4, Re5) with a legend indicating each line’s color. The lines show a negative peak around the middle of the wall normal distance and then increase towards zero. Panel C: A line graph shows R33 varying as a function of wall normal distance. The x-axis is labeled z/H and ranges from 0 to 1. The y-axis is labeled R33/(uτ)^2 and ranges from 0.05 to 0.25. The graph includes five lines representing different Reynolds numbers (Re1, Re2, Re3, Re4, Re5) with a legend indicating each line’s color. The lines show an increasing trend, peaking around the middle of the wall normal distance, and then slightly decreasing.

In order to provide first-order analysis on the reasonableness of the channel flow scaling results, the implied skin friction coefficient, $C_f$ is plotted, as well as the expected inverse $log(Re)$ scaling that is supported by Schultz & Flack (Reference Schultz and Flack2013). These results are shown in figure 24. As seen, the implied $C_f$ is decreasing for the first four Reynolds numbers, and is still reasonably following the overall trend. However, the last Reynolds number suggests an increase in $C_f$ , suggesting that the results are less trustworthy at the highest Reynolds number. However, it is very likely this is more of a wall modelling error as opposed to a subgrid stress model error, as in these really coarse simulations there are eight points across the channel and the wall model’s accuracy is suspect at these extreme Reynolds number cases ( $\textit{Re}_{\tau }$ of 500 000 000). This argument is further enhanced when the sigma model results are taken into account, as it exhibits similar trends as compared with the S4ND model. Therefore, the S4ND model likely extrapolates well to these extremely out of training set Reynolds numbers, and the error that is seen can be attributed to the non-data-driven wall model. Even if this largest Reynolds number data point is discarded, the S4ND model still obtains a 500 times friction Reynolds number scaling when the model is only trained on one friction Reynolds number. The largest Reynolds number data point is mainly kept to show that even if the results are less trustworthy, the subgrid stress model still remains stable even in this sort of highly extrapolatory situation.

Figure 24.

Implied wall skin friction coefficient versus a first-order trend, with the x-axis being friction Reynolds number. The trends show that the error is likely due to the wall model, as both the S4ND model (plotted in blue) and the sigma model (plotted in red) show similar trends and deviation from the expected trend (plotted in yellow (Schultz & Flack Reference Schultz and Flack2013)).

While another way to ensure that the neural networks are stable during Reynolds number scaling is to use clipping, where values of the SGS tensor are set to zero if the local dissipation is negative. While this is an effective method, it is ultimately unphysical. In all these cases, no clipping is utilized and the S4ND model is stable without clipping.