Turbulence Modeling in

Computational Fluid Dynamics (CFD)

Zhaoyuan Wang, Maysam
Mousaviraad, Jun Shao, Shanti
Bhushan, Tao Xing and Fred Stern

IIHR—Hydroscience & Engineering

C. Maxwell Stanley Hydraulics Laboratory
The University of Iowa

58:160 Intermediate Mechanics of Fluids

http://css.engineering.uiowa.edu/~me_160/
Oct. 9, 2020
 Outline
1. Characteristics of turbulence
2. Approaches to predicting turbulent flows
3. Reynolds averaging
4. RANS equations and unknowns
5. The Reynolds-Stress Equations
6. The closure problem of turbulence
7. Characteristics of wall-bound turbulent flows
8. Turbulence models and ranges of applicability
7.1. RANS
7.2. LES/DES
7.3. DNS
9. Example: diffuser
2
 Characteristics of turbulence
• Randomness and fluctuation: u = U + u

• Nonlinearity: Reynolds stresses from the

nonlinear convective terms
• Diffusion: enhanced diffusion of momentum,
energy etc.
• Vorticity/eddies/energy cascade: vortex
stretching
3
 Characteristics of turbulence

•Dissipation: occurs at smallest scales

•Three-dimensional: fluctuations are always 3D
•Coherent structures: responsible for a large part of
the mixing
•A broad range of length and time scales: making
DNS very difficult
4
Approaches to predicting turbulent flows
• AFD, EFD and CFD:
– AFD: No analytical solutions exist
– EFD: Expensive, time-consuming, and sometimes
impossible (e.g. fluctuating pressure within a flow)
– CFD: Promising, the need for turbulence modeling
• Another classification scheme for the approaches
– The use of correlations: CD = f ( Re )
– Integral equations: reduce PDE to ODE for simple cases
– One-point closure: RANS equations + turbulent models
– Two-point closure: rarely used, FFT of Two-point
equations
– LES: solve for large eddies while model small eddies
– DNS: solve NS equations directly without any model
5
Deeper insights on RANS/URANS/LES

RANS

URANS

LES

6
 Reynolds Averaging
• Time averaging: for stationary
turbulence
• Spatial averaging: for homogenous
turbulence
• Ensemble averaging: for any turbulence
• Phase averaging: for turbulence with
periodic motion

7
 RANS equations and unknowns
• RANS equation
U i U i P 

t
+ U j
x j
=− +
xi x j
(
2  S ji −  u j ui )
• RANS equation in conservative form
U i  P 

t
+
x j
(  )
U jU i + u j ui = − +
xi x j
( 2 S ji )

• Numbers of unknowns and equations

– Unknowns: 10 = P (1) + U (3) + ( −uiuj ) (6)
– Equations: 4 = Continuity (1) + Momentum (3)

8
 The Reynolds-Stress Equation
•Derivation: Taking moments of the NS equation.
Multiply the NS equation by a fluctuating property and time
average the product. Using this procedure, one can derive a
differential equation for the Reynolds-stress tensor.
 ij  ij U j U i ui u j
+Uk = − ik −  jk + 2
t xk xk xk xk xk
 ui p uj p     ij 
+ + +  + uiuj uk 
  x j  xi  xk  xk 
 
•NEW Equations: 6 = 6 equation for the Reynolds stress tensor
•NEW Unknowns: 22 = 6 + 6 + 10
 ui p uj p 
 +  → 6 unkowns
  x j  xi 
 
ui u j
2 → 6 unkowns uiuj uk → 10 unkowns
xk xk 9
 The closure problem of turbulence

• Because of the non-linearity of the Navier-Stokes

equation, as we take higher and higher moments,
we generate additional unknowns at each level.
• In essence, Reynolds averaging is a brutal
simplification that loses much of the information
contained in the Navier-Stokes equation.
• The function of turbulence modeling is to devise
approximations for the unknown correlations in
terms of flow properties that are known so that a
sufficient number of equations exist.
• In making such approximations, we close the
system.

10
 Characteristics of Wall-Bound Turbulent Flows
• The turbulent boundary layer (zero-
pressure gradient) has universal
velocity distribution near the wall
(inner-layer) (Clauser 1951)
u yu
= u+ =  = y+ : y+  5
u 
1
u+ = ln( y + ) + 5.1 : 104  y +  30
0.4
k
2
= 0.3 : 10 4
 y +
 30
u

• The velocity defect when plotted

vs. y/ collapse on a single curve
(outer-layer)
u − u y
= 9.6(1 − ) 2
u 
11
 Turbulent models (RANS)
• Boussinesq eddy-viscosity approximation
• Algebraic (zero-equation) models
– Mixing length
– Cebeci-Smith Model
– Baldwin-Lomax Model
• One-equation models
– Baldwin-Barth model
– Spalart-Allmaras model
• Two-equation models
– k- model
– k- model
• Four-equation (v2f) models
• Reynolds-stress (seven-equation) models
12
 Turbulent models (RANS)
• Boussinesq eddy-viscosity approximation

 U i U j  2
−uiuj =  T  +  −  ij k
 x xi
 j  3

• Dimensional analysis shows: T = C qL , where q is a

turbulence velocity scale and L is a turbulence length
1
scale. Usually q = 2k where k = uiui is the turbulent
2
kinetic energy. Models that do not provide a length
scale are called incomplete.

13
 Turbulent models (0-eqn RANS)
• Mixing length model:  T = lmix
2 dU
dy

• Assume lmix =  ( x ) for free shear flow, then

─ =0.180 for far wake
─ =0.071 for mixing layer
─ =0.098 for plane jet
─ =0.080 for round jet
• Comments:
─Reliable only for free shear flows with different  values

─Not applicable to wall-bounded flows

14
 Turbulent models (0-eqn RANS)
Cebeci-Smith Model (Two-layer model)
 Ti , y  ym

T =  Where ym is the smallest value of y for  T =  T
 To y  ym

i o

12
 U 2  V 2 
=  y 1 − e− y 
+
A0+
Inner layer:  Ti = lmix 
2
 +   lmix
 y   x    
Outer layer:  T = Uev* FKleb ( y; )
0
−1 2
 dP dx 
Closure coefficients:  = 0.40  = 0.0168 A+ = 26 1 + y
−1  U2 
  y 
6

FKleb ( y;  ) = 1 + 5.5   

 v* =  (1 − U U e ) dy
Comments:      0

Three key modifications to mixing length model;

Applicable to wall-bounded flows;
Applicable to 2D flows only;
Not reliable for separated flows;
 ,  v* ,U e difficult to determine in some cases;
15
 Turbulent models (0-eqn RANS)
Baldwin-Lomax Model (Two-layer model)
 Ti , y  ym

T =  Where ym is the smallest value of y for  T =  T
 To y  ym

i o

lmix =  y 1 − e− y 
+
A0+
Inner layer:  T = lmix
2
 i  
Outer layer: T = Ccp Fwake FKleb ( y; ymax CKleb )
0

Closure coefficients:  = 0.40  = 0.0168 A+ = 26 Ccp = 1.6 CKleb = 0.3 Cwk = 1

−1 12
  y  
6
 V U 2  W V  2  U W 2 
FKleb ( y; ymax CKleb ) = 1 + 5.5     =  −  + −  + −  
  ymax CKleb    x y   y z   z x  
Comments:
Applicable to 3D flows;
Not reliable for separated flows;
No need to determine  ,  v , U e
*

16
 Turbulent models (1-eqn RANS)
Baldwin-Barth model
Kinematic eddy viscosity :  T = C RT D1D2
Turbulence Reynolds number :
(
  RT ) + U  ( R ) = ( C (
 2  RT )− 1  T   RT ( )
 2 f 2 − C 1 )  RT P + ( +  T   )
T

t
j
xj xk xk   xk xk

Closure coefficients and auxiliary relations:

+ C
C 1 = 1.2 C 2 = 2.0 C = 0.09 A0+ = 26 A2 = 10  = 0.41
1
= ( C 1 − C 2 )
 2

 U U j  U 2 U U  − y + A2+
P =  T  i
+ 
i
− k k
 D1 = 1 − e
− y+ A0+ D2 = 1 − e
 x xi  x j 3 xk xk 
 j
C 1  C 1   1   1  D2 − y+ A0+ D +
A2+ 
f2 = + 1 −  + + D D 
1 2  D D +  +e + +1 e− y 
C 2  C 2    y
1 2
  D1D2  A0 A2  

17
 Turbulent models (1-eqn RANS)
Spalart-Allmaras model
Kinematic eddy viscosity :  T =  fv1
Eddy viscosity equation:
    1     cb1  
2

+U j = cb1S − cw1 f w   +  ( +  ) +
t x j  d   xk  xk   xk xk

Closure coefficients and auxiliary relations:

cb1 (1 + cb 2 )
cb1 = 0.1355 cb 2 = 0.622 cv1 = 7.1  = 2 3 cw1 = +
2 
 
cw2 = 0.3 cw2 = 2.0  = 0.41  = g = r + cw 2 ( r 6 − r ) r=
 S 2 d 2

3  16
 1 + cw6 3  
f v1 = fv 2 = 1 − fw = g  6 6  S=S+ fv2 S = 2ij ij
 +c
3 3
1 +  f v1  2d 2
v1  g + cw3 

18
 Turbulent models (1-eqn RANS)
Comments on one-equation models:
1.One-equation models based on turbulence kinetic energy are incomplete as
they relate the turbulence length scales to some typical flow dimension. They
are rarely used.
2.One-equation models based on an equation for the eddy viscosity are
complete such as Baldwin-Barth model and Spalart-Allmaras model.
3.They circumvent the need to specify a dissipation length by expressing the
decay, or dissipation, of the eddy viscosity in terms of spatial gradients.
4 Spalart-Allmaras model can predicts better results than Baldwin-Barth
model, and much better results for separated flow than Baldwin-Barth model
and algebraic models.
5 Also most of DES simulations are based on the Spalart-Allmaras model.

19
 Turbulent models (2-eqn RANS)
k- model:  T = C k 2 
k k U i   k 
+U j =  ij − + ( + T  k ) 
t x j x j x j  x j 

   U i 2    
+U j = C 1  ij − C 2 + ( + T   ) 
t x j k x j k x j  x j 
,
C 1 = 1.44 C 2 = 1.92 C = 0.09  k = 1.0   = 1.3  =  ( C k ) l = C k 3 2 

T = k 
,
k- model:
k k U i   k 
=  ij −  *k + ( +   T )
,
+U j *

t x j x j x j  x j 

,
   U i    
+U j =   ij −  2 + ( +  T ) 
t x j k x j x j  x j 

=
13
 = 0 f   * =  0* f  =
1 1
0 =
9 1 + 70 
* * = f =
25 2 2 125 1 + 80 
1, k  0
ij  jk S ki  1 k 
 = 0* =
9 f* = 1 + 680  k2 k =
( ) 1 + 400  2  k  0  3 x j x j
* 3
0 100
 k

 =  *k l = k1 2 
20
 Turbulent models (2-eqn RANS)
Comments on two-equation models:
1. Two-equation models are complete;
2. k- and k- models are the most widely used two-
equation models and a lot of versions exist. For example, a
,
popular variant of k- model introduced by Menter has
been used in our research code CFDSHIP-IOWA. There are
also a lot of low-Reynolds-number versions with different
damping functions.
3. k- model shows better results than k- model for
flows with adverse pressure gradient and separated flows
as well as better numerical stability.

21
 Turbulent models (4-eqn RANS)
v2f-k model:  T = C v 2T
Dv 2 v2   v 2  C1  v 2 2  P
= kf −  + ( + t )  L  f − f =  −  − C2 k
2 2

Dt k x j  x j  T  k 3  k

 D  2    t   
Dk   t  k  = C 1 P − C 2 +   +  
= P − +  +   k x j    x j 
x j    Dt k 
, Dt k  x j 

C1 = 0.4 C2 = 0.3 CL = 0.3 C = 70 C 2 = 1.9   = 1.0
,  k 3 2   3 14
 
4
   L = max C
C 1 = 1.3 + 0.25 1 + ( CL d 2 L )  k  
2
  =  
, C 
T max

 
, 6
 
 
L
    
, 1− n
v2f-k model:  T = C k n v2 
,
Dv 2    v 2  1  v2 2  1 v2 P
= kf − 6v 2 + ( + t )  L  − f = ( C1 − 1)  −  − 5
2 2
− C2 k
Dt k x j  x j  T  k 3 T k k
 
1− n
Dk U i    t  k  D  U i  2    t   
=  ij −  * k + 2 v
 +   =   ij −    +   +  
Dt x j x j   k  x j  Dt k x j k  x j     x j 
 

 = 0.45977  * = 0.09  = 3 40  k = 1.0   = 1.5

    k 3 2   3  
14
k L = max CL
T = max  , 6  , C   

   
      
22
 Turbulent models (4-eqn RANS)
Comments on four-equation models:
1. For two-equation models a major problem is that it is hard to
specify the proper conditions to be applied near walls.
2. Durbin suggested that the problem is that the Reynolds number is
low near a wall and that the impermeability condition (zero normal
, velocity) is far more important. That is the motivation for the equation
for the normal velocity fluctuation.
3. It was found that the model also required a damping function f,
hence the name v2f model.
4. They appear to give improved results at essentially the same cost
as the k- and k- models especially for separated flows. Hopefully the
v2f-k model can have better numerical stability than v2f-k model as
their counterparts behave in two-equation models.

23
 Turbulent models (7-eqn RANS)
Some of the most noteworthy types of applications for which models based
on the Boussinesq approximation fail are:
1. Flows with sudden changes in mean strain rate
2. Flow over curved surfaces
3. Flow in ducts with secondary motions
4. Flow in rotating fluids
5. Three-dimensional flows
6. Flows with boundary-layer separation

In Reynolds-stress models, the equations for the Reynolds stress

tensor are modeled and solved along with the -equation:

 ij  ij U j U i ui uj  ui p uj p     ij 

+ Uk = − ik −  jk + 2 + + +  + u  
i j k 
u u 
t xk xk xk xk xk   x j  xi  xk  xk 

24
 Turbulent models (7-eqn RANS)
These are some versions of Reynolds stress models:
1. LRR rapid pressure-strain model
2. Lumley pressure-strain model
3. SSG pressure-strain model
4. Wilcox stress- model

Comments on Reynolds stress models:

1.Reynolds stress models require the solution of seven additional
PDEs and those equation are even harder to solve than the two-equation
models.
2. Although Reynolds stress models have greater potential to
represent turbulent flow more correctly, their success so far has been
moderate.
3. There is a lot of current research in this field, and new models
are often proposed. Which model is best for which kind of flow is not clear
due to the fact that in many attempts to answer this question numerical
errors were too large to allow clear conclusions to be reached.

25
 Turbulent models (LES)

• Large scale motions are generally much more energetic than

the small scale ones.
• The size and strength of large scale motions make them to be
the most effective transporters of the conserved properties.
• LES treats the large eddies more exactly than the small ones
may make sense
• LES is 3D, time dependent and expensive but much less costly
than DNS.

26
 Turbulent models (LES, filtering)
• LES needs a velocity field that contains only the large scale
components of the total field, which is achieved by filtering
the velocity field (Leonard, 1974)

ui (x ) =  G (x −  )ui ( )d
• G(x-) is the filter kernel, is a localized function, which
includes a Gaussian, a box filter (a simple local average) and a
cutoff (a filter which eliminates all Fourier coefficients
belonging to wavenumbers above a cutoff)
• Each filter has a length scale associated with it, .
• Eddies of size large than  are large eddies while those
smaller than  are small eddies and need to be modeled.
1   − x − 2  sin (kc (x −  ))
 if x −     
12
G(x −  ) =
G(x − ) =  2 G ( x −  ) =  2  exp  
kc (x −  )
     2

0 otherwise   
kc =

Box or top-hat filter Gaussian filter Cutoff filter
27
 Turbulent models (LES,
Governing Equations)
• Filtered Navier-Stokes equations (constant density,
incompressible):
 ( u i )
=0
xi
(ui ) (ui u j ) p    u u j 
+ =− +   i + 
t x j xi x j   x j xi 
• Note that: ui u j  ui u j
Introducing Subgrid-scale  ijS = − (uiu j − uiu j )
Reynolds Stress

• The filter width > grid size h

• The models used to approximate the SGS Reynolds
stress are called subgrid-scale (SGS) or subfilter-scale
models.

28
Turbulent models (LES, Smagorinsky model)
• The earliest and most commonly used subgrid scale
model is one proposed by Smagorinsky (1963), which is
an eddy viscosity model.
• As the increased transport and dissipation are due to the
viscosity in laminar flow, it seems reasonable to assume
that
1 S  ui u j 
 −  kk  i , j = t 
S
+  = 2t Sij

 
i, j
3  x j xi 
Sij Strain rate of the large scale or resolved field
 t Eddy viscosity t = CS  S
2 2

CS  0.2 Model constant

• Drawbacks: 1. Cs is not constant
2. Changes of Cs are required in all shear flows
3. Need to be reduced near the wall.
4. Not accurate for complex and/or higher
Reynolds number flows.
29
Turbulent models (LES, Scale-similarity model
and Dynamic model)
• Dynamic model:
1. filtered LES solutions can be filtered again using a filter
broader than the previous filter to obtain a very large scale
field.
2. An effective subgrid-scale field can be obtained by
subtraction of the two fields.
3. model parameter can then be computed.
Advantages: 1. model parameter computed at every spatial
grid point and every time step from the results of LES
2. Self-consistent subgrid-scale model
3. Automatically change the parameter near the wall and in
shear flows
Disadvantages: backscatter (eddy viscosity<0) may cause
instability.

30
 Turbulent models (DES)
• Massively separated flows at high Re usually involve
both large and small scale vortical structures and very
thin turbulent boundary layer near the wall
• RANS approaches are efficient inside the boundary
layer but predict very excessive diffusion in the separated
regions
• LES is accurate in the separated regions but is unaffordable
for resolving thin near-wall turbulent boundary layers at
industrial Reynolds numbers
• Motivation for DES: combination of LES and RANS. RANS
inside attached boundary layer and LES in the separated
regions

31
 DES Formulation
• Modification to RANS models was straightforward by
substituting the length scale d w , which is the distance to the
closest wall, with the new DES length scale, l defined as:
l = min(dw , CDES )  = max( x ,  y ,  z )
• where C is the DES constant,  is the grid spacing and is
DES

based on the largest dimensions of the local grid cell, and  x ,  y ,  z

are the grid spacing in x, y and z coordinates respectively
• Applying the above modification will result in S-A Based,
standard k- (or k-) based and Menter’s SST based DES
models, etc.

=  *k =  k 3 2 / lk − lk − = k 1 2 (  * )
k
DRANS
32
D k
DES = k l = min(lk − , CDES )
l
32
 Resolved/Modeled/Total Reynolds
stress (DES)
TKE
Modeled

Resolved
Total

Modeled Reynolds stress: u 'v ' = − t  u + v 

 y x 
 33
 Resolved/Modeled/Total Reynolds
stress (URANS)
Modeled Total

Resolved

34
Deep insight into RANS/DES (Point oscillation on
free surface with power spectral analysis)
EFD DES URANS

-5/3

2HZ 35
 Direct Numerical Simulation (DNS)
• DNS is to solve the Navier-Stokes equation directly
without averaging or approximation other than numerical
discretizations whose errors can be estimated (V&V) and
controlled.
• The domain of DNS must be at least as large as the
physical domain or the largest turbulent eddy (scale L)
• The size of the grid must be no larger than a viscously
determined scale, Kolmogoroff scale, 
• The number of grid points in each direction must be at
least L/ 
• The computational cost is proportional to
Re 3L/ 4  (0.01Re )
34

• Provide detailed information on flow field

• Due to the computational cost, DNS is more likely to be a
research tool, not a design tool.

36
 Two-Phase Turbulence
• Most of the turbulence models (RANS and LES) are developed based on single-
phase flows. The effect of turbulence on the interface and the interface induced
turbulence are not considered.
• The eddy viscosity is found to be over-predicted in the gas area near the interface
(Wang et. al, 2010).
• A buoyancy term Gb is added to the TKE
Equation (Devolder et al., 2017/2018):

Eddy viscosity for KCS wave breaking. Left: V4.5; right: V5.5.

• This additional term suppresses the turbulence level at the air-water interface

37
 Examples (Diffuser)
• Asymmetric diffuser with
separation is a good test case
for turbulence models.
• A inlet channel was added at
the diffuser inlet to generate
fully developed velocity profile
• Boundary layer in the lower
diffuser wall will separate due
to the adverse pressure
gradient.
• Results shown next include
comparisons between V2f and
k-ε
• LES simulation of this geometry
can be found in:
M. Fatica, H. J. Kaltenbach, and
R. Mittal, “Validation of LES in a
Plain Asymmetric Diffuser”,
center for turbulence research,
annual research briefs, 1997

38
 Examples (Diffuser)
• Mean velocity
predicted by V2f
agreed very well
with EFD data,
particular the
separation
region is
captured.
• K-ε model fails
to predict the
separation
caused by
adverse
pressure
gradient.

39
 Examples (Diffuser)
v2f

k-ε

• TKE predicted by V2f agreed better with EFD

data than k- ε model, particular the
asymmetric distribution.
• Right column is for the skin friction coefficient 2x/H
on the lower wall, from which the separation
and reattachment point can be found. 40
 References

1. J. H. Ferziger, M. Peric, “Computational Methods for Fluid

Dynamics,” 3rd edition, Springer, 2002.
2. D. C. Wilcox, “Turbulence Modeling for CFD”, 1998.
3. S.B.Pope, “Turbulent Flows,” Cambridge, 2000
4. P. A. Durbin and B. A. Pettersson Reif, “Statistical Theory and
Modeling for Turbulent Flows,” John Wiley & Sons, LTD, 2001.
5. P.A. Davidson, “Turbulence: An Introduction for Scientists
and Engineers,” Oxford, 2004.

41