A Priori Assessment of Subgrid-Scale Models and Numerical Error in Forced Convective Flow at High Prandtl Numbers

An extensive a priori analysis has been carried out on data from Direct numerical simulation of fully developed heated turbulent pipe flow at high molecular Prandtl numbers Pr=10\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Pr=10$$\end{document}/20, testing three popular modelling candidates for subgrid-scale closure in Large-Eddy simulation (LES). Aside from assessing the models’ capabilities to describe quantitatively the unresolved turbulent fluxes, a special focus is also put on the role of the numerical error, which arises from the discretization of the filtered advective fluxes on a coarse LES grid. The present analysis extends here previous studies on subgrid-scale momentum transport in a isothermal mixing layer and channel flow carried out by Brandt (J Numer Methods Fluids 51: 635–657, 2006) and Vreman et al. (J Eng Math 29: 299–327, 1995), respectively, to the subgrid-scale transport of heat at high Prandtl numbers. The statistical dependence between the individual contributions (resolved, subgrid-scale, numerical discretization error) constituting the filtered advective flux divergence in the LES formulation is investigated as well, in terms of corresponding cross-correlations. The sensitivity of the tested sgs-models to a grid refinement is further examined performing also a posteriori LES, where the basically more sophisticated candidates turn out to be more demanding in terms of required grid resolution.


Introduction
The accurate modelling of the subgrid-scale (sgs) fluxes of momentum and heat near heated/cooled walls still represents a challenging task in Large-Eddy Simulations (LES). This is specially the case at molecular Prandtl numbers (Pr) far from unity, where the relevant length scales of the smallest turbulent thermal structures differ significantly from those of the smallest dynamic structures. This disparity strongly conflicts with the Reynolds analogy between the transport of momentum and heat, which is inherently assumed by most sgs-models for the turbulent heat flux. Using the popular Boussinesq-type diffusion Ansatz for the unresolved advective subgrid-scale fluxes, the Reynolds analogy relates 1 3 the eddy-viscosity directly to the thermal eddy-diffusivity through a subgrid-scale Prandtl number Pr sgs being of the order of unity. The near-wall variation of its counterpart, which is equivalently used in Reynolds Averaged Navier-Stokes (RANS), the turbulent Prandtl number Pr T , has been extensively investigated in various studies based on theoretical analysis and/or Direct Numerical Simulation (DNS), see, e.g. Antonia and Kim (1991), Kawamura et al. (2000), Piller (2005) and Reynolds (1975). Comparatively fewer studies looked at Pr sgs , considering mostly fluids with molecular Prandtl numbers near unity, like the work of Inagaki et al. (2012) and Moin et al. (1991). The adequate setting of this important model parameter for the considered high molecular Prandtl numbers shall be also addressed in the present work.
The principal capability of providing closure for the unresolved advective momentum fluxes, has been widely tested in many a priori LES for various subgrid-scale (sgs) models, see, e.g., Abbá et al. (2003), Jimenez et al. (2001), Moin et al. (1987) and Salvetti and Banerjee (1995). Most of these a priori tests evaluate the individual sgs models by comparing their predictions for the subgrid-scale fluxes against those obtained from an explicit filtering of the non-linear advective flux-terms, as provided by the fully resolved DNS solution. In fact, the subgrid-scale closure must rather describe the divergence of these fluxes, as appearing in the filtered equations of motion. Therefore, aside from the accurate modelling of the sgs-fluxes, their substitution into the filtered equations also includes a numerical error arising from the discretized representation of the divergence operator on the coarse LES grids. This often disregarded issue was addressed in only a few a priori studies, such as the works of Vreman et al. (1995) and Brandt (2006), who examined the contribution of the numerical error in the filtered momentum equation. For the case of a turbulent mixing layer (Vreman et al. 1995) demonstrated that up to a certain threshold of filter width the instantaneous numerical error statistically exceeds the sgs contribution to be provided by the model. This excess still appeared for filter-widths equal to the grid size of the coarse LES mesh, as typically assumed in LES with implicit grid-filtering, which makes the role of the actually applied sgs-modelling questionable. The a priori LES of Brandt (2006) analysed the individual contributions to the divergence of the momentum fluxes, as occurring in a set of explicitly filtered LES equations. The study demonstrated that applying the explicit filtering only to the non-linear advective flux terms instead of the full velocity field of the LES, produces plausible tendencies for increasing filter-width, in that the sgs-contribution increases, while the numerical error decreases. The present work shall extend the investigations of Brandt (2006) and Vreman et al. (1995) to the convective heat transport, considering fully developed turbulent pipe flow at molecular Prandtl numbers well beyond unity. The study shall further examine, how the invidual contributions occurring in the flux divergence terms are statistically cross-correlated. The effect of the resolution provided by the coarse LES grid shall be investiged as well. Finally, the findings from the a priori LES are backed up with results from a posteriori LES, which have been carried out for selected cases.

Filtered Governing Equations
For the considered incompressible flow with constant material properties, the spatially filtered non-dimensionalized conservation equations of mass, momentum and energy solved in LES read The pipe diameter D, the wall friction velocity w = √ w ∕ and temperature T = q w ∕ c p w have been used as length, velocity, and temperature scales, respectively, for non-dimensionalization. = (ū,v,w) T and represent the non-dimensional filtered velocity vector and temperature difference from the wall value, respectively. The overlined operator ∇ denotes the gradient vector computed on the LES grid. sgs and sgs represent the unresolved subgrid-scale (sgs) stress tensor and heat flux vector, respectively, arising from the spatial filtering of the non-linear advection terms. They are accordingly defined as For closing the set of LES equations the deviatoric part of sgs stress tensor and the sgs heat flux vector are modelled with the Boussinesq hypothesis as respectively. The subgrid-scale eddy viscosity sgs occurring in (6) is provided by the adopted sgs model. It is further used in (7) for modelling subgrid-scale thermal eddy diffusivity introducing the subgrid-scale Prandtl number, which is presently assumed as constant set to Pr sgs = 0.5.

A Priori LES
The a priori LES is based on fully resolved instantaneous velocity, U, and temperature, , fields obtained from DNS of fully developed heated turbulent pipe flow, considering four cases at wall friction Reynolds number Re = 360∕500 and Pr = 10∕20 . The DNS results are computed in cylindrical coordinates on a mesh of 256 × 512 × 1024 grid points in the radial r, azimuthal , and axial directions z, respectively. The total axial length of the computational domain is five pipe diameters. A 4th order accurate Finite-Volume scheme and a 2nd order explicit Adam-Bashforth scheme have been used for spatial and temporal discretization, respectively. The wall of the pipe is uniformly heated with a constant average heat flux ⟨q w ⟩ = const. , assuming no-slip wall boundary conditions for velocity, w = 0 , and Dirichlet-type thermal boundary conditions w = 0 . (1)

3
Periodic boundary conditions are applied in the axial and the azimuthal directions. The DNS data have been filtered with a spatial box filter, using a filter width ̄= 8 , which is eight times the cell size used in the DNS. This implies a spatial resolution for the a priori LES obtained on a mesh with 32 × 64 × 128 grid points in the radial r, azimuthal , and axial direction z, respectively. The present analysis splits the spatially filtered advective divergence terms, that arise from the filtering of the non-linear terms in the transport equations of momentum and heat, into three contributions, being the resolved component, the unresolved subgridscale component, and the numerical error due to the discretization on the coarse LES grid. Accordingly, the filtered divergence of the advective momentum and heat fluxes are decomposed into respectively. The vector and scalar representing the numerical error are estimated as proposed by Vreman et al. (1995), assuming the DNS solution on the fine DNS-grid as representative for the exact solution, such that where the operator ∇ refers to the fourth-order accurate discretization used on the fine DNS-grid, while ∇ refers to the second-order accurate discretization used on the coarse LES-grid. The subgrid-scale contribution to be delivered by the sgs-model is consistently computed on the coarse grid, applying alternatively the assessed candidates, namely the Standard Smagorinsky model with Van Driest-type wall dampening (SSM) (Smagorinsky 1963;Van Driest 1956), the Wall-adapting Local Eddy-viscosity model (WALE) (Nicoud and Ducros 2006), and the Coherent Structure Model (CSM) (Kobayashi 2005). The Standard Smagorinsky model (SSM) computes the eddy viscosity as using a constant model parameter set to C s = 0.1 , which is reduced to zero with decreasing non-dimensional wall distance y + applying a Van-Driest-type damping function. The Coherent Structure model (CSM) represents a Smagorinsky model using a non-constant the model parameter modelled as where the coherent structure function is computed from (1 − F cs )|F cs | 1.5 , 1 3 dependent on the resolved strain tensor S defined in (6) and the resolved vorticity tensor W = 1 2 (∇ − (∇ ) T ) . The Wall-adapting Local Eddy-viscosity (WALE) approach models the sgs eddy viscosity as dependent on the traceless part of the square of the velocity gradient tensor and model constant set to C w = 0.35 . The filter width is always defined as ̄= (̄1̄2̄3) 1∕3 . Neither of the assessed Boussinesq-based sgs-models requires any explicit test-filtering procedure for determining dynamically the model parameters, nor any subsequent statistical averaging of these parameters is needed for the sake of numerical stability in a posteriori LES. Since the coherent structure function F cs as well as the tensor ̄d go to zero near solid walls, CSM and WALE provide the additional benefit of not requiring any artificial wall damping functions to enforce a vanishing eddy viscosity in the viscous near-wall region. These features make the tested models attractive candidates in complex technical engineering applications far apart from simplified generic test flow configurations. The a priori representations of the eddy viscosity and thermal eddy diffusivity are computed from the filtered DNS-data by contracting Eqs. (6) and (7) with the filtered strain tensor S and temperature gradient ∇̄ , respectively , such that and The angular brackets always denote statistical averaging in the homogeneous spatial directions z and .

Results
The spatial resolution provided by the DNS for the smallest dynamics and thermal structures, measured in terms of the Kolmogorov and Bachelor length scales, K = ( 3 ∕ ) 1∕4 and B = K ∕Pr 1∕2 , respectively, is shown in Fig. 1. The depicted profiles for the spatial resolution in radial, azimuthal and axial directions refer to the most demanding case with Re = 500 and Pr = 20 . Due to the applied cylindrical mesh, the resolution is coarsest in the azimuthal direction near the wall. The observed maximum levels are well comparable to those presented in recent DNS by Lee et al. (2013), Nemati et al. (2016) and Zonta et al. (2012).
As follows from Eq. (7), the presently assessed sgs-models compute the thermal eddydiffusivity as directly related to the modelled eddy-viscosity through a subgrid-scale Prandtl number Pr sgs , which requires a suitable setting of this parameter. Figure 2 shows the radial variation of the local a priori value for Pr sgs = sgs ∕a sgs , as obtained from Eqs. (18) and (19) using the filtered DNS results. For all considered cases, Pr sgs remains evidently close 0.5 in the fully turbulent core region. It increases towards unity in the buffer layer and further approaches higher levels above unity very next to the wall ( y + < 1 ). As such the a priori results are well in line with the findings for channel flow presented in literature (Moin et al. 1991), and, as particularly suggested by the observations for the turbulent inner region, they support presently applied setting Pr sgs = 0.5.

P r sgs
Re τ =360 P r =20 10 −1 10 0 10 1 10 2 y + Re τ =500 P r =20 Fig. 2 Radial variation of the a priori Pr sgs = sgs ∕a sgs computed from the filtered DNS data at considered Re and Pr cases Figure 3 shows the radial variations of the statistically averaged budgets of the advective heat flux divergence for the four considered cases, as obtained from the decomposition (9). A pronounced negative peak always occurs near the upper limit of the viscous sublayer. The peak increases in magnitude for increasing Reynolds and Prandtl numbers, and it is located closer to the wall for the higher Prandtl number.
The numerical error evidently adds on average very little to the resolved contribution. The remaining considerable gap to the blue (DNS-based) target line is filled and partly overcompensated only by the Standard Smagorinsky subgrid-scale model (SSM), while the WALE model delivers a comparatively much smaller contribution and the CSM model produces even less. As follows from Eq. (7), the contribution from the subgrid-scale model is certainly also affected by the choice of the subgrid-scale Prandtl number. Being assumed as constant Pr sgs = 0.5, it falls below the a priori values in the near-wall region as already discussed above in Fig. 2. Accounting for this observed near-wall increase would effectively further reduce the contribution from the sgs-model. However, since the variation of Pr sgs affects only a narrow region very next to the wall, the effect such an ad-hoc modification on the advective heat flux budgets is still minor.
In contrast to the small average contribution observed in Fig. 3, the instantaneous fluctuating representations of the numerical discretization error are signficantly higher. This is clearly seen from Figs. 4 and 5, showing exemplarily rms value of the error contribution to the axial momentum transport w and the error contribution to the heat transport , respectively. The rms-values of the numerical error evidently exceed markedly their statistical means, as shown in Fig. 3 for the heat transport, and they tend to increase for higher Reynolds and Prandtl numbers. The also shown lower black lines represent the rms-values of the a priori sgs contributions w = ∇ ⋅ w sgs and = ∇ ⋅ sgs occurring in the budgets (8) and (9), respectively, where the a priori sgsflux vectors w sgs and sgs are obtained from Eqs. (4) and (5), respectively, using the filtered DNS data for velocity and temperature. The instantaneous contributions from the a priori sgs-terms evidently tend to stay considerably below those of the numerical error. This observation confirms the findings of Vreman et al. (1995). Their a priori analysis of a turbulent mixing layer unveiled the same dominance of the numerical error over the sgs-contribution, if the filter width is equal to the mesh size of the LESgrid, which is generally assumed in implicit LES and also applies in the present study.
The mutual statistical dependence of the individual contributions is further investigated in terms of the normalized cross correlations, as shown for the advective heat flux for the considered Reynolds and Prandtl numbers in Figs. 6, 7, 8 and 9. The correlation coefficient between the resolved contribution and the numerical error becomes strongly negative towards the wall. This indicates the dissipative nature of the numerical error, which effectively dampens the resolved turbulent fluctuations in this region. On the other hand, the resolved contribution mostly appears as positively correlated with the contribution from subgrid-scale model in the diffusive sub-and buffer layer y + < 10 , implying that the applied subgrid-scale model tends to enhance the resolved turbulent fluctuations in this region. This positive correlation is obviously most pronounced for the SSM, which also explains the overpredicted subgrid-scale contribution of SSM near the wall, as seen in the budgets in Fig. 3  Re τ =360 P r =10 β rms θ α rms θ Re τ =360 P r =20 10 0 10 1 10 2 y + 400 800 1200 Re τ =500 P r =10 10 0 10 1 10 2 y + Re τ =500 P r =20

Effect of Resolution
The effect of a different grid resolutions was examined in two ways. The first strategy consists in a homogeneous refinement of the filter width by a factor 2 in all directions leading to = 4 DNS . The corresponding LES grid contains a total number of 64 × 128 × 256 grid points in the radial r, azimuthal , and axial direction z, respectively. As such, this grid is 8 times larger in total size than the previous one associated with = 8 DNS . The second strategy consists in doubling the grid resolution in -direction, while decreasing the grid resolution in z-direction by a factor 2, keeping the radial resolution unchanged. As such, this grid contains a total number of 32 × 128 × 64 grid points, in the radial r, azimuthal , and axial direction z, respectively. This particular

Fig. 8
Cross correlation coefficient ⟨ ⟩ * between numerical error and resolved contribution (red line), numerical error and sgs-model contribution (blue line), resolved and sgs-model contribution (green line) for Re = 500 and Pr = 10 grid resolution, which will be referred to as ̄ * = 8 DNS in the further, provides an improved near-wall resolution of the azimuthal turbulent structures, but does not lead to a change in the local filter width , nor in the total mesh size with respect to the case = 8 DNS . Figure 10 shows the radial variation of the eddy-viscosities predicted by the sgsmodels on the different grids, compared against the corresponding a priori values computed from Eq. (18) using the filtered DNS data (indicated as "DNS" in the legend). Both the CSM and the WALE model produce evidently a markedly lower eddy-viscosity in the buffer layer and below ( y + < 10 ), while the wall-dampened prediction of SSM lies fairly close to the DNS-based a priori results. This deficit stays the same for all considered resolutions. It explains the insignificantly small sgs-contributions of CSM and WALE observed in Fig. 3 as well as in the Figs. 11 and 12, which exemplarily show the budgets of the advective flux divergence for the tested alternative grids at Re = 500 and  Pr = 10 . Using the homogeneously refined grid associated with a filter width = 4 DNS expectedly increases the resolved contribution to the advective flux divergence. As already indicated by the predicted eddy-viscosities shown in Fig. 11, owing to the reduced filter width, the contributions of the CSM and the WALE sgs-models remain consistently low, while the contribution of the SSM sgs-model is decreased in a way, such that the target line is closely reached, without the overshoot seen for SSM in Fig. 3 (right column). As can be seen from Fig. 12, using the azimuthally refined grid associated with * = 8 DNS does also notably improve the resolved component in the advection budgets, as compared to the LES grid with = 8 DNS shown in Fig. 3, third row. This effectively reduces as well the gap to the target line, while the sgs-model contributions essentially remain the same because of the unaltered filter width. Thus, adding to 10 0 10 1 10 2 y + (SSM) Fig. 12 Budgets of the averaged advective flux divergence for the heat transport comparing the WALE, CSM and SSM subgrid-scale model contributions for Re = 500, Pr = 10 for * = 8 DNS a higher resolved base contribution, the budgets with the CSM and WALE models come finally closer to the target line, whereas the SSM produces a considerable overshoot.

A Posteriori LES
The aforementioned results in the a priori analysis for = 8 DNS and * = 8 DNS were further assessed with a posteriori LES. The governing set of LES equations were solved using a finite volume scheme with 2nd order accuracy in time and space. WALE, CSM and SSM were applied as sgs-models. Two grid resolutions were used, corresponding exactly to those associated with the previously defined filter widths = 8 DNS and * = 8 DNS . Figure 13 exemplarily shows the predicted averaged velocity and temperature for Re = 500, Pr = 10 . The improved direct resolution of advection provided by the azimuthal refinement, which was already observed in the a priori analysis, evidently produces better predictions for all considered sgs-models. WALE and CSM clearly outperform the SSM approach. This is also seen in the turbulent fluxes for momentum and heat shown in Figs. 14 and 15, respectively. In the case with * = 8 DNS , the resolved contributions for CSM and WALE are almost on top of the DNS profiles in the buffer layer, and no significant sgs-contribution is observed here. On the other hand, the SSM appears to dissipate a great part of the turbulent motion, as indicated by the reduced resolved fluxes near the wall. Due to the higher predicted eddy-viscosity, SSM compensates for this loss with a higher sgs contribution, which partly closes the considerable gap to the DNS results for the total turbulent momentum flux. On the other hand, the higher eddy diffusivity, being directly related to the eddy viscosity through Pr sgs , evidently leads to an overprediction for the total turbulent heat flux near the wall. This overprediction is also well in line with the overshoot observed in the a priori budgets for SSM in Fig. 12, and it suggests a higher setting of Pr sgs for SSM.

Conclusions
The performance of three Boussinesq-type sgs-models was tested in a priori and a posteriori LES at molecular Prandtl numbers well beyond unity. Among the tested approaches, the CSM and the WALE model turned out to yield a much lower eddy-viscosity/diffusivity in the near wall region, as compared to SSM. This particular feature, which translates into very low contributions to the advective flux budgets, makes both approaches basically more sensitive to the spatial resolution of the LES grid. The comparatively more rapid decrease in eddy-viscosity/diffusivity predicted by CSM and WALE near the wall requires a suitably refined LES grid to resolve sufficiently the relevant turbulent structures, so that only a very small sgs-contribution needs to be supplied by the model. In contrast, producing generally a higher level of eddy-viscosity, the wall-dampened SSM appears to compensate better for the unresolved sgs-structures typically met with fairly coarse LES grids. However, as shown in the present a priori and a posteriori LES on refined grids, this benefit may turn into a deficit, where the higher levels of eddy-viscosity for SSM produce too much dissipation. While the contribution of the numerical error to the advective flux divergence is on average very small, it exhibits a very high instantaneous local variation. Its rms values are shown to strongly exceed the rms-values of the contributions from the a priori sgs-components obtained from the filtered DNS data. The observed increasing tendency with Reynolds and molecular Prandtl number underlines the importance of this contribution when considering high Prandtl number flow. The analysis of the statistical cross-correlations of the invidual contributions to the advective flux divergence unveils a strong negative correlation between the numerical error component and the resolved component near the wall. This indicates an effective dampening of the resolved turbulent motion by the numerical error.