









Study with the several resources on Docsity
Earn points by helping other students or get them with a premium plan
Prepare for your exams
Study with the several resources on Docsity
Earn points to download
Earn points by helping other students or get them with a premium plan
Although there is a strong confidence in our ability to solve the Navier Stokes equations on modern computers, that confidence is limited to laminar – not turbulent – flows. Unfortunately, most flows of engineering interest are turbulent.
Typology: Study notes
1 / 16
This page cannot be seen from the preview
Don't miss anything!










Jacaranda (Engineering) 3333 Mail Code Phone: 818.677.
Email: [email protected] 8348 Fax: 818.677.
College of Engineering and Computer Science Mechanical Engineering Department Computational Fluid Dynamics Notes
Larry Caretto January 30, 2010
Introduction
Although there is a strong confidence in our ability to solve the Navier Stokes equations on
modern computers, that confidence is limited to laminar – not turbulent – flows. Unfortunately,
most flows of engineering interest are turbulent. Although there have been some solutions of the
Navier Stokes equations directly for turbulent flows (a process known as DNS for direct numerical
simulation), these computations are not practical for engineering analysis. They require
extensive computer resources, are applicable only to simple geometries, and do not directly
provide the average quantities of interest in engineering design.
Instead of direct simulation, CFD applications to turbulent flow use models. These models range
from simple algebraic models to those which require the solution of one or more partial differential
equations. In general, the more complex the flow, the more complex of a turbulence model is
required. Particular problems occur in modeling turbulence with large amounts of swirl as occurs
in combustors and turbomachinery.
These notes begin with a simple discussion of turbulence and move on to the formal statistical
analyses that are used in turbulence. We will see that the formal statistical analysis leads to a
closure problem that requires additional assumptions to model physical processes.
The most important idea in these notes is that the choice of wall boundary conditions that you use
for turbulent flows determines limits on the grid spacing. The choice to (a) model the laminar
sublayer at a wall or to (b) use wall functions sets conditions on the grid size near the wall that
must be followed to ensure a correct solution.
What is turbulence?
The transition between laminar and turbulent flow, in any geometry, is characterized by the
Reynolds number, Re, based on the appropriate characteristic velocity, U, and the appropriate
length parameter, L, for the flow.
Re = [7-1]
As usual, ν is the kinematic velocity. When the Reynolds number is below a critical level, the flow
is laminar. Increases in the Reynolds number above this critical level lead to a transition region,
and further increases lead to a fully turbulent region.
For free convection flows, the borders between laminar, transition and turbulent flow are
determined by the Rayleigh number,
3 g TL Ra
= , where T (^) P
is the coefficient of
thermal expansion, which equals 1/T for an ideal gas, g is the acceleration of gravity, ΔT is a
characteristic temperature difference, and α is the thermal diffusivity.
Page 7.2 L. S. Caretto, January 30, 2010 Turbulence Modeling
Turbulence is characterized by random behavior. The value of the flow properties (velocity
components, temperature, species concentrations) at any point in the flow vary in a random
manner over time. This random variation may be about a constant (steady-state) average, or it
may be about some overall average trend in time.
Turbulence is a three-dimensional phenomenon. A fully developed laminar flow in a cylindrical
pipe varies only in the radial direction. (There is no variation in the angular coordinate or down
the length of the pipe after the flow becomes fully developed.) However, if such a flow were
turbulent, there would be variation in the flow in all three coordinate directions.
Turbulent flows are characterized by turbulent eddies. In
pictures of turbulent flows we see underlying structures that
move (translate) and rotate in the mean flow. These
underlying structures are called turbulent eddies. When we
examine a visualization of a turbulent flow we see that these
eddies have many different length scales.
The largest eddies
have a characteristic velocity and length scale that are of the
same order of magnitude as the characteristic velocity and
length scale (U and L) for the main flow. In these large eddies
(with large Reynolds numbers) the inertia effects dominate and
viscous effects are negligible. [Recall that the Reynolds
number is a measure of the ratio of inertia forces (ρU
2 ) to
viscous forces (μU/L).]
The main flow (at some mean velocity like U) supplies energy
to the large eddies. This increases their rotation rate and
decreases their size. The larger eddies supply smaller eddies with energy. This creates the
following energy cascade.
In this way, the turbulence produces increased energy dissipation (and hence increased pressure
loss) for the flow. The turbulence structures provide a mechanism by which energy is
transformed from the main flow kinetic energy into viscous dissipation.
** This mechanism is not
present in laminar flows. Although turbulence increases pressure drop because of this energy
transfer to viscous forces, it also produces increased mixing rates that are important in heat
transfer, mixing processes, and combustion.
Although smoking is not permitted in campus buildings, all of us have probably noticed the
smoke from a burning cigarette, which provides a good visualization of the development and
structure of turbulent flows. This flow usually starts as a laminar flow next to the cigarette, then
goes through a transition region to a turbulent flow. In the turbulent flow there are many eddies
with varying characteristic sizes. These eddy structures persist for a while then break up. This is
illustrated in the photograph shown here taken from the web site shown below:
Humphrey_Bogart_by_Karsh_(Library_and_Archives_Canada).jpg
** Lewis F. Richardson (1881-1953), who wrote an Article entitled “Numerical Prediction of
Weather Process: in 1922, wrote a poem that characterizes this process.
Big whorls have little whorls That feed on their velocity And little whorls have lesser whorls And so on to viscosity.
Laminar
Turbulent
Page 7.4 L. S. Caretto, January 30, 2010 Turbulence Modeling
The overall result of equation [7-5] is usually summarized as follows.
Here we see that the average of the product consists of two terms. The first is the product of the
averages. The second is the average of the product of the fluctuations. Only if these fluctuations
were perfectly anticorrelated would this average of the product be zero. It is this basic result that
provides an introduction to turbulence quantities in the time-averaged Navier-Stokes equations.
Equation [7-6] tells us that the average of the product of two fluctuating quantities, which may be
the same, will typically not vanish in a turbulent flow. We can use this result to define the root
mean square average as our typical measure of the magnitude of a turbulent fluctuation. We
define this quantity as follows.
Δ
t
rms dt t (^) 0
2 2 ( ')
ϕ (ϕ') ϕ [7-7]
The velocity fluctuations in a turbulent flow, using Cartesian coordinates and the usual velocity
components, are u’, v’, and w’. The kinetic energy of the turbulence, k, is given by the sum of the
squares of these velocity fluctuations.
2 2 2 ( ') ( ') ( ') 2
k = u + v + w [7-8]
We see that this is a typical fluid mechanics “kinetic energy” term, which is actually a kinetic
energy per unit mass.
The time average of the spatial derivative of a flow quantity, φ, is found as follows.
i
t
i
t
i i x
dt x t
dt x t x ∂
Δ Δ
0 0
[7-9]
Thus the average of the derivative (with respect to the coordinate x (^) i ) is simply the derivative of the
average. This result applies to derivatives of any order.
We are now ready to derive the time-average Navier-Stokes equations. We start with equation
[1-81]; this general transport equation is copied below.
ρφ ρ ϕ (ϕ ) ϕ ( ϕ ) S x x x
u
t
c i i i
i
∂
[1-81]
To simplify the derivation, we will assume that we have constant properties and a zero source
term. Further assume that we have a steady problem so that the time derivative is zero. We will
divide the equation by the product, ρc, and write our initial equation as follows.
i i i
i x x x
u
∂
[7-10]
Turbulence Modeling L. S. Caretto, January 30, 2010 Page 7.
In this equation we have defined γ
(φ) as follows
ϕ ϕ
( ) ( ) Γ = [7-11]
Recall that in the general equation c = 1 unless φ is the temperature in the energy equation. In
most cases, then, the denominator in equation [7-11] is simply the density, ρ. If φ represents the
temperature in the energy equation then c = c (^) p or c (^) v, depending on the equation. In the
momentum equation, γ
(φ) would be the kinematic viscosity, ν = μ/ρ; in the energy equation with
temperature, it would be the thermal diffusivity, α = k/ρc (^) p. (Note that for the momentum equation,
Γ
(φ) is the dynamic viscosity, μ, and for the energy equation, with temperature as the variable, Γ
(φ)
is the thermal conductivity.)
We want to obtain the time average of the entire transport equation in [7-10]. We do this by
substituting both sides of this equation into the definition of the time average from equation [7-3].
Δ Δ
⎥ ⎦
t
i i
t
i
i dt t x x
dt x
u
t (^) 0
( )
0
[7-12]
Both sides of equation [7-12] are averages of spatial derivatives. We can use the result of
equation [7-9] that the average of (spatial) derivatives are simply the derivatives of averages to
rewrite equation [7-12] as follows.
i i i
i x x x
u
∂
[7-13]
We can use the result of equation [7-6] for the average of a product to rewrite equation [7-13] as
follows.
i i i
i i x x x
u u
∂
[7-14]
If we compare equations [7-10] and [7-14], we see that the time averaged equation in [7-14] has
two differences. The first is that the average quantities appear in [7-14] in place of the
instantaneous quantities in [7-10]. The second change is the addition of the correlation term
involving the average of two fluctuating quantities. This is the mathematical term that leads to
problems in the computations of turbulent flows using the time-averaged Navier Stokes
equations.
When we use the general balance equation to represent the continuity equation, we take φ = 1
and Γ = 0. In this case the fluctuation term vanishes and we have the constant-property, steady-
state equation for turbulent flows as follows.
i
i
x
u [7-15]
In principle, we could take the time average of equation [7-14] to get a differential equation for the
correlation terms. However, this would introduce new, higher order correlations. This creates
Turbulence Modeling L. S. Caretto, January 30, 2010 Page 7.
= − k i j x
u
x
u k x
u
x
u u u ij i
j
j
i ij t i
j
j
t i
[7-19]
Here we have defined μt and νt as the turbulent dynamic and kinematic viscosities. Again, we are
hiding our ignorance of turbulent flows in these terms. We are left with equations that are similar
to the laminar flow equations, but we have to find a way to predict μt and νt. We see that the term
called the Reynolds stress terms. We see that they are symmetric, so there are only six different
Reynolds stress terms. We will show below that advanced turbulence models seek to predict
each Reynolds stress.
The mixing-length model of simple turbulent flows
Considerations of turbulent flows started with Reynolds’s work in the late 1800s and the
contemporary assumption by Boussinesq that the Reynolds stresses could be related to the
gradients of the mean flow velocity components. Subsequently in the early 1900s Prandtl
developed the model of the boundary layer in fluid flows. Von Karman later expanded this work.
This work developed Prandtl’s mixing-length theory as a useful approach for modeling simple
turbulent flows.
The mixing-length model remains a useful model for simple two-dimensional flows in jets and
near-wall regions. In addition, this basic analysis is used to develop the wall functions that are
used as boundary conditions in CFD turbulent flow models. In the mixing-length model, there is a
viscous (laminar) sublayer close to the wall. Beyond this sublayer, there is a transition region into
a fully turbulent boundary layer. The turbulent boundary layer near the wall is characterized by
the following quantities in a two dimensional flow where the x direction is the predominant flow
direction, with flow velocity u, and the y direction is the direction perpendicular to the flow. For
flows near walls, the point y = 0 is the wall.
τ τ
τ
u y
y
y y u
u u = w^ u = = =
[7-20]
Here, τw is the wall shear stress and the quantity uτ is called the friction velocity. The law of the
wall postulates that there is a relationship between the dimensionless velocity and distance
variables defined in equation [7-20].
u = f y [7-21]
For flows near a wall, this relationship has different forms. The region closest to the wall (y
< 5)
is the laminar sublayer in which the normal equations for viscous stresses are valid. Beyond this
region the values of u
are given by the following equations.
ln 500
ln( ) ln( )
max A y
y u
y
Ey B
y
y y
u
[7-22]
Page 7.8 L. S. Caretto, January 30, 2010 Turbulence Modeling
In equation [7-22], we have introduced the constants κ, known as the von Karman constant,
which has a value of 0.4, B = 5.5, E = 9.8, and A. The values of these constants are for smooth
walls. We have also introduced the boundary-layer thickness, δ, for the velocity profile away from
the boundary layer.
The results above have not made explicit use of the mixing-length; instead, it has shown the final
results. The mixing-length model is a way of determining a characteristic length, l, for the
turbulence. If we take a simple dimensional analysis of the turbulent kinematic viscosity, νt, with
dimensions of length squared over time, we see that this variable has the same dimensions as
the product of a velocity and a length. Thus we can write the turbulent kinematic viscosity as a
product of a velocity and a length. If we take the velocity as U and the length scale as l, we can
write:
νt = Cμ U l [7-23]
where = Cμ is an empirical constant.
The k-ε model of turbulent flow
The most popular model for modeling turbulence in common flows of engineering interest is
called the k-ε model. This model uses the kinetic energy of turbulence, defined in equation [7-8]
and the turbulence dissipation rate, ε. This dissipation rate measures the rate of energy transfer
from kinetic energy in the smallest eddies when they do work against the viscous forces.
Formally the dissipation is defined in terms of the deformation (or strain) rates, eij which are
defined below.
ij ij ij i
j
j
i ij sothat e x
u
x
u
The formal definition of the dissipation is shown in equation [7-25]. The dimensions of dissipation
are energy divided by (mass times time). (As noted in the paragraph following equation [7-8], we
conventionally define the energy per unit mass as the “kinetic energy”; that is being done here.)
The SI units for ε are m
2 /s
3 .
We are using the summation convention here. The implied summation over the two repeated
indices, i and j, gives nine terms in this equation. Also, ν in this equation is the actually kinematic
viscosity, not the turbulent kinematic viscosity.
The basic approach of equation [7-23] in which the turbulent velocity is treated as the product of a
length scale times a velocity scale is also used here. In the k-ε model the velocity scale is the
square root of the kinetic energy of turbulence; i.e., U = k
1/
. The length scale is equal to k
3/ /ε.
Substituting these terms into equation [7-23] gives the following equation for the turbulent
viscosity.
(^3222) 2
(^1) k C
k C
k t =^ C^ k = t = [7-26]
Page 7.10 L. S. Caretto, January 30, 2010 Turbulence Modeling
The first two terms on the right-hand side can be combined to give the equation for the kinetic
energy of turbulence.
k k i
t
i i
i P x
k
x x
uk
t
k [7-31]
The laminar viscosity in the diffusion term is much smaller than the turbulent viscosity (divided by
the Prandtl number for k) and is usually not included in the computations. Note that this equation
has the same form as the general CFD equation (transient – convection – diffusion – source) so
we perform the numerical analysis of this equation by the usual CFD algorithms for other
variables.
Although we defined the turbulent dissipation, ε, in equation [7-25], we do not use this equation to
compute ε. Instead, we have to solve another partial differential equation to compute this
quantity. The final equation differential equation for the dissipation, with all the modeling
assumptions already made, is shown below.
k
k
x x x
u
t
k i
t
i i
i
2
1 2
ε ε ε
[7-32]
The first three terms in this equation are the usual terms in the general balance equation,
representing the transient term, the convection term and the diffusion term. For the turbulent
flow, we use the mean velocities in the convection term and the turbulent transport coefficient in
the diffusion term. The effective transport coefficient is expressed as the ratio of the viscosity to
the Prandtl number for dissipation, σε. The last two terms in equation [7-32], with the empirical
constants Cε 1 and Cε 2 , represent the production and destruction of dissipation. This equation is
very similar to the corresponding equation [7-31] for the kinetic energy of turbulence.
The various empirical constants used in the standard k-ε model are summarized below.
These notes have described the original k-ε model which still has wide use in engineering
applications. Other versions of this model known as the renormalizable group (RNG) and the
realizable k-ε models use alternative derivations of the model to be more consistent with the basic
physics of turbulence. These models have different differential equations for the production and
dissipation terms and different equations for computing the turbulent viscosity from k and ε. The
empirical constants used in these models are also different. These models are recommended in
flows that high strain rates ( i.e. large velocity gradients).
Boundary conditions and wall functions
In applying turbulence models like the k-ε model it is necessary to specify boundary conditions for
the dependent variables, such as k and ε, at a variety of physical boundaries encountered in the
flow. The most complex application of boundary conditions is at solid walls. At large Reynolds
numbers, the laminar sublayer is so small that it is not possible to use a grid that is fine enough to
resolve this region. As a consequence, the boundary conditions at solid walls are usually
handled by the mixing length relations from equation [7-22].
The starting point for the derivation of wall functions is the part of equation [7-22] which is valid for
y
values between 30 and 500. Here we use the second form of that equation.
Turbulence Modeling L. S. Caretto, January 30, 2010 Page 7.
τ
τ
Ey Eyu
u
u u ln
ln( ) 1 [7-34]
In the general consideration of wall functions, u
is regarded as a normalized velocity component
parallel to the wall and y
is regarded as a normalized distance perpendicular to the wall. To use
this equation we need a value for the wall shear stress, τw, which is used in the normalizing
and y
.) This parameter is usually calculated from the equilibrium assumption that the production
and dissipation of turbulence are equal. In this case one can derive the following relationship
between uτ and the turbulent kinetic energy.
u C^4 k
1 τ =^ μ [7-35]
We can combine equations [7-34] and [7-35] to get an expression for the velocity at the first node
in from the wall.
μ μ
τ τ
EyC k C k
Eyu u u
4
1
4
1 ln ln [7-36]
In this equation Cμ, and E are known dimensionless constants, 1 ν is the (laminar) kinematic
viscosity, and k and y are the values, respectively, of turbulent kinetic energy and distance from
the wall at the first node. Other wall functions are used for k, ε, and other transported variables
such as temperature and species concentrations.
At outlets, the gradients of turbulent kinetic energy and dissipation are assumed to be zero in the
direction of the flow. At inlets it is necessary to specify profiles of k and ε. Ideally these should
come from measurements on flows similar to the one being modeled. If such data are not
available, the inlet values of k and ε can be estimated form the following equations.
l
l
k k Cuinlet L C^4
2 3
Here L is a characteristic length, which may be taken as the hydraulic diameter of an inlet.
(Recall that the hydraulic diameter is 4A/P, where A is the area and P is the perimeter. For
circular flow passages this is the diameter.) The constant C is a small number in the range of
0.01 to 0.05. Another way to specify k is to specify the turbulence intensity. This is defined as
the ratio of the kinetic energy of turbulence to the mean velocity squared. This is equivalent to
allowing the user to specify C in the initialization equation for k above
The wall functions proposed above are limited to “large” Reynolds numbers. Separate wall
functions are required for lower Reynolds numbers. The renormalizable group (RNG) and the
realizable k-ε models, mentioned above as alternative k-ε models have different initial condition
specifications and different wall functions.
1 E = 9.8 and for the k-ε model, Cμ = 0.09.
Turbulence Modeling L. S. Caretto, January 30, 2010 Page 7.
Recall that we defined the general coefficient γ
(φ) = Γ
(φ) /(ρc) in equation [7-11]. For momentum
transport (where c = 1), γ is the kinematic viscosity, ν, and Γ is the dynamic viscosity, μ. Although
the Prandtl number is formally defined as μc (^) p/k = ν/α, in turbulence models, the turbulent
“Prandtl” number for a particular flow variable, φ, is defined as the ratio ν/γ
(φ) and is given the
symbol σ (^) φ. (Here I put quotation marks around Prandtl to indicate that the name Prandtl number
is applied to ratios other than the original definition of the Prandtl number as ν/α. With this
definition of the general turbulent Prandtl number we can rewrite equation [7-17] as follows.
turb i i
turb
i lam
i x x x
u
∂
∂ϕ
∂
σ
ν
σ
∂
∂ ϕ
ϕ , ϕ,
[7-38]
As noted above, the turbulent transport property is usually much larger than the laminar transport
property and is often ignored in computations. The laminar “Prandtl” number is a fluid property
that can be found from data supplied with the CFD codes or entered as input by the user from
property tables. The turbulent “Prandtl” number is usually set as a default value in the CFD code.
These default values are usually satisfactory for most applications. Users can override these
default values if they have other data available to them to justify this.
Other turbulence models
Several other models that solve two partial differential equations have been proposed, and each
has their proponents. However none of these seems to have received better acceptance than the
k-ε model.
For complex flows, typically those with large amounts of swirl, the Reynolds Stress model is used.
This model solves six partial differential equations, one for each independent Reynolds stress
term. In addition, the dissipation equation [7-32] must be solved to obtain terms used in the
differential equations for the Reynolds stresses. The use of the Reynolds stress model has not
been completely successful. For some flows this model gives much better agreement with
experimental data as compared to the k-ε model. For other flows, little improvement is obtained
by going from the k-ε model to the Reynolds stress model.
The general equation for the Reynolds stress model is shown below.
ij ij
i ij j
j i
i j
k
i j i j t
x
u u u x
u u u x
u u
x x
uu u
t
u u − +Π +Ω ⎟
l
l l
l l l l
l [7-39]
Transient Convection Diffusion Production Dissipation Pressure Strain Rotation
The first five terms appearing in full in the above equation, reading from left to right, are the usual
transient, convection, diffusion, production and dissipation terms. The remaining terms, the
pressure-strain term, Πij , and the rotational term Ωij are described below. The δij term is the
Kronecker delta. All other terms in the above equation have been defined previously. Note that
we have a repeated sum over the repeated indexl. The i and j indices represent a particular
Reynolds stress and we have to solve equation [7-39] for the six unique Reynolds stresses.
The pressure strain interaction terms are important ones in this equation, but they are also the
most difficult ones to model accurately. They come about from pressure fluctuations that occur
when two eddies interact with each other and when an eddy in one flow region interacts with the
main flow in another region. The simplest form of this term is
Page 7.14 L. S. Caretto, January 30, 2010 Turbulence Modeling
ij k ij
i j
j ij i j ij i P x
u u u x
u u u k C u u k
l
l l
l [7-40]
In this equation, C 1 = 1.8, C 2 = 0.6, and the Pk term is the total production of kinetic energy
defined in equation [7-29]. The rotational term uses the local rotation vector, ω, whose
components are ωk. It also uses the symbol εijk which is zero if any two indices are the same and
+1 if i = 1, j = 2, and k = 3. For all other combinations of ijk, εijk = -1 if the indices ijk are an odd
permutation of ijk or +1 if the indices are an even permutation of ijk. (In an odd or even
permutation an odd or even number of moves are required to get 123 into the given set of
indices. For example 321 is an odd permutation since it only takes one move, swapping the 3
with the 1. The order 231 is an even permutation since we need two moves. We first swap the 1
and the 3 (giving 321) then we swap the 2 and the 3 to get the 231 result. Other orders of swaps
are possible, but it will always take an even number of swaps to permute 123 into 231. With
these definitions, the rotation term is
This term has an implied summation over the two repeated indices k and l.
In the Reynolds stress model, the kinetic energy of turbulence is found from the sum of the three
Reynolds stresses with the same index.
2 3
2 2
2 u 1 u u u u k
j j + + = = [7-42]
Because k can be found from this equation, it is not necessary to solve a separate equation for k.
However, it is still necessary to solve an equation for the dissipation, ε. Thus, the Reynolds
stress model requires the solution of seven partial differential equations to compute the
turbulence properties.
An approximate form of the Reynolds stress model, known as the algebraic Reynolds stress
model, has been used. In this model the source terms for production and destruction of the
Reynolds stress terms are assumed to be equal in magnitude. This means that it is not
necessary to solve differential equations. Instead, a system of six simultaneous linear equations
is solved at each node for the Reynolds stresses at that node. It is still necessary to solve a
differential equation for the kinetic energy and the dissipation.
A new model called the v2f (or v
2 -f) model has been used to allow improved results (as compared
to the k-ε model) for low Reynolds number turbulent flows.
The Spalart-Allmaras model is a one-equation model that solves a transport equation for the
turbulent viscosity. It was designed specifically for aerospace applications involving flows along a
surface. The detached eddy simulation (DES) model is a modified version of the Spalart-
Allmaras model that is designed to be a less expensive alternative to LES.
These notes have not discussed turbulence models for compressible flows. In such models it is
necessary to consider density fluctuations. One approach for doing this is the use of Favre
averages in which the instantaneous flow quantity is split into two parts as follows.
Page 7.16 L. S. Caretto, January 30, 2010 Turbulence Modeling
external flows the Spalart-Allmaras model should be used. This model may be used with wall
functions or resolution of the viscous sublayer. Chapter 12 of the Fluent Users’ manual gives a
full discussion of turbulence models and guidance for when to apply a particular model and the
types of boundary conditions that can be used with a particular model.
Some references on turbulence and turbulence models
As noted above, Chapter 12 in the user’s guide for Fluent discusses the basic turbulence models
and provides guidance for their use in Fluent. This manual is available from the help menu in
Fluent.
The two main references used for this set of notes were the chapters on turbulence modeling in
the following texts:
J. H. Ferzinger and M. Perić, Computational Methods for Fluid Dynamics , (third edition) Springer,
H. K. Versteeg and W. Malalasekera, An Introduction to Computational Fluid Dynamics. The
Finite Volume Method ” (Second edition), Prearson, Prentice-Hall, 2007. (Chapter three on
turbulence and its modeling.)
The Ferzinger and Perić text has more discussion of advanced methods such as DNS and LES.
The Versteeg and Malalasekera text has more background on the fundamentals of turbulence.
The two classical references in this field are:
H. Tennekes and J. L. Lumley, A First Course in Turbulence , The MIT Press, 1972.
B. E. Launder and D. B. Spaulding, Mathematical Models of Turbulence , Academic Press, 1972.
(also, “The Numerical Computation of Turbulent Flows,” Computer Methods in Applied Mechanics
and Engineering , 3 :269, 1974.)
I have not reviewed the text listed below, but the first edition of this text has some good online
reviews.
D. C. Wilcox, Turbulence Modeling for CFD , (second edition) DCW Industries, Inc., 1998. See
the book web site: http://www.dcwindustries.com/books/0963605151.htm