1

Recent Developments in Eddy Viscosity Modelling of

Review Turbulence
Kazuhiko Suga

Abstract

Eddy viscosity modelling is still a "standard" approach for industrial CFD applications for turbulent flows
despite its serious deficiencies. Thus, a number of research studies, including the author's, have been
recently made to improve models of this kind. This article reviews these efforts and suggests a future
direction for tackling turbulent flows of industrial importance.

Keywords CFD, Modelling turbulence, Eddy viscosity, Near-wall model,

Nonlinear stress-strain relation

(EVM), Reynolds Stress Model (RSM), Large Eddy

1. Introduction
Simulation (LES) and Direct Numerical Simulation
Turbulence is an irregular motion in fluid flows. (DNS). Amongst these, only DNS closely simulates
The various flow quantities thus show random the actual physics of turbulence, while industrially
variation with time and space and only statistically interesting schemes are at present the simplest
averaged values are distinctively discerned. Since EVM's because the others require far more
exactly dealing with turbulence with mathematics computer resources than are available for routine
has been one of the most notoriously thorny work.
problems of classical physics, several well known Although the eddy viscosity concept assumes a
scientists, who had interests in turbulence, did not crude relation between turbulent quantities, models
dare to pursue its physics. For example, W. based on this concept such as the k-ε two-equation
Heisenberg, the great Nobel Prize laureate, did his model have made many successful predictions in
doctoral research on turbulent flows but later many flow fields with numerical stability. The k-ε
changed his major to quantum theory. A. Einstein, EVM is thus the main scheme for routine work in
the greatest Nobel Prize physicist in this century, present industrial laboratories. Nevertheless, over
knew the difficulty of turbulence and thus did not go the wide range of flow predictions made over the
into its physics. R. Feynman, another famous Nobel last two decades, it is now recognized that
Prize laureate, noted that turbulence was the most conventional EVM's have severe defects in many
important unsolved problem of classical physics. complicated flow fields. Consequently, research
However, the development of modern computer studies aimed at extending the applicability of
technology drastically changed this situation and has EVM's have been highly demanded by industry. In
provided some opportunities even for industrial fact, the continuous research efforts have been
engineers to challenge turbulent flows using significantly extending the performance of EVM's.
Computational Fluid Dynamics (CFD). The major Therefore, this article particularly focuses on these
CFD treatments of turbulence can be classified into recent achievements.
four types of approaches: Eddy Viscosity Model The following section §2 surveys the main

R&D Review of Toyota CRDL Vol. 33 No. 1 ( 1998.3 )

2

historical establishments related to eddy viscosity transport equations for the turbulent kinetic energy
modelling, section §3 summarizes the recent efforts, k and its dissipation rate ε are solved with
then finally, section §4 concludes and suggests a approximations.
future direction for the treatment of industrially
Dk = Dk – ui uj ∂Ui – ν ∂uj ∂uj
important turbulent flows. ·······················(6)
Dt ∂x j ∂x i ∂x i









2. Historical foundations Pk ε

In this section, firstly, major establishments in

Dε = Dε + cε 1Pk – cε 2ε
eddy viscosity modelling are surveyed, then a more ···························(7)
Dt k/ε
extended modelling concept that forms algebraic
Reynolds stress models (ASM's) is summarized where D φ is the diffusive transport term of the
because it is mathematically very close to EVM's variable φ which is normally modelled as:
and links some recent nonlinear eddy viscosity ∂ ∂φ
approaches. Dφ = ν + νt . ························(8)
∂x j σφ ∂x j
2. 1 Eddy viscosity models
The EVM's are based on an algebraic expression The coefficients cε1 and cε2 were given by referring
which represents the Reynolds stresses appearing in to the measured rate of decay in grid turbulence and
the ensemble averaged Navier-Stokes equations as the local equilibrium turbulence, respectively. These
unknown properties. The ensemble averaged forms standard values 3) are listed on Table 1. This
of the transport equations for mean velocity of standard version does not have any near-wall
incompressible flows with constant properties can be dependence upon molecular viscosity, so that wall
written as: functions4,5) are employed in place of the no-slip
Continuity, wall boundary condition.
∂Ui To take account of the viscous effects, Jones &
=0 ······································(1) Launder2) (JL) first devised a low-Reynolds-number
∂x i
(LRN) version of the k-ε EVM. They implemented
Momentum, the near-wall viscous effects by damping νt toward a
DUi = – 1 ∂P + ∂ ν ∂Ui – ∂ ( u u ) wall with introduction of a damping function, fµ, as:
i j
Dt ρ ∂x i ∂x j ∂x j ∂x j
{

2
Reynolds stress ν t = cµf µ k . ······································(9)
ε
·······························(2)
This fµ function was designed to reduce its value
where ρ and ν are, respectively, the density and the from unity toward a wall. They also modified the ε
kinematic viscosity of the fluid of interest. equation with the introduction of the other damping
Following Boussinesq1), the Reynolds stresses are functions, f1 and f2, as:
represented by the eddy viscosity νt and the strain
Dε = Dε + cε 1f 1Pk – cε 2f 2ε + Pε 3 ···············(10)
tensor Sij as:
Dt k/ε
ui uj = 2 kδ ij – ν t S ij . ~
3
································(3) where ε is the isotropic part of ε defined as
. The reason why they chose ~ε
1/2 2
ε ≡ ε – 2ν ( ∂k∂x )
j
The strain tensor is defined as: rather than ε itself is that ~
ε vanishes to exactly 0 at a
∂Ui ∂Uj
wall boundary. This simple boundary condition
S ij ≡ + . ································(4) makes numerical solutions more stable. The
∂x j ∂x i

The k-ε EVM2) takes the eddy viscosity as:

Table 1 Empirical constants for the k-ε EVM.
2
ν t = cµ k ································(5)
ε
cµ cε1 cε2 σk σε
where cµ is a constant value given by referring to 0.09 1.44 1.92 1.0 1.3
local equilibrium shear layers. To obtain νt, the

R&D Review of Toyota CRDL Vol. 33 No. 1 ( 1998.3 )

 3

gradient production term, Pε3, was also modelled Eq.(15), one may rewrite Eq.(14) as:
using a gradient diffusion hypothesis as:
uu
2 2 Pij + Πij – εij = ki j ( Pk – ε ) ···············(16)
∂ Ui
Pε 3 = 2 νν t . ·······························(11) Then, one just needs models for the terms Πij and εij.
∂x j ∂x k There are many established models for Πij such as
Many versions followed this original work. Some Launder, Reece & Rodi model10), and each model
researchers such as Wilcox6), for example, chose a forms a different version of the ASM. However, the
substitute quantity for ε . This model solves a basic model11) may be
1
modelled equation for ω ( ≡ ε / k ): Πij = –c1εaij – c2 ( Pij – δijPkk )
3
····································(17)
Dω = Dω + cω1Pk ω – cω2ω 2. ······················(12) 2
Dt k where the anisotropic stress tensor aij = uiuj / k – δij
3
and the values of 1.8 and 0.6 are normally used for
The coefficients such as cω1 and cω2 were tuned by the coefficients c 1 and c 2, respectively. In high
referring to the k- ε model equations since the ω Reynolds number isotropic flows, the following
equation was derived by manipulating the k and ε treatment is normally applied.
equations.
Patel et al. 7) concluded, however, following εij = 23 εδij ···································(18)
systematic comparisons between eight LRN models,
that an amended version of the JL model by Launder Consequently, the algebraic expression for the
& Sharma8) (LS) was one of the most successful for a Reynolds stress is obtained as:
number of straight thin shear flows. The LS model
uses the following damping functions of a turbulent
~ ( 1 – c2 ) ( Pij – 23 δ ij Pk )
Reynolds number Rt ( ≡ k2 / ( νε ). ui uj = k 2 δ ij + .
3
( c1 – 1 + Pk / ε ) ε

 fµ = exp { –3.4
} ·································(19)
2
 ( 1 + Rt / 50 )
 f1 = 1.0 Since Pij and Pk consist of the Reynolds stresses

 f2 = 1.0 – 0.3exp ( –Rt2) and the mean velocity gradients, this algebraic form
is implicit in terms of the Reynolds stress.
····································(13)
Therefore, ASM's need to solve the transport
2. 2 Algebraic Reynolds Stress Models equations of k and ε with successive matrix
The transport equation of the Reynolds stress uiuj is inversions of the implicit algebraic equation set for
the Reynolds stresses.
Dui uj ∂Uj ∂Ui
= Dij – ui uk + uj uk + Πij – εij
Dt ∂x k ∂x k 3. Toward a new standard











Pij In this section, the recent novel attempts to

·····························(14)
improve the EVM's are discussed especially
where Π ij and ε ij are, respectively, the pressure focusing on wall detecting parameters and nonlinear
correlation and the dissipation rate terms of uiuj. constitutive relations used in the models. However,
The terms that contain the gradient of uiuj are the due to page limitations, unfortunately, many of their
transport terms (i.e., the convection and the diffusion model equations are not described, therefore, the
terms). In the ASM scheme, the transport of the referenced papers should be consulted for more
Reynolds stresses is approximated in terms of that of details.
turbulence energy k to reduce the differential 3. 1 Near-wall modelling
equations for uiuj to a set of algebraic ones. This Since the near-wall variation in the local turbulent
scheme was firstly introduced by Rodi9) as: Reynolds number R t significantly changes
depending on the bulk Reynolds number Re as
ui uj
Tij = Tk ·······························(15) shown in Fig. 1‡, finding another near-wall detector
k
which has general near-wall characteristics has been
where Tφ is namely the net transport (convection the main concern of modelling near-wall turbulence.
minus diffusion) of φ. Since Tk =Pk – ε, by adapting ‡
The results of DNS can be treated as almost exact physics.

R&D Review of Toyota CRDL Vol. 33 No. 1 ( 1998.3 )

4

In their review, Patel et al. 7) emphasized the al. 21) using their DNS results suggested that
necessity to have a reasonable near-wall f µ constructing a universal model depended on
distribution because none of their cited models identifying dimensionless parameters such as the
agreed with the data deduced from several different normalized strain invariant:
experiments (Fig. 2). This triggered many
S≡τ S ij S ij / 2 ····························(22)
researchers to modify the LRN models, suggesting
ways of improving the EVM's. Moreover, since the where τ is the characteristic time scale normally given
emergence of the DNS12-14), a lot of attention has as the ratio of k and ε, k/ε. In response, in order to
been given to the near-wall asymptotic behaviour of obtain the substitute parameter for the wall distance,
each turbulent quantity because the DNS provided Yang & Shih22) (YS) devised the parameter R which
reliable data for every process including consisted of the strain invariant as:
unmeasurable correlations. k 2 / ( νε )
R≡ k = = Rt / S .
To obtain a reasonable near-wall distribution of ν S ij S ij / 2 ( k / ε ) S ij S ij / 2
the fµ damping function, many recent versions of the ····························(23)
LRN k-ε EVM15-17) have implemented the effects of
the dimensionless wall distance (wall unit): The use of this in the dumping function fµ led to
y+ ≡ uτ y/ν ···································(20) good predictive performance in wall shear flows
where u τ and y are the friction velocity and the with zero or favorable pressure gradients.
distance from a wall, respectively. However, one Due to its general characteristics in shear flows as
can easily see that none of these LRN EVM's is shown in Fig. 3, the strain invariant has been
useful to apply for a flow with a recirculation. In recently employed as a near-wall parameter in
particular, the use of uτ is not suitable for such a several other proposals such as Cotton & Ismael23)
flow case because it becomes zero at a reattaching
point. In this case, the wall unit was sometimes
replaced with R y ( ≡ k y / ν ). Abe et al.18), however, 1
replaced u τ with the Kolmogorov velocity scale, 0.8
(νε)1/4, and devised the parameter:
0.6
y* ≡ ( νε )1/4 y/ν ··································(21) fµ
to damp the eddy viscosity in order to obtain 0.4

reasonable predictions of backward facing step 0.2

flows. (They later extended the model to a nonlinear 0
k-ε model19)) The parameter y* was also used in 10 20 30 40 50 60 70 80 90 100
y+
Kawamura & Kawashima's (KK) LRN k-ε EVM20).
Nevertheless, the use of the wall distance y limits
the model's applicability when considering flow in Fig. 2 Experimentally suggested distribution
7).
more complex geometry. The discussion by Lee et of the damping function, fµ

500 20
Re=2800 Re=2800
400 Re=7000 Re=7000
15
300
Rt S 10
200
5
100

0 0
0 20 40 60 80 100 0 20 40 60 80 100
wall
+ +
y y

Fig. 1 Turbulent Reynolds number in plane Fig. 3 Normalized strain invariant in plane
channel flows by DNS12,13). channel flows by DNS12,13).

R&D Review of Toyota CRDL Vol. 33 No. 1 ( 1998.3 )

 5

and the nonlinear k- ε EVM of Craft, Launder & necessarily the wall-boundary. To support this, Fig. 5
Suga24,25). Cotton & Ismael later proposed a k-ε-S shows a similar damping profile of the eddy
model26) coupling with a transport equation for the viscosity near the free surface (y/ δ ) of an open
strain invariant: channel flow. Hence, directly implementing this
effect by the use of υ2 in the damping model has a
DS = DS + 0.5 τS ij2 – S .
Dt τ
····························(24) physically correct reason. However, the k- ε - υ 2
model is only applicable in a flow parallel to a wall
The transport effects of S gave some reliability in because υ is not always normal to a wall in
predicting buoyant flows. complicated geometry. Furthermore, the use of υ2
Since the desirable variation in the strain invariant alone in a scalar variable ν t leads to severe
as a near-wall parameter is relatively limited near fundamental inconsistencies since υ2 is a component
walls, Craft, Launder & Suga27-29) further introduced of the Reynolds stress tensor and should not appear
the stress invariant A2 ( ≡ aijaij) as another near-wall in any scalar value.
detector into their nonlinear k- ε EVM. Stress Thus, it is necessary for a more general eddy
anisotropy is high near a wall and its measure is viscosity formula to have a physically and
represented by A2 as shown in Fig. 4, hence, A2 can mathematically correct damping parameter toward
be a near-wall parameter. The value of A 2 was wall or shear-free boundaries. Accordingly, the
obtained by solving its transport equation: author noticed the flatness parameter of the
Reynolds stress tensor, A (≡ 1– 89 (aijaij - aijajkaki)), as
DA 2 = – 2 A 2 D + 2 aij D – 2 A 2 P
k ij k a damping parameter. He thus extended his work on
Dt k k k
a a a the k- ε -A 2 model 27-29) to a k- ε -A model 34) by
+ 2 ij Pij + 2 ij Πij + 2 A 2 ε – 2 ij εij
k k k k substituting the following A-transport equation for
····························(25) the A2 equation.
DA
with proper models for Πij and εij. = – 9 ( 32 A3Dkk + 2aijDij – 3ajkakiDij )
Dt 8k
Durbin 33) introduced the Reynolds stress
component normal to a wall, υ 2 , as a damping – 9 ( 32 A3Pkk + 2aijPij – 3ajkakiPij )
8k
parameter for the eddy viscosity of the k-ε EVM as:
νt = cµυ2τ ·································(26) – 9 ( 32 A3Πkk + 2aijΠij – 3ajkakiΠij )
8k
The values of υ are obtained by solving its
2

modelled transport equation. In fact, the near-wall + 9 ( 32 A3εkk + 2aijεij – 3ajkakiεij )

8k
damping in the eddy viscosity comes from the ·································(27)
blocking effect on the fluctuating velocity Because A is a scalar and vanishes at the wall and
component normal to a wall-boundary by the shear-free boundaries as shown in Fig. 6, its
existence of the boundary. In this sense, the
boundary which gives the blocking effect is not
0.12
0.1
0.08
νt/(δuτ)

2
Re=2800 0.06
A2 Re=7000
1.5 0.04
0.02
A,A2 1
A 0
0 0.2 0.4 0.6 0.8 1
0.5 wall symmetry plane
y/δ shear-free boundary

0
0 20 40 60 80 100
Fig. 5 Eddy viscosity:
y+
, Open-channel30) ;
, Couette-Poiseuille32) ;
Fig. 4 Stress invariants in plane channel flows , Channel (Re =2800)12) ;
by DNS12,13). , Channel (Re =7000)13).

R&D Review of Toyota CRDL Vol. 33 No. 1 ( 1998.3 )

6

introduction into the damping function of the eddy In fact, the cubic c4~c7 terms have sensitivity to
viscosity allows one to form a physically and swirl and streamline curvature. They afterwards
mathematically correct model. modified their cubic NLEVM coupling with the
3. 2 Nonlinear eddy viscosity modelling effects of A 2 to correctly mimic near-wall
Another important topic in the recent eddy turbulence27-29).
viscosity modelling is a nonlinear extension of the Pope 36) showed that the generalized nonlinear
(linear) stress-strain relation, Eq. (3). This approach stress-strain relation was mathematically equivalent
forms a nonlinear eddy viscosity model (NLEVM). to an explicit form of the ASM. He generalized the
Note that a sort of the NLEVM is sometimes called nonlinear constitutive relation using the Cayley-
an explicit ASM due to its optimization process for Hamilton theorem and solved a matrix obtained by
the coefficients. substituting the Reynolds stresses in Eq. (19) with
The original linear stress-strain relation does not the constitutive relation. Although he outlined the
produce meaningful differences between the normal procedure to obtain the coefficients, he was not able
stresses. For example, in shear flows where only S12 to provide the coefficients generally due to the
is nonzero, Eq. (3) leads to isotropic turbulence as: complexity of the algebra. In fact, the generalized
u1u1 = u2u2 = u3u3 =
2
k constitutive relation includes up to fifth-order
3 ·························(28)
products of strain and vorticity tensors. Recently,
while the values of the normal stresses are very following Pope's methodology, Taulbee42) and Gatski
different from one another in actual flow cases. & Speziale43) proposed elaborate coefficients for the
Thus, the linear model lacks the capability of three-dimensional flows. Their NLEVM's (explicit
predicting anisotropic turbulence in many ASM's) thus include up to fifth-order terms.
industrially important flows such as turbulence- However, the roles and necessity of fourth- and
driven secondary flows, swirling flows, etc. fifth-order terms have never been clarified.
Although the ideas of NLEVM themselves Very recently, the author pointed out an inherent
emerged back in the 70's 35,36), until recently, the defect in the stress-strain relation and tried to
models of this type were not widely explored. Many remove it. In shear-free turbulence appearing, for
attempts at developing and using such schemes have example, near the free surface of an open channel
been recently made 37-41) . They all introduced flow, all strain and vorticity tensor components
quadratic terms into Eq. (3) as:
vanish. The linear and nonlinear stress-strain
1
aij = – cµτSij + c1τ2 ( SikSkj – S S δ ) relations thus always return isotropic turbulence
3 kl kl ij
+ c2τ ( ΩikSkj + ΩjkSki ) there (e.g., all the terms on the right hand side of Eq.
2

1
+ c3τ ( ΩikΩjk – ΩlkΩlkδij )
2 (30) become 0 in this case) while the actual
3
·························(29) turbulence is significantly anisotropic. Therefore,
the author introduced the following additional term
where the vorticity tensor, Ωij ≡ ∂Ui / ∂xj – ∂Uj / ∂xi.
The quadratic c1~c3 terms produce discrepancies
between the normal stresses. These quadratic 1
NLEVM's thus successfully reproduced turbulence 0.8
driven secondary flows, however, they did not have 0.6
sensitivity to streamline curvature (including swirl). A
0.4
Therefore, in order to capture the streamline
24)
curvature effects, Craft, Launder & Suga further 0.2
introduced cubic terms as: 0
0 0.2 0.4 0.6 0.8 1
aij = – cµτSij + c1τ ( SikSkj – 1 SklSklδij )
2
3 y/δ
+ c2τ2 ( ΩikSkj + ΩjkSki )
+ c3τ2 ( ΩikΩjk – 1 ΩlkΩlkδij )
3 Fig. 6 Stress flatness parameter:
+ c4τ3 ( SkiΩlj + SkjΩli ) Skl
, Open-channel31);
+ c5τ3 ( ΩilΩlmSmj + SilΩlmΩmj– 2 SlmΩmnΩnlδij )
3 , Couette-Poiseuille32) ;
+ c6τ3SijSklSkl + c7τ3SijΩklΩkl. , Channel (Re =2800)12) ;
·························(30) , Channel (Re =7000)13).

R&D Review of Toyota CRDL Vol. 33 No. 1 ( 1998.3 )

 7

Aij composed of the gradient of the stress flatness stress distributions with the DNS13) data near the
parameter A into the cubic stress-strain relation: Eq. wall. All the models reproduce the DNS results
(30). quite well, though the profile by the W92 model
∂ Ak ∂ Ak ∂ Ak ∂ Ak distinctively deviates from the data in the region y+< 20.
A ij = ca' τ 2 – 1 δ ij
∂x i ∂x j 3 ∂x k ∂x k The predicted turbulence energy distributions
·························(31) shown in Fig. 8, however, display an interesting fact.
The recently proposed nonlinear SA and ARG
Since the distribution of A has a steep gradient near models predict very similar profiles to that of the
the shear-free boundary as shown in Fig. 6, this rather dated NY model and they are poorer than that
additional term does produce anisotropy of of the 24-year-old LS model. Except for them, the
turbulence there. The author showed its usefulness recent versions of EVM's have shown quite
for capturing shear-free turbulence combining it with successful performance. In fact, many of them are
the k-ε-A three equation NLEVM34) in excellent agreement with the DNS.
3. 3 Comparisons of model performance Since modelling the ε (or ω) equation is much
This subsection displays the near-wall more difficult than modelling the k equation, thus
performance of typical linear and nonlinear EVM's many of the predicted ε distributions poorly accord
listed in Tables 2 and 3. All the models listed, with the DNS data as seen in Fig. 9. Nonetheless,
except for the W92 44), the SA 45) and the ARG 46) the result of the nonlinear CLS model shows quite
models, have already been discussed or referred to. excellent agreement with the DNS, and those of the
The W92 model is the latest version of the linear k-ω linear KK and the nonlinear AKN models are also
model. Among the nonlinear EVM's, the SA and the fairly acceptable.
ARG models are, respectively, a k-ε model and a k-ω Another important feature of a LRN model is grid
model based on the nonlinear stress-strain model of dependency on the predicted results. Fig. 10 shows
Gatski & Speziale43). the grid dependency on the predictive performance
Fig. 7 compares the predicted turbulent shear of the mean velocity in the pipe flow measured by
Linear EVM's
1
Table 2 Linear EVM's.
0.8
Model transport near-wall
0.6 Kim(1989),
−uv+

variables parameters LS,

W92,
k, ε~
8)
LS: Launder & Sharma(1974) Rt 0.4 YS,
RM,
W92: Wilcox(1992) 44)
k, ω Rt KK,
0.2
YS: Yang & Shih(1993)22) k, ε Rt, S CI,

RM: Rodi & Mansour(1993)17) k, ε Rt, y+ 0

KK: Kawamura & Kawashima(1994)20) k, ε~ Rt, y* 0 20 40 60 80 100
~ y+
k, ε, S
26)
CI: Cotton & Ismael(1994) Rt, S
NLEVM's
1

Table 3 Nonlinear EVM's. 0.8

Model transport near-wall 0.6 Kim(1989),

−uv+

CLS,
variables parameters ARG,
0.4 SA,
k, ε
38)
NY: Nisizima & Yoshizawa(1987) Rt, y+ AKN,
MK,
MK: Myong & Kasagi(1990)40) k, ε Rt, y+ 0.2 NY,
k, ε
45)
SA: Speziale & Abid(1995) Rt, Ry
0
ARG: Abid, Rumsey & Gatski(1995)46) k, ω Rt 0 20 40 60 80 100
AKN: Abe, Kondoh & Nagano(1995)19) k, ε Rt, y* y+
~
CLS: Craft, Launder & Suga(1995)29) k, ε, A2 Rt, S, A2
~
Suga(1997)34) k, ε, A Rt, A
Fig. 7 Turbulent shear stress.

R&D Review of Toyota CRDL Vol. 33 No. 1 ( 1998.3 )

8

Laufer47). The solid lines noted as 100% are the apparently shown. As mentioned in §3.2 and clearly
results with a fine enough grid whose first grid node shown in Fig. 14, however, the author's k- ε -A
is located just under unity of the wall unit (y1+< 1.0). NLEVM can capture stress anisotropy near the free
The lines noted as x% are the results using a grid surface while the CLS model cannot. This model
whose grid node density normal to the wall is x% of performance of the k-ε-A model is believed to be
that of the fine enough grid. Obviously, the LS very useful if the model is used to calculate heat and
model is very sensitive to the grid density and many mass transfer through a shear-free interface which is
of the other LRN models need at least a 50% grid one of the key phenomena of the environmental
node density of the fine enough grid. (The first grid issues.
node's y+ of this 50% grid is about 2.0: y1+~2.0.) 4. Conclusions
The nonlinear CLS, MK and NY models, however,
show equivalent performance even with the 40% The following aspects may be summarized
grid distributed from y1+~4.0. On the whole, it can through this review covering the recent research
be said the CLS model shows the best performance works on the eddy viscosity modelling of turbulence.
in the models compared in Fig. 10 in terms of the 1. Until recently, the wall distance y was often
predictive accuracy and the grid sensitivity. used in the low-Reynolds-number eddy viscosity
Fig. 11-13 show the predicted near-wall turbulent models, while the use of y limited the model's
intensities by the NLEVM's compared with the applicability.
DNS13) data. The CLS model clearly demonstrates 2. Many researchers have started to find general
the best performance while the other models do not local parameters for detecting wall effects. The
successfully mimic the stress anisotropy. proposed near-wall invariant parameters so far are
The near-wall performance of the author's k-ε-A the strain invariant S, the stress invariant A2 and the
three equation NLEVM is comparable to that of the stress flatness parameter A. Although they require
CLS (k- ε -A 2 ) model though it has not been
Linear EVM's
Linear EVM's 0.3
5 Kim(1989),
0.25 LS,
W92,
4 0.2 YS,
RM,
KK,
3 ε+ 0.15 CI,
+
k Kim(1989),
LS, 0.1
2
W92,
YS, 0.05
1 RM,
KK, 0
CI, 0 20 40 60 80 100
0
0 20 40 60 80 100 y+
y+

NLEVM's NLEVM's
5 0.3
Kim(1989),
0.25 CLS,
4 ARG,
0.2 SA,
AKN,
3 MK,
k+ ε+ 0.15 NY,
Kim(1989),
2 CLS, 0.1
ARG,
SA,
1 AKN, 0.05
MK,
0 NY, 0
0 20 40 60 80 100 0 20 40 60 80 100
y+ y+

Fig. 8 Turbulence energy. Fig. 9 Turbulence dissipation rate.

R&D Review of Toyota CRDL Vol. 33 No. 1 ( 1998.3 )

 9

solving their transport equations, the models streamline curvature and swirl effects.
including their effects showed encouraging results. 4. The combined effects of the new local near-wall
3. The use of nonlinear terms in the stress-strain parameters and nonlinear stress-strain relations have
relation is essential to predict complex strain fields. significantly extended performance of the eddy
Moreover, the cubic terms are necessary to mimic viscosity scheme. In particular, the use of A has

Re=45000; Linear EVM’s Re=45000; NLEVM’s

25 25

20 20
CI CLS
15 15

10 10
KK AKN
5 5

0 0
RM ARG
5 5

0 0
U+

U+
YS SA
5 5

0 0
W92 MK
5 5

0 0
LS NY
5 5

0 0

5 5

0 0
1 10 100 1000 1 10 100 1000
y+ y+

Fig. 10 Mean velocity distributions:

Symbols, Expt.47); , 100%; , 50%; , 40%; ··········, 35%.

Channel, Re=7000; CLS & MK Channel, Re=7000; AKN & NY

3 3

2 2
v’+ w’+ u’+
v’+ w’+ u’+

1 1

0 0
0 20 40 60 80 100 0 20 40 60 80 100
y+ y+

Fig. 11 Turbulent intensities: Fig. 12 Turbulent intensities:

Symbols, DNS13); Symbols, DNS13);
29)
, CLS ; , MK40). , AKN19); , NY38).

R&D Review of Toyota CRDL Vol. 33 No. 1 ( 1998.3 )

10

extended the applicability of the eddy viscosity

References
model toward capturing shear-free turbulence which
is very important for environmentally oriented 1) Boussinesq, J. : "Théorie de l'ecoulement tourbillant",
issues. Mém. Prés. par div. savant à l´acad. Sci. Paris,
23(1877), 46
5. Overall, it may be concluded that the nonlinear 2) Jones, W. P. and Launder, B. E. : "The Prediction of
k-ε-A three equation model is the most promising Laminarization with a Two-equation Model of
scheme in the eddy viscosity models of turbulence. Turbulence", Int. J. Heat Mass Transfer, 15(1972),
In the very near future, the author believes that the 301~314
3) Launder, B. E. and Spalding, D. B. : "The Numerical
standard k-ε model in industrial applications will be
Computation of Turbulent Flows", Comp. Methods
replaced with the recently developed low-Reynolds- Appl. Mech. Engng., 3(1974), 269~289
number nonlinear eddy viscosity models. 4) Wolfshtein, M. W. : "The Velocity and Temperature
To achieve the full potential of the eddy viscosity Distribution in one Dimensional Flow with
modelling, however, further attention is expected for Turbulence Augmentation and Pressure Gradient",
Int. J. Heat Mass Transfer, 12(1969), 301
optimizing the transport equation for ε. Since the rapid 5) Patankar, S. V. and Spalding, D. B. : "Heat and Mass
development of DNS will be providing much more Transfer in Boundary Layers", 2nd Edition, (1970),
detailed data for the modelling, the existing too Intertext Books, London
empirically modelled ε equations will soon be replaced. 6) Wilcox, D. C. : "Reassessment of the Scale-Determining
Equation for Advanced Turbulence Models", AIAA
J., 26-11(1988), 1299~1310
7) Patel, V. C., Rodi, W. and Scheuerer, G. : "Turbulence
Channel, Re=7000; ARG & SA Models for Near-Wall and Low-Reynolds Number
3 Flows", AIAA J., 23(1985), 1308~1319
8) Launder, B. E. and Sharma, B. I. : "Application of the
Energy-Dissipation Model of Turbulence to the
2 Calculation of Flow Near a Spinning Disc", Letters in
v’+ w’+ u’+

Heat Mass Transfer, 1(1974), 131~138

9) Rodi, W. : "The Prediction of Free Turbulent
Boundary Layers by Use of a Two-Equation Model of
1 Turbulence", Ph. D. Thesis, Univ. of London, (1972)
10) Launder, B. E., Reece, G. J. and Rodi, W. : "Progress
in the Development of a Reynolds Stress Turbulence
0
Closure", J. Fluid Mech., 68(1975), 537~566
0 20 40 60 80 100 11) Launder, B. E. : "An Introduction to Single-Point
y+ Closure Methodology, An Introduction to the
Modelling of Turbulence", Lecture Series 1987-06,
Fig. 13 Turbulent intensities: VKI, Rhode-Saint-Genese, Belgium, (1987)
Symbols, DNS ;
13) 12) Kim, J., Moin, P. and Moser, R. : "Turbulence
46) Statistics in Fully Developed Channel Flow at Low
, ARG ; , SA45).
Reynolds Number", J. Fluid Mech., 177(1987),
133~166
13) Kim, J. : Personal communication, (1989)
3 14) Spalart, P. R. : "Direct Simulation of a Turbulent
2.5 Boundary Layer up to Rθ =1410", J. Fluid Mech.,
187(1988), 61~98
2
v’+,w’+,u’+

15) Myong, H. K. and Kasagi, N. : "A New Approach to

1.5 the Improvement of k-ε Turbulence Model for Wall-
Bounded Shear Flows", JSME Int. J., Ser. II,
1
33(1990), 63~72
0.5 16) Nagano, Y. and Tagawa, M. : "An Improved k-ε
Model for Boundary Layer Flows", ASME J. Fluids
0
0 0.2 0.4 0.6 0.8 1 Engng., 112(1990), 33~39
y/δ 17) Rodi, W. and Mansour, N. : "Low Reynolds Number
k-ε Modelling with the Aid of Direct Simulation
Data", J. Fluid Mech., 250(1993), 509~529
Fig. 14 Turbulent intensities in an open channel flow: 18) Abe, K., Kondo, T. and Nagano, Y. : "A New
symbols, DNS31) ; Turbulence Model for Predicting Fluid Flow and Heat
Transfer in Separating and Reattaching Flows-I. Flow
, k-ε-A model .
29) 34)
, CLS ; Field Calculations", Int. J. Heat Mass Transfer,

R&D Review of Toyota CRDL Vol. 33 No. 1 ( 1998.3 )

 11

37-1(1994), 139~151 33) Durbin, P. A. : "Near-Wall Turbulence Closure

19) Abe, K., Kondoh, T. and Nagano, Y. : "A Low- Modelling without 'Damping Functions'", Theoret.
Reynolds Number k-ε Model Reflecting Comput. Fluid Dynamics, 3(1991), 1~13
Characteristics of a Reynolds Stress Model", 34) Suga, K. : "Nonliner Eddy Viscosity Modelling with
Trans. JSME Ser.B, 61-585(1995), 1714~1721 a Transport Equation for Lumley's Stress Flatness
20) Kawamura, H. and Kawashima, N. : "A Proposal of Parameter, Proc. 11th Symp. on Turbulent Shear
k-ε Model with Relevance to the Near-Wall Flows, (1997), 13-18~13-23
Turbulence", Proc. Int. Symp. Turbulence and Heat 35) Lumley, J. L. : "Toward a Turbulent Constitutive
Transfer, (1994) Relation", J. Fluid Mech., 41-2(1970), 413~434
21) Lee, M. J., Kim, J. and Moin, P. : "Structure of 36) Pope, S. : "A more General Effective-Viscosity
Turbulence at High Shear Rate", J. Fluid Mech., Hypothesis", J. Fluid Mech., 72(1975), 331~340
216(1990), 561~583 37) Speziale, C. G. : "On Nonlinear k-l and k-ε Models of
22) Yang, Z. and Shih, T. H. : "A Galilean and Tensorial Turbulence", J. Fluid Mech., 178(1987), 459~475
Invariant k-ε Model for Near Wall Turbulence", 38) Nisizima, S. and Yoshizawa, A. : "Turbulent Channel
NASA Technical Memorandum, 106263, (1993) and Couette Flows Using an Anisotropic k-ε Model",
23) Cotton, M. A. and Ismael, J. O. : "Development of a AIAA J., 25-3(1987), 414~420
Two-Equation Turbulence Model with Reference to a 39) Rubinstein, R. and Barton, J. M. : "Non-linear
Strain Parameter", Proc. 5th Int. Symp. Refined Flow Reynolds Stress Models and Renormalization
Modelling and Turbulence Measurements, (1993), Group", Phys. Fluids, A 2(1990), 1472~1476
117~124 40) Myong, H. K. and Kasagi, N. : "Prediction of
24) Craft, T. J., Launder, B. E. and Suga, K. : "Extending Anisotropy of the Near Wall Turbulence with an
the Applicability of Eddy Viscosity Models through Anisotropic Low-Reynolds Number k-ε Turbulence
the Use of Deformation Invariants and Non-linear Model", ASME J. Fluid Eng., 112(1990), 521~524
Elements", Proc. 5th Int. Symp. Refined Flow 41) Shih, T. H., Zhu, J. and Lumley, J. L. : "A Realizable
Modelling and Turbulence Measurements, (1993), Reynolds Stress Algebraic Equation Model", NASA
125~132 Technical Memorandum, 105993, (1993)
25) Craft, T. J., Launder, B. E. and Suga, K. : 42) Taulbee, D. B. : "An Improved Algebraic Reynolds
"Development and Application of a Cubic Eddy- Stress Model and Corresponding Nonlinear Stress
Viscosity Model of Turbulence", Int. J. Heat Fluid Model", Phys. Fluids A4-11(1992), 2555~2561
Flow, 17(1996), 108~115 43) Gatski, T. B. and Speziale, C. G. : "On Explicit
26) Cotton, M. A. and Ismael, J. O. : "A 3-equation Algebraic Stress Models for Complex Turbulent
Transported Strain Parameter Model of Turbulence Flows", J. Fluid Mech., 254(1993), 59~78
and Its Application to Turbulent Mixed Convection", 44) Wilcox, D. C. : "The Remarkable Ability of
Proc. Int. Symp. Turbulence, Heat and Mass Transfer, Turbulence Model Equations to Describe Transition",
(1994) Proc. 5th Symp. on Numerical and Physical Aspects
27) Suga, K. : "Development and Application of a Non- of Aerodynamic Flows, (1992)
linear Eddy Viscosity Model Sensitized to Stress and 45) Speziale, C. G. and Abid, R. : "Near-wall Integration
Strain Invariants, Ph. D. Thesis, UMIST, Manchester, of Reynolds Stress Turbulence Closure with No Wall
U.K., (1995) Damping", AIAA J., 33-10(1995), 1974~1977
28) Craft, T. J., Launder, B. E. and Suga, K. : "A Non- 46) Abid, R., Rumsey, C. and Gatski, T. : "Prediction of
Linear Eddy-Viscosity Model Including Sensitivity to Nonequilibrium Turbulent Flows with Explicit
Stress Anisotropy", Proc. Tenth Symp. on Turbulent Algebraic Stress Models", AIAA J., 33-11(1995),
Shear Flows, (1995), 23-19~24 2026~2031
29) Craft, T. J., Launder, B. E. and Suga, K. : "Prediction 47) Laufer, J. : "The Structure of Turbulence in Fully
of Turbulent Transitional Phenomena with a Developed Pipe Flow", NACA Rep., 1174, (1954)
Nonlinear Eddy-Viscosity Model", Int. J. Heat Fluid
Flow, 18(1997), 15~28 Kazuhiko Suga was born in 1959. He
30) Nezu, I. and Rodi, W. : "Open-Channel Flow received Ph. D. degree in mechanical
Measurements with a Laser-Doppler Anemometer", engineering in 1996 from the University
J. Hydraulic Eng., 112(1986), 335~355 of Manchester Institute of Science and
31) Lombardi, P., De Angelis, V. and Banerjee, S. : Technology. He belongs to the Heat
"Direct Numerical Simulation of Near-Interface Transfer Lab., where he works on
Turbulence in Coupled Gas-Liquid Flow", Phys. computational fluid dynamics and
Fluids, 8-6(1996), 1643~1665 turbulence modeling. He is a member
32) Kuroda, A., Kasagi, N. and Hirata, M. : "Direct of J.S.M.E. and the Heat Transfer Soc.
Numerical Simulation of Turbulent Plane Couette- of Jpn.
Poiseuille Flows: Effect of Mean Shear Rate on the
Near Wall Turbulence Structures", Turbulent Shear
Flows 9, Ed. by Durst, F., et al., (1995), 241~257,
Springer

R&D Review of Toyota CRDL Vol. 33 No. 1 ( 1998.3 )