T h e o r y f o r t h e T u r b u l e n t F l o w I n t e r f a c e s

The Single-Phase Flow, Turbulent Flow Interfaces theory is described in this section:

•Turbulence Modeling

•The k- Turbulence Model

•The k- Turbulence Model

•The Low Reynolds Number k- Turbulence Model

•The Spalart-Allmaras Turbulence Model

•Inlet Values for the Turbulence Length Scale and Intensity

•Pseudo Time Stepping for Turbulent Flow Models

•References for the Single-Phase Flow, Turbulent Flow Interfaces

Turbulence Modeling

Turbulence is a property of the flow field and it is mainly characterized by a wide range of flow scales: the largest occurring scales, which depend on the geometry, the smallest, quickly fluctuating scales and all the scales there in between. The tendency for an isothermal flow to become turbulent is measured by the Reynolds number

where is the dynamic viscosity, the density, and U and L are velocity and length scales of the flow, respectively. Flows with high Reynolds numbers tend to become turbulent and this is the case for most engineering applications.

The Navier-Stokes equations can be used for turbulent flow simulations, although this would require a large number of elements, to capture the wide range of scales in the flow. An alternative approach is to divide the flow in large, resolved scales and small, unresolved scales. The small scales are then modeled using a turbulence model with the intension that the model is numerically less expensive than resolving all present scales. Different turbulence models invoke different assumptions on the unresolved scales resulting in different degree of accuracy for different flow cases.

T H E O R Y F O R T H E TU R B U L E N T F L O W I N T E R F A C E S | 143

This module includes Reynolds-averaged Navier-Stokes (RANS) models which is the model type most commonly used for industrial flow applications

R E Y N O L D S - A V E R A G E D N A V I E R - S T O K E S ( R A N S ) E Q U A T I O N S

The information below assumes that the flow fluid is incompressible and Newtonian in which case the Navier-Stokes equations take the form:

u = 0

Once the flow has become turbulent, all quantities fluctuate in time and space. It is seldom worth the extreme computational cost to obtain detailed information about the fluctuations. An averaged representation often provides sufficient information about the flow.

The Reynolds-averaged representation of turbulent flows divides the flow quantities into an averaged value and a fluctuating part,

= +

where can represent any scalar quantity of the flow. In general, the mean value can vary in space and time. This is exemplified in Figure 4-5, which shows time averaging of one component of the velocity vector for nonstationary turbulence. The unfiltered flow has a time scale t1. After a time filter with width t2 t1 has been applied, there is a fluctuating part, u i, and an average part, Ui. Because the flow field also varies on a time scale longer than t2, Ui is still time dependent but is much smoother than the unfiltered velocity ui.

144 | C H A P T E R 4 : S I N G L E - P H A S E F L O W B R A N C H

Figure 4-5: The unfiltered velocity component ui, with a time scale t1, and the averaged velocity component, Ui, with time scale t2.

Decomposition of flow fields into an averaged part and a fluctuating part, followed by insertion into the Navier-Stokes equation, then averaging, gives the Reynolds-averaged Navier-Stokes (RANS) equations:

U = 0

where U is the averaged velocity field and is the outer vector product. A comparison with Equation 4-11 indicates that the only difference is the appearance of the last term on the left-hand side of Equation 4-12. This term represents interaction between the fluctuating velocities and is called the Reynolds stress tensor. This means that to obtain the mean flow characteristics, information about the small-scale structure of the flow is needed. In this case, that information is the correlation between fluctuations in different directions.

E D D Y V I S C O S I T Y

The most common by far way to model turbulence is by assuming the turbulence to be of a purely diffusive nature. The deviating part of the Reynolds stress is then expressed by

u' u' – -- trace u' u' I = – T U + U T

3

where T is the eddy viscosity, also known as the turbulent viscosity. The spherical part can be written

T H E O R Y F O R T H E TU R B U L E N T F L O W I N T E R F A C E S | 145

where k is the turbulent kinetic energy. In simulations of incompressible flows, this term is included in the pressure, but when the absolute pressure level is of importance (in compressible flows, for example) this term must be explicitly included

T U R B U L E N T C O M P R E S S I B L E F L O W

If the Reynolds average is applied to the compressible form of the Navier-Stokes, terms of the form u appear and need to be modeled. To avoid this, a density based average, known as the Favre average, is introduced:

Using Equation 4-14 and Equation 4-15 along with some modeling assumption for compressible flows (Ref. 7), Equation 4-13 can be written in the form

The Favre-averaged Reynolds stress tensor is modeled using the same argument as for incompressible flows:

146 | C H A P T E R 4 : S I N G L E - P H A S E F L O W B R A N C H

where k is the turbulent kinetic energy. Comparing Equation 4-16 to its incompressible counterpart (Equation 4-12), it can be seen that apart from the term

– 2 3 k ij , the compressible and incompressible formulations are exactly the same,

except that the free variables are ˜ instead of = . More information on ui Ui ui

modeling of turbulent compressible flows can be found in Ref. 1 and Ref. 7.

The turbulent transport equations are used in their fully compressible formulations (Ref. 8).

The k- Turbulence Model

The k- model is one of the most used turbulence models for industrial applications. This module includes the standard k- model (Ref. 1). This model introduces two additional transport equations and two dependent variables: the turbulent kinetic energy, k, and the dissipation rate of turbulence energy, . Turbulent viscosity is modeled by

T H E O R Y F O R T H E TU R B U L E N T F L O W I N T E R F A C E S | 147

The model constants in Equation 4-17, Equation 4-18, and Equation 4-20 are determined from experimental data (Ref. 1); the values are listed in Table 4-8.

TABLE 4-8: MODEL CONSTANTS

M I X I N G L E N G T H L I M I T

Equation 4-18 and Equation 4-20 cannot be implemented directly as written. There is for example nothing that prevents division by zero. The equations are instead implemented as suggested in Ref. 9. The implementation includes an upper limit on the mixing length, lmixlim :

The mixing length is used to calculated the turbulent viscosity. lmixlim should not be active in a converged solution but is merely a tool to obtain convergence.

R E A L I Z A B I L I T Y C O N S T R A I N T S

The eddy-viscosity model of the Reynolds stress tensor can be written

2

uiuj = – 2 TSij + -3- k ij

where ij is the Kronecker delta and Sij is the strain rate tensor. The diagonal elements of the Reynolds stress tensor must be nonnegative, but calculating T from Equation 4-17 does not guarantee this. To assert that uiui 0 i , the turbulent viscosity is subjected to a realizability constraint. The constraint for 2D and 2D axisymmetry is:

and for 3D and 2D axisymmetry with swirl flow it reads:

148 | C H A P T E R 4 : S I N G L E - P H A S E F L O W B R A N C H

Combining equation Equation 4-22 with Equation 4-17 and the definition of the mixing length gives a limit on the mixing length scale:

Equivalently, combining Equation 4-23 with Equation 4-17 and Equation 4-21 gives:

Observe that this means that this means that there are two limitations on lmix: the realizability constraint and the imposed limit via Equation 4-21.

The effect of not applying a realizability constraint is typically excessive turbulence production. The effect is most clearly visible in stagnation points. To avoid such artifacts, the realizability constraint is always applied for the RANS models. More details can be found in Ref. 4, Ref. 5, and Ref. 6.

M O D E L L I M I T A T I O N S

The k- turbulence model relies on several assumptions, the most important of which are that the Reynolds number is high enough and that the turbulence is in equilibrium in boundary layers, which means that production equal dissipation. These assumptions limit the accuracy of the model because they are not always true. It does not, for example, respond correctly to flows with adverse pressure gradients which can result in under predicting the spatial extension of recirculation zones (Ref. 1). Furthermore, in the description of rotating flows, the model often shows poor agreement with experimental data (Ref. 2). In most cases, the limited accuracy is a fair trade-off for the amount of computational resources saved compared to more complicated turbulence models.

W A L L F U N C T I O N S

The flow close to a solid wall is for a turbulent flow is very different compared to the free stream. This means that the assumptions used to derive the k- model are not valid close to walls. While it is possible to modify the k- model so that it describes the flow in wall regions (see The Low Reynolds Number k- Turbulence Model), this is not always desirable because of the very high resolution requirements that follows. Instead,

T H E O R Y F O R T H E TU R B U L E N T F L O W I N T E R F A C E S | 149

analytical expressions are used to describe the flow at the walls. These expressions are known as wall functions.

The wall functions in COMSOL are such that the computational domain is assumed to start a distance w from the wall (see Figure 4-6).

Mesh cells

w

Solid wall

Figure 4-6: For wall functions, the computational domain starts a distance w from the wall.

The distance w is automatically computed so that w+ = u w , where

u C 1/4 k is the friction velocity, becomes 11.06. This corresponds to the distance from the wall where the logarithmic layer meets the viscous sublayer (or rather would meet if there was not a buffer layer in between). w is limited from below so that it never becomes smaller than half of the height of the boundary mesh cell. This means that w+ can become higher than 11.06 if the mesh is relatively coarse.

The boundary conditions for the velocities are a no-penetration condition u n = 0 and a shear stress condition

where = u + u T is the viscous stress tensor and

150 | C H A P T E R 4 : S I N G L E - P H A S E F L O W B R A N C H

where in turn v is the von Kárman constant (default value 0.41) and B is a constant that per default is set to 5.2.

The turbulent kinetic energy is subject to a homogeneous Neumann condition n k = 0 and the boundary condition for reads:

C3/ 4k3/ 2

= -----------------------

v w

See Ref. 9 and Ref. 10 for further details.

I N I T I A L V A L U E S

The default initial values for a stationary simulation are (Ref. 9),

The k- Turbulence Model

The k- model solves for the turbulent kinetic energy, k, and for the dissipation per unit turbulent kinetic energy, . is also commonly know as the specific dissipation rate. The CFD Module has the Wilcox revised k- model (Ref. 1)

T H E O R Y F O R T H E TU R B U L E N T F L O W I N T E R F A C E S | 151

where in turn ij is the mean rotation-rate tensor

Pk is given by Equation 4-19. The following auxiliary relations for the dissipation, , and the turbulent mixing length, l , are also be used:

152 | C H A P T E R 4 : S I N G L E - P H A S E F L O W B R A N C H