Turbulence: a necessity!
- Sir Horace Lamb, 1932, in an address to the British Association for the Advancement of Science
"I am an old man now, and when I die and go to heaven there are two matters on which I hope for enlightenment. One is quantum electrodynamics, and the other is the turbulent motion of fluids. And about the former I am rather optimistic."
Despite such as remarkable and intriguing characteristics of "turbulent motion", there is no formal definition of turbulence and it has been characterized intuitively using terms like turbulent kinetic energy, turbulent dissipation rate and so on. One of the best known and a bit poetic definition is by Richardson which explain the energy cascade from mean flow to viscous dissipation at walls.
Big whorls have little whorls,
which feed on their velocity;
And little whorls have lesser whorls,
And so on to viscosity.
Before moving on the necessity of turbulence modeling, few words on the technologies involved in CFD simulations. Any numerical simulation program has to deal with governing physics and its numerical applications, conversion of partial differential equations into equivalent matrix form and using computer to let users visualize the inputs and outputs. The picture below summarizes the inter-dependence of technology and applied mathematics. "Turbulence modeling" is one of the most critical steps in overall simulation process: be it development of the program itself or application to real-life problems.
This page deals with the concept of turbulence and how it is mathematically represented in CFD simulation programs. Walls are main source of vorticity and turbulence and its presence gives rise to turbulent momentum and thermal boundary layers: accurate prediction of frictional drag for external flows and pressure drop for internal duct flows depends on accuracy of local wall shear stress predictions.
In this context, the steep variations of field varaibls (velocity, temperature ad pressure) are in the very near-wall regions shown by following picture. The zone near the walls are known as boundary layer and has further been divided into sub-layers. The turbulence and boundary layer are two very closely related topics. While overall flow is divided into two zones: boundary layer and free-stream, boundary layers themselves are divided into 4 zones: viscous sub-layer, buffer-layer, log-law region and outer layer.
Use of very fine mesh to resolve these steep profiles is in most of the applications computationally too expensive for application of CFD tools for industrial processes. Hence special near-wall treatments are necessary since governing equations cannot be integrated down to wall and thus the concept of wall functions and near wall treament came into existence.
But turbulence is not an undesirable thing under all circumstances. The image describes the features of "turbulent motions" an well as few of the benefits too.
The complicated motion of turbulence is the result of the non-linear advection that creates interactions between different scales of motion. The large scales are governed by the geometrical scales and the smallest scale, known as Kolmogorov scale, is constrained by the viscosity of the fluid.
- Though all fluid flows are governed by Navier-Stokes equations, the wide range of length and time scales poses difficulties in treating turbulence, both analytically and numerically. Statistically large number of length scales also cause non-repeatibility where calculations tend to be sensitive to initial conditions and turbulence parameter values at the inlet boundary.
- The presence of too may flow scales computational resources too high for most industrial applications. To address this, the governing equations
are either time/spatially filtered to remove all/some of the turbulent scales.
- This approach results in a tractable set of equations for numerical calculations, closure problem arises from the correlated velocity fluctuations. So, a turbulence model is needed to represent the effect of the small scales on the mean flow.
- There are no turbulence model that fit all the requirements. k-ε model is the workhorse of the industry but needs to be used with caution. k-ω SST generally gives good results in most of the cases where k-ε fall short. Yet one must investigate the effect of different turbulence models for the particular case being investigated.
Last but not the least, an excerpts from Academic Resource Center of Illinois Institute of Technology: "The phenomenon of turbulence, caused by the convective terms, is considered the last unsolved problem of classical mechanics. We know more about quantum particles and supernova than we do about the swirling of creamer in a steaming cup of coffee!"
Key Parameters for Specification of Turbulence
- Energy cascade Cascade refers to a "transfer process". In turbulence flows, it refers to the transfer of kinetic energy from large scales of motion, through successively smaller scales, ending with viscous dissipation and conversion to heat.
- Ergodicity and Homogeneous Turbulence Ergodiity refers to the case when time averages and ensemble averages are equivalent. Homogeneous turbulence is when statistical properties do not change with spatial translation (turbulent behviour do not change with position).
- Auto-correlations, cross-correlations, intermittency and stochastic (variable whose auto-correlation decays to zero exponentially) are few other terms used to define turbulence and its attributes.
- Vorticity: This is neither a physical property of the fluid nor a directly measurable physical charateristics (like vvelocity of pressure) of the flow field. It is a mathematical definition where vorticity ω = ∇ × V.
- Turbulent (or Eddy) Viscosity - μ_{t}
- Analogous to molecular viscosity, turbulence viscosity is related to fluctuating component of the velocity. It is the constant of proportionality between turbulent (Reynolds) stresses and mean (large-scale) strain rate, analogous to laminar viscosity in Newton’s law of viscosity.
- μ_{t} = C_{μ} * k^{2} / ε where C_{μ} = 0.09, an empirical coefficient.
- Turbulent Kinetic Energy - k [m^{2}/s^{2}]
- It is specific (per unit mass) kinetic energy of eddies = 1/2 * v'^{2}
- This is trace of Reynolds stress tensor, computed by summing the diagonals of the matrix of components. Thus, TKE,
- It defines velocity scale of the turbulence (eddies). Remember Ryenolds number is ratio of length scales * velocity scale / kinematic viscosity
- Since u', v' and w' are not modeled in CFD except DNS approach, how do we estimate u' and k? The answer is Boussinesq hypothesis: which relates the Reynolds stresses to the mean velocity gradients. This hypotehsis is used in Spallarat-Allmaras model, k-ε and k-ω models.
- Boussinesq hypothesis: Turbulent viscosity is expressed based on analogy with molecular transport of momentum [viscous stresses are linearly proportional to mean strain rates]- hence assumptions is valid at molecular level, not necessarily always valid at macroscopic level.
- Turbulent Intensity - I
- The turbulence intensity I, is defined as the ratio of the root-mean-square of the velocity fluctuations u', to the mean flow velocity, u_{mean}.
- A turbulence intensity of 1% or less is generally considered low and turbulence intensities greater than 10% are considered high.
- Turbulent Eddy Dissipation Rate - ε = dk/dt [m^{2}/s^{3}]
- It is rate of conversion of k into heat or thermal energy per unit mass.
- Along with k, it defines length scale of the flows (turbulence eddies)
- Turbulent Eddy Dissipation Frequency - ω [m^{2}/s]
- It is rate of conversion of k into thermal energy per unit volume.
- It can be defined implicitly by k and ε
- Flow scales
- Length scale, η = (ν^{3} / ε)^{0.25}
- Time scale, &tau = (ν / ε)^{0.50}
- Velocity scale, v = (ν . ε)^{0.25}
- Summary: RANS methods
- Step-1: Mean flow velosity and pressure equations are solved.
- Step-2: TKE (k) and TED (ε) equations are solved.
- Step-3: Turbulent viscosity is solved using equations described above.
- Step-4: Reynolds stress tensor is modeled (estimated) using Boussinesq hypothesis.
Following plot described the production and dissipation of TKE. Note that k and ε are in almost equalibrium far away from the wall.
Key variable of turbulence modeling is Y-plus (y^{+}).
The limitations imposed on this parameter governs the size of mesh near the walls. This refers to two parameters:
- The height of the 1st layer of the mesh near the walls
- Number of layers in the boundary layer (which is not known a priori)!
Further deliberations need to be done along with the wall functions
- "High-Re number wall treatment" refers to case where boundary layers are modeled in terms of empirical (and well-established) correlations known as wall functions. Wall laws used to calculate wall shear stress.
- "Low-Re number wall treatment" refers to the case when boundary layers are resolved and not modeled. The Re-number here refers to Turbulent Reynolds Number defined as Re_{TURB} = k^{2}ν/ε.
- Note that it is related to the way boundary layers are resolved and has nothing to do with the "MEAN FLOW REYNOLDS NUMBER".
- A "low Reynolds number" turbulence model doesn't use wall functions, that is does not involve any assumptions about the near-wal variation of velocity.
- Examples include k-ω SST or low Reynolds number k-ε (a variant of k-ε that does not use wall functions and integrates right through the boundary layer).
- When "Low-Re number wall treatment" is used, following requirements on mesh are imposed:
- The (y^{+}) should be of the order of 1.
- Under no circumstances, the y^{+} should increase above 5.
- The number of elements in the boundary layer should be ~ 10. This ensures that the mesh corresponding to y^{+} falls in the log-law region.
- For all "high-Re number wall treatments" such as (a) scalable wall function [exploits a limiter in the (y^{+}) calculations such that y^{+} = min [y^{+}, 11.36] in CFX, (b) standard wall functions in FLUENT and (c) standard wall function in STAR-CCM+.
- The centroid of cell adjacent to walls should lie in logarithmic region: y^{+} >= 30 in FLUENT and STAR CCM+ because they use cell-centred schemes (staggered arrangement).
- The nodes of 1st layer of elements should fall in y^{+} >= 30 in ANSYS CFX. This software uses vertex-centred scheme (collocated arrangement).
- A scalable wall function is another way of saying wall surface coincides with the edge of the viscous sublayer if the first grid point is too close to the wall. Still other explanation is that scalable wall functions effectively shifts the near-wall mesh point to y^{+} =11 regardless of the close it is to the wall. Thus, scalable wall functions can be assumed less sensitive to mesh as compared to standard or non-equilibrium wall functions but that does not imply they are necessarily more accurate.
- What happens if more than 1 layers of mesh are inside viscous layer - how does scalable wall function treats the second node falling inside viscous layer?
- No noticeable impact on result is likely even if the y^{+} falls up to the level of 11.36, in any "high-Re number wall treatments"
- It is not possible to maintain y^{+} value in a very narrow band. It is advised to target a y^{+} band of 30 ~ 50 for higher accuracy. A range of 30 ~ 100 on y^{+} yields result which is acceptable for design iterations.
- The y^{+} value should not be allowed to cross the upper limit of 300.
- It is not possible to maintain this limit (y^{+} > 11.36) in few areas of the fluid flow such as separation and stagnation (recirculation) zones. This is numerically acceptable and no impact of results should be anticipated because boundary layers are neither orderly nor as steep as normal wall-bounded turbulent flows with separations.
- There are "numerical improvisations" where the wall treatment toggles automatically between "Low-Re number wall treatment" and "high-Re number wall treatments". This is known as "automatic near-wall treatment" where wall treatment automatically switches from wall-functions to a low-Re formulation as the y^{+} value changes based on local mesh refinement.
- Enhancced Wall Treatment or Enhanced Wall Functions (EWF): Blended law-of-the-wall are applied for momentum, energy, species and turbulent dissipation frequency. The blending function for momentum is of the form u^{+} = e^{Γ} . u^{+}_{LAM} + 1/e^{Γ} . u^{+}_{TURB}.
- Enhanced Wall Treatment (EWT) and Scalable Wall Functions are expected to be insensitive to y+ among near wall modeling approaches.
- In other words, as the grid is refined or coarsed in certain range of y^{+}, the results will not show strong dependence on grid resolution.
- EWT is also a "viscous sub-layer resolving approach". Thus, when it is used with a mesh that is fine enough to capture the flow behavior in the viscous sublayer, the accuracy of results can be be expected to be better than with scalable wall functions.
- In a case where the nodes of the first layer of elements near the wall are always in the log-layer, no siginificant difference between EWT and scalable wall functions would be observed.
- Realizable k-ε Model: Stanadard k-ε model cannot universally ensure positive value of normal stresses and Schwarz’s inequality of shear stresses. Here equation for k remains same as standard k-ε model but the ε-equation is based on a transport equation for the mean-square vorticity fluctuation and C_{μ} is a function of mean velocity field and turbulence. It addressed the anamoly of standard k-ε model in round jet (jet emanating from a circular hole) and plane jet (jet emanating from a rectangaulr opening of high aspect ratio) and predicts the spreading rate correctly for both types of jets.
- Law of the Wall:
- Friction Velocity: U_{τ} = [τ_{W}/ρ]^{0.5}, τ_{W} = wall shear stress, ρ = density of fluid
- Dimensionless Velocity: U^{+} = U_{τ} / U_{Y1}, U_{Y1} = Velocity at first node near wall and parallel to the wall
- Logarithmic law of the wall valid in the range 30 ≤ y+ ≤ 130: U^{+} = ln[E . y^{+} / κ], E = surface roughness parameter, equal to 8.6 for smooth walls and κ is von Karman constant = 0.41
The term "Isotropic" refers to both Isotropic diffusion and Isotropic dissipation whereas in reality near-wall flows are anisotropic due to presence of walls. Imagine the flow perpendicular to the wall and near the wall - the fluctations in direction perpendicular to the wall will slow-down faster in the direction parallel to it. In the following diagram, one-, two- and four-equation models refer to the number of additional partial differential equations required (beyond the equations for mean quantities u, v, w, p, T) to complete a RANS model.
Some characteristics of turbulence models
k-ε model
- Advantages:
- Good results for many industrial applications dealing with wall-bounded flows with low pressure gradients.
- Stable and numerically robust
- Computationally less expensive, no sensitivity of free-steam conditions
- Realizable version of k-ε model is good for complex flows with large strain rates, recirculation, rotation, separation, strong pressure gradients.
- Disadvantages:
- Limited ability to predict secondary flow characteristics such separation and reattachment (such as flow over aerofoil with non-zero angle of attack, cross flow over a cylinder).
- Not accurate in flows with high streamline curvature and sudden changes in the mean strain rate (such as back-facing step), strong swirling flow (such as cyclones, hydrocyclones and stirred tanks).
- Like any eddy-viscosity model, fails to predict the cases where turbulence transport or non-equilibrium effects are important in flows. It is based on isotropic turbulence.
Spalart-Allmaras model
- Key characteristics:
- In standard form, it is a Low-Re number model and hence no wall function is used. That further imposes a restrictions to have mesh fine enough so that y^{+} is < 1.0 everywhere.
- Since this is a Low-Re number model, it can be used with "Automatic Wall Treatment" or "All y^{+} treatement" methods of turbulence modeling.
- This models is a good choice for applications in which the boundary layers are largely attached and separation is not present or mild separation is expected. Typical examples would be flow over a wing, ship-hulls, missiles, fuselage or other aerospace external-flow applications.
- The Spalart-Allmaras model for RANS equations is not recommended for flows dominated by free-shear layers (such as jets), flows where complex recirculation occurs (especially with heat transfer) and natural convection.
- The Spalart-Allmaras models with DES can be used for flows dominated by free-shear layers (such as jets), flows where complex recirculation occurs (especially with heat transfer) and natural convection.
Reynolds Stress Model (RSM) or Second Moment Closure Methods
- Key characteristics:
- Reynolds stresses are not modeled as Boussinesq Hypothesis. However, modeling is still required for many terms in the transport equations.
- Recommended for complex 3-D turbulent flows with large streamline curvature and swirl, but the model is computationally intensive, difficult to converge than eddy-viscosity models such as k-ε or Spalart-Allmaras models.
- Anisotropy of turbulence is accounted for, quadratic pressure-strain option improves performance for many basic shear flows
- Most suitable for curved ducts say U or S-bends, rotating flow passages, combustors with large inlet swirl and cyclone separators.
Standard k-ω Model (SKO)
- Key characteristics:
- Specific dissipation rate ω = k/ε solved instead of ε
- Demonstrates better performance for wall bounded and low-Re flows and potential to handle transitional flows (though tend to predict a transition early).
- Suitable for complex boundary layer flows with adverse pressure gradient and separation (external aerodynamics and turbomachinery).
- Separation is typically predicted to be higher and earlier than experimentally observed values.
Shear Stress Transport (SST) k-ω Model
- Key characteristics:
- Specific dissipation rate ω = k/ε solved in inner layer (log-layer and viscous sub-layer) and transitions to a k-ε model away from the wall (but not same as standard k-ε equations).
- The boundary conditions for SST model are the same as the k-ω model and is relatively less sensitive to the free stream value of ω.
- Suitable for complex boundary layer flows with adverse pressure gradient and separation (external aerodynamics and turbomachinery).
References
- User manuals of commercial programs like FLUENT, CFX, STAR-CCM+
- Biswas Gautam, Som S. K., Chakraborty Suman: Introduction to Fluid Mechanics and Fluid Machines, Mc Graw Hill Education, Third Edition, 2012
- J. H. Ferziger and M. Peric’: Computational Methods for Fluid Dynamics, 4th Edition, Springer, Berlin, 2002