Most of nowadays CFD simulations are conducted with the Reynolds Averaging approach. Reynolds Averaged Navier-Stokes (RANS) simulation, the “working horse” of industrial CFD is based on the Reynolds decomposition according to which a flow variable is decomposed into mean and fluctuating quantities. When the decomposition is applied to Navier-Stokes equation an extra term known as the Reynolds Stress Tensor arises and a modelling methodology is needed to close the equations. The “closure problem” is apparent as higher and higher moments of the set of equations may be taken, more unknown terms arise and the number of equations never suffices. This is of course an obvious consequence to the fact that taking these higher moments is simply a mathematical endeavor and has no physical contribution what so ever.
Levels of modeling are related to the number of differential equations added to Reynolds Averaged Navier-Stokes equations in order to “close” them.
0-equation (algebraic) models are the simplest form of turbulence models, a turbulence length scale is specified in advance through experimenting. 0-equations models are very limited in applications as they fail to take into account history effects, assuming turbulence is dissipated where it’s generated, a direct consequence of their algebraic nature.
1-equation and 2-equations models, incorporate a differential transport equation for the turbulent velocity scale (or the related the turbulent kinetic energy) and in the case of 2-equation models another transport equation for the length scale (or time scale), subsequently invoking the “Boussinesq Hypothesis” relating an eddy-viscosity analog to its kinetic gasses theory derived counterpart (albeit flow dependent and not a flow property) and relating it to the Reynolds stress through the mean strain.
In this sense 2-equation models can be viewed as “closed” because unlike 0-equation and 1-equation models (with exception of 1-equations transport for the eddy viscosity itself as described in a former post: Understanding The Spalart-Allmaras Turbulence Model) these models possess sufficient equations for constructing the eddy viscosity with no direct use for experimental results.
2-equations models do however contain many assumptions along the way for achieving the final form of the transport equations and as such are calibrated to work well only according to well-known features of the applications they are designed to solve. Nonetheless although their inherent limitations, today industry need for rapid answers dictates CFD simulations to be mainly conducted by 2-equations models whose strength has proven itself for wall bounded attached flows at high Reynolds number (thin boundary layers) due to calibration according to the law-of-the-wall.
A drawback evident in almost all eddy-viscosity models is the inability to inherently account for rotation and curvature. This drawback is resulted from relating the Reynolds stress to the mean flow strain and in fact is the major difference between such a modeling approach and a full Reynolds-stress model (RSM). The RSM approach accounts for the important effect of the transport of the principal turbulent shear-stress. On the other hand, RSM simulations are not computationally cost-effective, in as much that one does get an improved physical fidelity that is worth the time and computational resources consumed, not only that, they often do not converge.
Shortcomings In Wall Treatment for the Standard k-ε Turbulence Model
Although being perhaps the most popular and known turbulence model, the standard k-ε turbulence model carries along some harsh shortcomings which are important to acknowledge.
First, it is important to note that the model is essentially a high Reynolds model, meaning the law of the wall must be employed and provide velocity “boundary conditions” away from solid boundaries (what is termed “wall-functions”). From a mathematical standpoint, even if one could impose Dirichlet conditions for ε on solid boundary, after meshing it would still be difficult to numerically approach the problem due to what is termed in numerical analysis as stiffness of the numerical problem, partially related to the high gradients.
In order to integrate the equations through the viscous/laminar sublayer a “Low Reynolds” approach must be employed. This is achieved as additional highly non-linear damping functions are needed to be added to low-Reynolds formulations (low as in entering the viscous/laminar sublayer) to be able to integrate through the laminar sublayer (y+<10). This again produces numerical stiffness and in case is problematic to handle in view of linear numerical algorithms.
Furthermore the low Reynolds methodology should not be confused with transitional Reynolds modeling, as sometimes many practitioners rely on low Reynolds models to achieve a measure of transition from laminar to turbulence prediction. As the low Reynolds methodology is devised to handle the near wall viscous/laminar sublayer, there is no reason to expect it also satisfactory transition predictions, especially noting that there are many mechanisms for transition onset. Perhaps the only predictions which shall be close to satisfactory by low Reynolds model (maybe “pseudo-transition” behavior) are the ones related to bypass transition of which high levels of turbulence in the free stream occur and transition is dominated by diffusion effects (a brief description of the mechanism in my former post: A Forest of Hairpins – on the quest for turbulence coherent structures )
Another major drawback is the model lack of sensitivity to adverse pressure-gradient. It was observed that under such conditions it overestimates the shear stress and by that delays separation. Menter’s k-ω SST alleviates this drawback through the Shear Stress Transport (SST) concept. This drawback is evident in almost all eddy-viscosity models, relating the Reynolds stress to the mean flow strain and in fact is the major difference between such a modeling approach and a full Reynolds-stress model (RSM). The RSM approach accounts for the important effect of the transport of the principal turbulent shear-stress. The ingenious idea of Menter to include it in the revised k-ω model (termed the Baseline (BSL) model) is related observed success in implementing what is termed as the Bradshaw’s assumption, that the shear-stress in the boundary layer is proportional to the turbulent kinetic energy (a brief evaluation of the model in former post: Understanding The k-ω SST Model ).
Remedies in ANSYS Fluent for Near-Wall Treatment
Upon approaching a wall bounded turbulence flow RANS simulation a decision for the level of fidelity by which the turbulent boundary layer needs to be resolved needs to be made. The apparent question at hand is to whether or not resolving the viscous sublayer (y+<5-10) is mandatory to capture the physical phenomena at hand.
Commonly encountered situations where the flow demands integrating through the viscous sublayer include: heat transfer, calculation of aerodynamic drag, flows with adverse pressure gradients, etc’…
Integrating through the viscous sublayer would generally impose resolving the mesh to y+<1.
Typical situations for which the flow does not demand integrating the solution through the viscous sublayer are cases where wall-bounded effects are secondary, or where the flow separates to a dominant shear flow due to geometry such as the case of bluff body surrounding flow. For such cases the first cell away from the wall shall be placed at a y+ of 30-300.
As mentioned above unlike the standard ε-equation, the ω-equation can be integrated through the viscous sublayer without the need for a two-layer approach. This feature can be utilized for a y+-insensitive wall treatment by blending the viscous sublayer formulation and the logarithmic layer formulation based on y+. This formulation is the default for all ω. Such is also the case for the one equation eddy-viscosity transport Spalart-Allmaras turbulence model. While the k-ω model does that by taking an “elliptic” near wall behavior (partial differential wise), meaning that it has an inherent nature of being able to “communicate” with the wall and actually has Dirchlet (as in no-slip in this case) boundary conditions, The Spalart-Allmaras does that by indirectly incorporating the distance from the wall. See “Understanding The k-ω SST Model” and “Understanding The Spalart-Allmaras Turbulence Model” for in-depth review.
As the above also suggests, this is not the case for the k-ε model and traditionally, there are two approaches to modeling the near-wall region. In the first, the viscous sublayer and buffer layer (viscosity-affected inner region) are not resolved. Instead, semi-empirical formulas called “wall functions” are used to bridge the viscosity-affected region between the wall and the fully-turbulent region. The use of wall functions obviates the need to modify the turbulence models to account for the presence of the wall.
In the second approach, the turbulence models are modified to enable the viscosity-affected region to be resolved all the way to the wall, this is termed “Low-Reynolds Modeling”.
The main shortcoming of standard wall functions treatment is that the numerical results deteriorate under refinement of the grid in the wall normal direction. y+ values less than O(10) will gradually result in unbounded errors in wall shear stress and wall heat transfer.
On the other hand, the traditional low Reynolds approach demands for the first grid point away from the wall to be located at y+<1. For cases of which y+>1, the low-Re formulation results in inaccurate wall values for the wall shear stress (and heat transfer). This has is a major drawback on simulations demanding an integration through the viscous sublayer by such a restriction which depicts itself in highly unwanted computational demand that stems from the fact that for most industrial uses, complex geometries are involved.
Moreover, even when a y+<1 is achieved, lets say for a geometry which allows such control of the mesh, the methodology by which it is achieved induces numerical outcomes which are numerically stiff as the ε-equation can’t be integrated through the viscous sublayer to the wall and additional highly non-linear damping functions are needed to construct low-Reynolds formulations (low as in entering the viscous/laminar sublayer) that are able to bridge that gap.
To alleviate the above restriction ANSYS Fluent proposes a few remedies that circumvent the above presented shortcomings.
In what follows, it should be noted that in ANSYS Fluent, the law-of-the-wall for mean velocity and temperature are based on the wall unit, y*, rather than:
This means that as the law-of-the-wall for mean velocity yields:
where U* is dimensionless velocity defined as:
The dimensionless distance from the wall, y*, may be defined as:
It is further important to note that in ANSYS Fluent the log-law is employed when y*>11.225. For values below 11.225 a laminar stress-strain relationship is employed (i.e. y*=U*).
ANSYS Fluent proposes the following near wall
- Standard wall functions:
The standard wall functions in ANSYS Fluent are based on the Launder and Spalding route.
Wall functions are applicable for a suitable range of y* and by that to the flow’s Reynolds. As the purpose is to allow not integrating through the viscous sublayer, the lower limit of y* is 11 hence standard wall functions should not be used below that limit as the solution’s accuracy might deteriorate in an uncontrolled manner.
The upper limit, as the range of the integral layer in the law-of-the-wall depends on the Reynolds number on a much more physical level. This conceptual understanding of the law-of -the wall is important and should be confronted upon approaching a wall-bounded turbulence problem. For very high Reynolds numbers (such as submarines, ships, airplanes, etc’…) the logarithmic layer can extend to values as high as several thousand, while for low Reynolds applications (such as UAVs, turbine blades, etc’…) the upper limit for the integral layer may be as low as 100. Moreover, for the later, it is much less profound, meaning that the level by which the law-of the-wall and the integral layer could be recognized is unclear.
The last notion implies that the wall functions approach (especially in the standard form) is most suitable for high Reynolds flows of which integrating through the viscous sublayer has only a slight impact on the desired outcome.In k-ε models, the turbulence kinetic energy is solved in the entire domain under the boundary condition zero normal turbulence kinetic energy flux.
The source terms for the kinetic turbulence are its production Gk ,and its dissipation rate, ε. Those are computed by the assumption that they are equal in the near wall control volume – this is termed the local equilibrium hypothesis.Following the above explanation, the production of the turbulence kinetic energy is based on the log-law and may be computed as:and the dissipation equation not solved at the wall-adjacent cells, but instead is computed using:
The standard wall functions work reasonably well for a broad range of wall-bounded flows – the flows inside the range of applicability of constant-shear and local equilibrium assumptions that do not incorporate physical phenomena of interest such as flows involving separation, reattachment, and impingement where the mean flow and turbulence are subjected to pressure gradients and rapid changes.
Scalable wall functions:
Industrial problem often contain complex geometries making the mission of defining consistently the first grid point away from the wall at y*>11 quite a conundrum. To avoid deterioration of the solution by standard wall functions in situations where it’s unavoidable for the first grid point to be located at y*<11, ANSYS Fluent proposes wall functions that produce consistent results for arbitrary grid refinement by forcing the usage of the log law in conjunction with the standard wall functions approach:
- Non-equilibrium wall functions:
As mentioned above standard wall functions will not work well outside their range of applicability (i.e. constant-shear and local equilibrium assumptions).
For flows involving separation, reattachment, and impingement where the mean flow and turbulence are subjected to pressure gradients and rapid changes, improvements can be obtained, particularly in the prediction of wall shear and heat transfer characteristics. This is achieved in ANSYS Fluent by means of sensitizing the log-law for mean velocity to pressure gradient effects and by the use of the two-layer-based concept to compute the reciprocal relations between turbulence kinetic energy production, Gk, and its dissipation ε, by non-equilibrium means. The log-law mean velocity sensitized mean velocity to pressure gradients is achieved as follows:
the physical viscous sublayer thickness is computed from:
The sensitized mean velocity may be drawn from:
The wall-neighboring cells are assumed to consist of a viscous sublayer and a fully turbulent layer in order to employ the two-layer concept in computing the reciprocal relations between turbulence kinetic energy production, Gk, and its dissipation ε, by non-equilibrium means to solve the turbulence kinetic energy equation at the wall-neighboring cells.
To be able to that the following profile assumptions for turbulence quantities are made:
What is done next is the computation of cell-averaged production of turbulence kinetic energy and cell-averaged dissipation of turbulence kinetic energy from the volume average of Gk and ε of the wall-adjacent cells. These proportions between viscous sublayer and the fully turbulent layer which varies widely from cell to cell in highly non-equilibrium flows actively determine
The turbulence kinetic energy budget for the wall-neighboring cells.
Therefore to some extent the non-equilibrium formulation for the wall functions takes into account the effect of pressure gradients on the distortion of the velocity profiles and by that account for some non-equilibrium effects that are not accounted for in the standard formulation due to the local equilibrium assumption.
Enhanced Wall Treatment ε-Equation (EWT-ε):
If the near-wall mesh is fine enough to be able to resolve the viscous sublayer (typically with the first near-wall node placed at y+<1), then the enhanced wall treatment will be identical to the traditional two-layer zonal model.
However, as noted above the restriction that the near-wall mesh must be sufficiently fine everywhere be limiting from a computational resources standpoint.
The best scenario by concept is to have a near wall treatment suitable for viscous sublayer integration when the mesh is at y+<1 and to use wall functions formulation when the first grid away from the wall falls in the log-layer and the viscous sublayer solution is of minor consequence. Moreover, the scenario should prevent excessive error to occur for the intermediate meshes where the first near-wall node is placed neither in the fully turbulent region, where the wall functions are suitable, nor in the direct vicinity of the wall at y+<1, where the low-Reynold approach is adequate.
This is achieved in ANSYS Fluent by incorporating the two-layer model with enhanced wall functions.
The separation of the two regions is conducted via a wall-distance-based, turbulent Reynolds number:
Where y is the distance to the nearest wall and defined uniquely for complex shaped domains and independent on mesh topology by applying the following on the union of all boundaries involved:
For regions where the above defined Reynolds number is above 200 the original standard k-ε is employed.
By the same token if the above Reynolds number is below 200, a one equation for the transport of turbulence kinetic energy is employed (Wolfstein’s k-equation), while the turbulent viscosity is computed from:
The length scale in the above:
The above intention is two blend the enhanced wall treatment eddy viscosity with that of high Reynolds at the outer region by:
the blending functions is such that it is 1 away from the wall and 0 in their viscinity:
In the viscosity affected region the 1-equation for the turbulence kinetic energy is completed by algebraically computing the dissipation, ε, as:
The length scale in the above (Chen & Patel):
The constants in the above are achieved by calibration to be:
In order for the method’s applicability to be extended to the entire near wall region Fluent also blends the different near wall regions to a single wall-law formulation. One outcomes is that this approach allows the fully turbulent law to be easily modified and extended to take into account other effects such as pressure gradients and variable properties due to compressibility or heat transfer.
Menter-Lechner ε-Equation (ML-ε):
In the above methodology it was shown how the two layer model may serve as a basis for a y+ insensitive model.
We also recall a definition of a wall-distance-based, turbulent Reynolds number to clearly separate two regions:
The use of such a definition hides some shortcomings within. One of its major ones is that regions with small turbulence kinetic energy might also have a turbulent Reynolds number below 200 and as such treated with a near-wall formulation even though they are actually away from the wall.
The Menter-Lechner ε-Equation (ML-ε) approach offers a y+ insensitive which is a not based on a two layer model assuming a sufficient resolution of the boundary layer, predicts in a y+ independent manner the wall shear stress and wall heat flux. Such a formulation switches gradually from wall functions to a low-Reynolds formulation when the mesh is refined.
The Menter-Lechner ε-Equation (ML-ε) concept is achieved by adding a a source term to the transport equation of the turbulence kinetic energy that accounts for near-wall effects:
This source term (marked in red) is active in the viscous sublayer and becomes zero in the logarithmic layer.
This closes the near wall treatment for what generally be referred to as ε-equation models.
To conclude, RANS shall take part in near wall modeling for long time to come. Large-Eddy Simulation (LES) allows high fidelity simulations for many applications but has severe limitations in the near wall regions, as the computational effort required to reliably model the innermost portion of the boundary layer (sometimes constituting more than 90% of the mesh) where turbulence length scale becomes very small is far from the resources available to the industry. Anecdotally, best estimates speculate that a full LES simulation for a complete airborne vehicle at a reasonably high Reynolds number will not be possible until approximately 2050…
Nevertheless even with the enormous reduction in nowadays computing resources limitations encountered with RANS it should be noted that the practitioner should always strive to gain more knowledge on the physical phenomena at hand before clearly preferring one methodology over the other.
Ref: ANSYS Fluent Theory Guide