next up [*] [*]
Next: 4.2.2 ReNormalisation Group k Up: 4.2 Turbulence Previous: 4.2 Turbulence

4.2.1 Standard k - epsilon Turbulence Model

At high Reynolds numbers the rate of dissipation of kinetic energy $\epsilon $ is equal to the viscosity multiplied by the fluctuating vorticity. An exact transport equation for the fluctuating vorticity, and thus the dissipation rate, can be derived from the Navier Stokes equation. The k - epsilon model consists of the turbulent kinetic energy equation

 \begin{displaymath}\nfour \nfour \nfour \nfour \nfour \nfour
\dxdt{k}{t} + \tex...
... - \rho \epsilon
\nomsort{gz}{$G$ }{Turbulent generation rate}
\end{displaymath} (4.2.13)

and the dissipation rate equation

 \begin{displaymath}\begin{split}
\dxdt{\epsilon}{t} + \text{div}(\rho \underlin...
...}{k} -
C_{2\epsilon} \rho \frac{\epsilon^2}{k} &
\end{split} \end{displaymath} (4.2.14)

where G represents the turbulent generation rate which is equal to

 \begin{displaymath}\begin{split}
G\;=\; 2 \left( \left[ \dxdt{u}{x} \right]^2 +...
...)^2
&+ \left( \dxdt{w}{y} + \dxdt{v}{z} \right)^2
\end{split} \end{displaymath} (4.2.15)

In the implementation of this model the Kolmogorov - Prandtl expression for the turbulent viscosity is used

 \begin{displaymath}\nu_t \; =\; C_\mu \frac{k^2}{\epsilon}
\end{displaymath} (4.2.16)

together with equation (4.2.11) for the eddy diffusivity. In the above $C_\mu$, $\sigma_k$, $\sigma_\epsilon$, $C_{1\epsilon}$ and $C_{2\epsilon}$ are all taken to be constants and are given respectively the values[96] 0.09, 1.0, 1.3, 1.44 and 1.92.

The methods used in the discretisation of the above equations are the same as those defined in previous chapters with the term on the right of each of (4.2.13) and (4.2.14) being treated as source terms. The continuity equation remains unchanged whilst the only differences in the momentum and enthalpy equations is the change in the diffusion coefficient. The laminar viscosity $\rho \nu_{lam}$ is replaced by $\rho \nu_{lam} + \rho \nu_t$ and the conductivity k in the enthalpy equation is replaced by $k + (\rho \nu_t / \sigma_t$).

Handling of the boundary conditions concerned with turbulent flows is the same as for laminar flows except in the case of walls. The k - epsilon model provides accurate solutions only for fully turbulent flows. In the part of the flow near to walls there exist regions in which the local Reynolds number of turbulence, $k^2/\nu \epsilon$, is so small that the viscous effects become more significant than the turbulent ones. In this viscous sublayer very steep gradients occur so for accurate modelling many grid points would be required in this region. Fortunately it is not necessary to discretise the k and $\epsilon $ conservation equations over this region as there exist suitable laws which relate the wall conditions to values of the dependent variables just outside the viscous sublayer. In the region just outside the sublayer[97], such that

 \begin{displaymath}\text{30} \; < \; \frac{y U_\tau} {\nu_{lam} } \; <
\; \text{100}
\nomsort{uz}{$U_\tau$ }{Resultant friction velocity}
\end{displaymath} (4.2.17)

where y is the perpendicular distance to the wall and $U_\tau$ is the resultant friction velocity, the generation of turbulent kinetic energy is balanced by its dissipation. In this region

 \begin{displaymath}\frac{\tau}{\rho k} \; = \; \frac{U^2_\tau}{k} \; = \;
\sqrt...
...all shear stress}
\nomsort{s}{$s$ }{Wall skin friction factor}
\end{displaymath} (4.2.18)

where $\tau$ is the wall shear stress. The wall skin friction factor, s, which is defined by

 \begin{displaymath}s \; = \; \frac{\tau}{\rho V^2}
\end{displaymath} (4.2.19)

in which V is the fluid speed, is determined from the formula

 \begin{displaymath}\sqrt{s} \; = \; \frac{\kappa}{ \text{ln} ( 1.01 +
E \, Re \...
...Wall roughness parameter}
\nomsort{Re}{$Re$ }{Reynolds number}
\end{displaymath} (4.2.20)

where $\kappa$, the Kaman constant, is equal to 0.435, E, the wall roughness parameter, usually being set to 9.0, which represents smooth walls, and Re, the local Reynolds number, which is defined by the formula

 \begin{displaymath}Re \; = \; \frac {Vy} {\nu_{lam}}
\end{displaymath} (4.2.21)

For a node positioned just outside the viscous sublayer the effective viscosity is set to the maximum of the laminar dynamic viscosity $\rho \nu_{lam}$ and sVy. From equation (4.2.18) the kinetic energy of turbulence at this point is given by

 \begin{displaymath}k \; = \; \frac{\tau}{\rho \sqrt{C_\mu}} \; = \;
\frac{s V^2}{\sqrt{C_\mu}}
\end{displaymath} (4.2.22)

As previously stated in this region the production of turbulent kinetic energy is equal to its dissipation. This sets

 \begin{displaymath}\epsilon \; = \; U^2_\tau \dxdt{V}{y}
\end{displaymath} (4.2.23)

The value of the gradient term is calculated by substituting (4.2.19) into equation (4.2.20) to remove s. Differentiation of the resulting equation, with respect to y, gives a formula for $\partial V/ \partial y$

\begin{displaymath}\begin{split}
\dxdt{V}{y} \; &= \; \frac{1}{\kappa} \sqrt{\f...
...}{\kappa} \sqrt{\frac{\tau}{\rho}}
\, \frac{1}{y}
\end{split} \end{displaymath} (4.2.24)

Substitution of this into equation (4.2.23) and using (4.2.18) gives

 \begin{displaymath}\epsilon \; = \; \frac{0.1643}{\kappa} \frac{k^{1.5}}{y}
\end{displaymath} (4.2.25)

The source terms in the conservation equations (4.2.13) and (4.2.14) include the quantity G which is the rate of generation of turbulent kinetic energy the formula for which is given in equation (4.2.15). The calculation of G requires the values of the derivatives of the three Cartesian velocity components with respect to the three Cartesian directions. The simplest method for the estimation of any of these is to use the following formulae

 \begin{displaymath}\int_V \dxdt{u_i}{x_j} \; \approx \; V \dxdt{u_i}{x_j} \Bigg\vert_P
\end{displaymath} (4.2.26)

but using the divergence theorem

 \begin{displaymath}\int_V \dxdt{u_i}{x_j} \; = \; \sum_F \int_F u_i n_j dS
\end{displaymath} (4.2.27)

By making the usual approximations to the right hand side of (4.2.27) and equating to (4.2.26) the following estimate for the derivative's value at P is obtained.

 \begin{displaymath}\dxdt{u_i}{x_j} \Bigg\vert_P \; \approx \; \frac{1}{V}
\sum_F A_F n_j (u_i)_F
\end{displaymath} (4.2.28)

This formula only requires the estimation of a representative value of ui on each face of the control volume and so can use the methods discussed previously, for example those included in the section on the cross product diffusion terms. The only care that is needed is in the handling of walls. At a wall there is no generation of turbulence so for near wall elements care is required in the calculation of the differentials.

The source terms, particularly those in the dissipation rate equation, are non linear. Three methods of linearisation are available within the code in line with those used by PHOENICS[60]. Assuming that the linearised form of the sources can be written in the form SC - $S_L \phi$ then the following gives the formulae for the coefficients for the various methods.

METHOD 1  
$\quad$ k source terms :- $\qquad$ $S_C = \nu_t \rho G$
  $S_L = C_D \nu_t \rho / ( C_\mu l^2)$
$\quad \epsilon$ source terms :- $\qquad$ $S_C = C_D \rho \nu_t C_{1\epsilon} G \nu_t / ( C_\mu l^2)$
  $S_L = C_{2\epsilon}C_D \rho \nu_t / (C_\mu l^2)$
METHOD 2  
$\quad$ k source terms :- $\qquad$ $S_C = (C_{2\epsilon}-1.0) C_\mu C_D k^2 \rho / \nu_t +
1.5 \rho \nu_t G$
  $S_L = 0.5 \nu_t \rho G / k +
C_{2\epsilon} C_\mu C_D k \rho / nu_t$
$\quad \epsilon$ source terms :- $\qquad$ $S_C = (C_\mu C_D \rho k / \nu_t )( C_{1\epsilon} \nu_t G +
(C_{2\epsilon}-1.0) \epsilon)$
  $S_L = (2.0C_{2\epsilon}-1.0) C_\mu C_D \rho k / \nu_t $
METHOD 3  
$\quad$ k source terms :- $\qquad$ $S_C = \nu_t \rho G$
  $S_L = C_D C_\mu k \rho / \nu_t$
$\quad \epsilon$ source terms :- $\qquad$ $S_C = C_\mu C_D \rho k C_{1\epsilon} G $
  $S_L = C_{2\epsilon} C_\mu C_D \rho k / \nu_t $
The earlier the method in the above list the quicker convergence will be obtained whilst the later methods offer greater numerical stability.


next up [*] [*]
Next: 4.2.2 ReNormalisation Group k Up: 4.2 Turbulence Previous: 4.2 Turbulence

2002-04-23