Numerical Investigation of Flow Past a Circular Cylinder and in a Staggered Tube Bundle Using Various Turbulence Models

(1)

LAPPEENRANTA UNIVERSITY OF TECHNOLOGY Faculty of Technology

Department of Mathematics and Physics

Numerical Investigation of Flow Past a Circular Cylinder and in a Staggered Tube

Bundle Using Various Turbulence Models

The topic of this Master’s thesis was approved by the departmental council of the De- partment of Mathematics and Physics on27^th May, 2010.

Supervisors: Professor Heikki Haario

Associate Professor Teemu Turunen-Saaresti

Examiners: Professor Heikki Haario

Associate Professor Teemu Turunen-Saaresti

In Lappeenranta August 24, 2010

Yogini Patel

Teknologiapuistonkatu 2 B 28 53850 Lappeenranta

Phone: +358466175120

Email: Yogini.Patel@lut.fi & yogi.msu@gmail.com

(2)

Abstract

Lappeenranta University of Technology Department of Mathematics and Physics

Yogini Patel

Numerical Investigation of Flow Past a Circular Cylinder and in a Staggered Tube Bundle Using Various Turbulence Models

Master’s thesis 2010

87 pages, 48 figures, 13 tables

Key words: Turbulent flow, RANS models, Large eddy simulation, Circular cylinder, Vortex shedding, Tube bundle.

Transitional flow past a three-dimensional circular cylinder is a widely studied phenomenon since this problem is of interest with respect to many technical applications.

In the present work, the numerical simulation of flow past a circular cylinder, performed by using a commercial CFD code (ANSYS Fluent 12.1) with large eddy simulation (LES) and RANS (κ−and Shear-Stress Transport (SST)κ−ωmodel) approaches. The turbulent flow forRe_D = 1000 & 3900 is simulated to investigate the force coefficient, Strouhal number, flow separation angle, pressure distribution on cylinder and the complex three dimensional vortex shedding of the cylinder wake region. The numerical results extracted from these simulations have good agreement with the experimental data (Zdravkovich, 1997). Moreover, grid refinement and time-step influence have been examined.

Numerical calculations of turbulent cross-flow in a staggered tube bundle continues to attract interest due to its importance in the engineering application as well as the fact that this complex flow represents a challenging problem for CFD. In the present work a time dependent simulation using κ−,κ−ω and SST models are performed in two dimensional for a subcritical flow through a staggered tube bundle. The predicted turbulence statistics (mean and r.m.s velocities) have good agreement with the experimental data (S. Balabani, 1996). Turbulent quantities such as turbulent kinetic energy and dissipation rate are predicted using RANS models and compared with each other. The sensitivity of grid and time-step size have been analyzed. Model constants sensitivity study have been carried out by adoptingκ−model. It has been observed that model constants are very sensitive to turbulence statistics and turbulent quantities.

(3)

Acknowledgements

In the first place I would like to record my gratitude to the Department of Mathematics from Lappeenranta University of Technology for financing of my studies.

I would like to express my deepest sense of gratitude to my supervisor, Professor Heikki Haario, for his patient guidance, encouragement, understanding and excellent advice throughout this thesis as well as providing me all facilities during my thesis.

My sincere appreciation goes to Associate Professor Teemu Turunen-Saaresti for his constant support, great advice and comments, directing and assistance during the work of this thesis, which appointed him a backbone of this work. I attribute the level of my thesis to his encouragement and effort and without him this thesis, too, would not have been completed or written.

Most important thanks here goes to my best friend, Gitesh. I don’t know what my life would be without you. You are like the sunshine, always giving me the feeling of warmth, hope and peace. It is so wonderful to have you beside me, in the past, present, and future. So Thank you so much for your faithful love and endless help. I could say, without you, this thesis wouldn’t exist.

I am also grateful to all my friends in Lappeenranta who helped me to have enjoyable and memorial stay.

Where would I be without my family? None of this would have been possible without the love and patience of my family. My parents deserve special mention for their inseparable support and prayers. My dad, Natvarlal, in the first place is the person who put the fundament my learning character, showing me the joy of intellectual pursuit ever since I was a child. My mom, Kantaben, is the one who sincerely raised me with her caring and unconditional love. I convey special acknowledgement to my elder brother, Alpeshbhai, and his wife, Deepali, for their support and love. I would like to express the dearest thanks to my niece, Aastha, for her loving support. I am extending my sincere gratitude to my kind grandparents for their blessings. I dedicate this thesis to my family and Gitesh, the most special person in my life.

I am ever grateful to Lord Shiva, the Creator and the Guardian, and to whom I owe my very existence.

Lappeenranta, August 24, 2010 Yogini Patel

(4)

List of Tables

1 Simulation settings for flow past a circular cylinder case with RANS models 39 2 Simulation settings for flow past a circular cylinder case with LES model . 39

3 Experimental and computational results forRe_D = 1000 . . . 41

4 Experimental and computational results forRe_D = 3900 . . . 45

5 Details of grids used in mesh-independence tests . . . 50

6 Drag coefficient of the flow past a circular cylinder using the RANS models 51 7 Computed flow parameters in comparison with experimental results using LES . . . 51

8 CPU times details used by each model for different grids. . . 55

9 Effect of time-step size on theC_dand S_t using LES. . . 56

10 Simulation settings of flow past in a staggered tube bundle case . . . 62

11 Details of grids used in mesh-independence tests and theiry⁺ values . . . 68

12 Averaged difference in the prediction of turbulent kinetic energy and dissipation rate with respect to originalC₁ . . . 79

13 Averaged difference in the prediction of turbulent kinetic energy and dissipation rate with respect to originalC₂ . . . 82

(7)

List of Figures

1 (a) Fluid element for conservation laws (b) mass flows in and out of fluid

element . . . 4

2 (a) Stress components on three faces of fluid element (b) stress components in the x-direction . . . 6

3 Overview of the CFD . . . 10

4 Turbulence models classification . . . 14

5 Extend of modelling for certain types of turbulent models [21] . . . 18

6 Subdivisions of the Near-Wall Region [8] . . . 31

7 Vortex shedding in the wake region of the flow past a circular cylinder [31]. 34 8 Diagram of forces acting around a circular cylinder. . . 35

9 Computational geometry and boundary conditions. . . 37

10 (a) Computational domain, (b) grid around the cylinder, and (c) 3D closer view of the cylinder. . . 38

11 Time histories of (a) drag coefficient,C_d, and (b) lift coefficient,C_l, for LES. 41 12 Pressure coefficient distribution around circular cylinder atRe_D = 1000. . 42

13 Iso-surfaces of x-vorticity produced from ske and SST models withRe_D = 1000. . . 43

14 Iso-surfaces of instantaneous x-vorticity withRe_D = 1000. . . 43

15 Contours of velocity magnitude and velocity vector withRe_D = 1000. . . 44

16 C_l history for κ− and SST withRe_D = 3900. . . 45

17 C_d and C_l histories for LES withRe_D = 3900. . . 46

18 Drag coefficient of the flow past a circular cylinder compared to experimental data. . . 47

19 Pressure coefficient distribution around circular cylinder atReD = 3900. . 47

20 Iso-surfaces of x-vorticity produced from ske and SST models withReD = 3900. . . 48

(8)

21 Iso-surfaces of instantaneous x-vorticity withRe_D = 3900. . . 49 22 C_d and C_l histories with different grids for LES. . . 52 23 Iso-surface of the magnitude of instantaneous vorticity for LES with (a)

Grid A, (b) Grid B, and (c) Grid C. . . 53 24 Contour of x-velocity for LES with (a) Grid A, (b) Grid B, and (c) Grid C. 54 25 C_d and C_l histories with different timestep size. . . 56 26 Strouhal number of the flow past a circular cylinder compared to experi-

mental data. . . 57 27 Iso-surface of instantaneous x-vorticity for LES using different time-step. . 57 28 (a) Cross-sectional view of the tube bundle, and (b) locations at which

results are presented. . . 60 29 (a) Boundary conditions, (b) computational grid, and (C) closer view of

the tube surface. . . 61 30 Contours of U-component of velocity (a-c) and vorticity magnitude (d-f). . 63 31 Contours of velocity vectors of SST model. . . 64 32 Comparison between profiles of predicted streamwise mean velocity with

experimental values at selected axial locations. . . 65 33 Comparison between profiles of predicted transverse mean velocity with

experimental values at selected axial locations. . . 66 34 Comparison between profiles of predicted streamwise fluctuations of the

velocity component with experimental values at selected axial locations. . 67 35 Comparison the predicted dissipation rate at selectedx/d locations. . . 68 36 Contours of turbulent kinetic energy with (a) Grid A, (b) Grid B, (c) Grid

C, and (d) Grid D using SST model. . . 69 37 Contours of U-component of velocity with (a) Grid A, (b) Grid B, (c)

Grid C, and (d) Grid D using SST model. . . 69 38 Profiles of grid-independence test of streamwise mean velocity at x/d =

1.25 and 5.45. . . 71

(9)

39 Profiles of grid-independence test of transverse mean velocity at x/d = 1.25 and 5.45. . . 72 40 Profiles of grid-independence test of turbulent kinetic energy atx/d= 1.25

and5.45. . . 73 41 Time-step effect on streamwise mean velocity atx/d= 3.35 and5.45. . . . 74 42 Time-step effect on transverse mean velocity atx/d= 3.35 and5.45. . . . 75 43 Comparison between profiles of predicted streamwise (a-d) and transverse

(e-f) mean velocities with experimental values at selected axial locations. . 77 44 Comparison the predicted turbulent kinetic energy at selectedx/dlocations. 78 45 Comparison the predicted dissipation rate at selectedx/d locations. . . 78 46 Comparison between profiles of predicted streamwise (a-d) and transverse

(e-f) mean velocity with experimental values at selected axial locations. . 80 47 Comparison the predicted turbulent kinetic energy at selectedx/dlocations. 81 48 Comparison the predicted dissipation rate at selectedx/d locations. . . 81

(10)

Nomenclature

A projected area [m²]

B additive constant C1, C2 model constants C_d drag coefficient C_l lift coefficient C_p pressure coefficient CSGS sub-grid scale constant

d tube diameter [m]

D diameter of the cylinder [m]

f force [N]

fs Strouhal frequency [Hz]

k von karamn’s constant N number of nodes

p pressure [N m⁻²]

Re Reynolds number

S source term

S_L longitudinal pitch-to-diameter ratio St Strouhal number

ST transverse pitch-to-diameter ratio

t time [s]

u velocity [ms⁻¹]

u⁰ fluctuation velocity [ms⁻¹]

(11)

U x-direction mean velocity component [ms⁻¹]

U∞ approach velocity [ms⁻¹]

V y-direction mean velocity component [ms⁻¹]

x streamwise coordinate [m]

y transverse coordinate [m]

y⁺ non-dimensional normal distance from wall

Greek Letters

∆ filter cutoff width [m]

turbulent dissipation rate [m²s⁻³]

ε_k dissipation term in the turbulent kinetic energy budget [m²s⁻³]

θ separation angle [degree]

` turbulent length scale [m]

κ turbulent kinetic energy [m²s⁻²]

λ viscosity [m²s⁻¹]

µ dynamic viscosity [N sm⁻²]

µt turbulent viscosity [N sm⁻²]

ρ density [kgm⁻³]

σ Prandtl number

τ viscous stress [N m⁻²]

ϑ velocity scale [ms⁻¹]

φ filtered function

ω specific dissipation [s⁻¹]

(12)

Abbreviations

CFD Computational Fluid Dynamics DES Detached Eddy Simulation DNS Direct Numerical Simulation DSM Dynamic Smagorinsky Model EVM Eddy-viscosity Models FFT Fast Fourier Transform FVM Finite Volume Method LDA Laser Doppler Anemometry LES Large Eddy Simulation PDE Partial Differential Equation RANS Reynolds-averaged Navier-Stokes RSM Reynolds Stress Model

SGS Sub-grid Scale Stresses ske standard κ−

SST Shear-Stress Transport TKE Turbulent Kinetic Energy

Subscript or Superscript

d, l referring drag and lift respectively i, j referring i- and j-directions respectively

x, y, z referring to x-, y- and z-directions respectively

(13)

1 Introduction

The topic of this thesis is the numerical investigation of flow past a circular cylinder and in a staggered tube bundle using various turbulence models. Computational Fluid Dynamics (CFD) calculates numerical solutions using the equations governing fluid flow.

One of the classical problems in fluid mechanics is the determination of the flow field past a bluff body represented by a circular or rectangular cylinder. The flow past circular cylinders has been extensively studied due to its importance in many practical applications, such as heat exchangers, chimneys, hydrodynamic loading on ocean marine piles and offshore platform risers and support legs [12]. In scientific terms, the flow around circular cylinders includes a variety of fluid dynamics phenomena, such as separation, vortex shedding and the transition to turbulence. The mechanisms of vortex shedding and its suppression have significant effects on the various fluid-mechanical properties of practical interest: flow-induced forces such as drag and lift forces and pressure coefficient.

Furthermore, the predictions of turbulent cross-flow in a staggered tube bundle continues to attract interest due to its importance in the engineering application as well as the fact that this complex flow represents a challenging problem for CFD. Cross-flow in tube bundles has wide practical applications in the design of heat exchangers, in flow across overhead cables, and in cooling systems for nuclear power plants [13]. In addition to the complexity arising from the flow instabilities in the tube bundle, one must also consider whether the flow is turbulent or laminar. Flows in tube bundles are usually subcritical (mixed, transition to turbulence occurs after separation) or critical (predominantly turbulent, only part of the boundary layer developing on the tube surface is laminar). In critical flows, transition to turbulence occurs before separation and turbulence is promi- nent in the rest of the boundary layer and in the flow inside the bundle (Zukauskas, 1989).

The combination of the flow instabilities and the transitional phenomena present in the boundary layers makes this type of flow difficult to model numerically. To appropriately simulate these flow characteristics, various turbulence models in CFD are used.

1.1 Objectives of the thesis

The objectives of the thesis have been set as follows:

• To investigate the flow past a circular cylinder, which is sensitive to changes of Reynolds number. Two Reynolds numbers (1000 and 3900) have been tested using unsteady turbulence models such as the RANS turbulence models namely, the κ− model and the Shear-Stress Transport (SST) κ−ω model, and the Large Eddy Simulation (LES). The study includes the simulation of vortex shedding phenomenon, force coefficient, Strouhal number and pressure distribution of the

(14)

flow. The numerical results extracted from CFD simulations are compared with each other (κ−model, SST model and LES model) and with the experimental data in order to determine the relative performance of these turbulence models and to find the best model for the flow of interest. Furthermore, grid-independence tests and time integration are performed to investigate the influence of grid refinement and time-step size effect respectively on the solution.

• In a staggered tube bundle case, time dependent calculations of the subcritical cross flow through a tube bundle have been performed in two dimensions with RANS models. RANS turbulence models like, the κ−model, the κ−ω model and the SST model have been chosen to test the suitability and the applicability of the models on the flow past a staggered tube bundle. For models comparison pur- pose, the streamwise and transverse, mean and r.m.s velocities are compared with experimental data at different location. Turbulent quantities such as turbulent kinetic energy and dissipation rate are calculated using RANS models and compared with each other. Moreover, grid and time-step size sensitivity are performed using above mentioned turbulence models.

1.2 Thesis structure

This chapter consists of the objectives and the methodology of the thesis work. The reader is then introduced to the main content of the following chapters.

Chapter 2, discusses the theoretical background of basic equations describing fluid motion. It explains how CFD formulates those equations. By using those equations, the Navier-Stokes equations are derived. In the middle part of this chapter discusses the principles of the CFD. The last part describes the definition of turbulence and its rich history.

Chapter 3 gives the brief details of turbulence models. The classification of the models are discussed based on filtering and time averaging. Different turbulence models such as the large eddy simulation (LES) and the Reynolds Averaged Navier-Stokes (RANS) models are explained with the suitability of each model in the applications of the flow past a circular cylinder and in a staggered tube bundle. The end part of this chapter contains the information about law of the wall.

Chapter 4 presents work done on the flow past a circular cylinder using RANS and LES methods at two different Reynolds numbers. The chapter focuses mainly on the verification and validation of RANS models and LES on the flow past a circular cylinder.

Pressure distribution, as well as the comparison of the Strouhal number and the drag coefficient of the flow from the prediction of RANS models and LES are compared to

(15)

experimental data. Grid and time-step size effect cases are discussed at the end also.

Chapter 5 contains work done on the flow past in a staggered tube bundle. The beginning part of this chapter displays the information about the staggered tube bundle.

The chapter focuses on the applicability and suitability of RANS models with present application and compared with experimental results. At the end, results obtained from CFD simulations for grid and time-step size influences are studied.

Chapter 6 draws conclusions of the whole thesis. This focuses on the objectives of the work and how they are achieved throughout the thesis.

(16)

2 Theoretical background

This chapter is a review of general theory of the governing equations for fluid flow. The governing equations of fluid flow are called the Navier-Stokes equations. In this section, concisely we will discuss the principles of the CFD with its components. Moreover, fundamental description of turbulence and its history will demonstrate.

2.1 Governing equations of fluid flow

In mid 18th century, the French engineer Claude Navier and the Irish mathematician George Stokes derived the well-known equations of fluid motion, known as the Navier- Stokes equations. These equations have been derived based on the fundamental governing equations of fluid dynamics, called the continuity, the momentum and the energy equations, which represent the conservation laws of physics [15].

2.1.1 Continuity equation

Figure 1: (a) Fluid element for conservation laws (b) mass flows in and out of fluid element

In Figure 1a, the six faces are labeled N, S, E, W, T and B, which stands for North, South, East, West, Top and Bottom. The centre of the element is located at position (x, y, z). The derivation of the mass conservation equation is to write down a mass balance for the element:

Rate of increase of mass in fluid element = Net rate of flow of mass into fluid element [1]

The rate of increase of mass in the fluid element is

(17)

∂(ρδxδyδz)

∂t = ∂ρ

∂tδxδyδz (1)

Next we need to account flow rate across a face of the element, which is given by the product of density, area and the velocity component normal to the face. From Figure 1b it can be seen that the net rate of flow of mass into the element across its boundaries is given by

ρu−∂(ρu)

∂x 1 2δx

δyδz−

ρu+∂(ρu)

∂x 1 2δx

δyδz

+

ρv−∂(ρv)

∂y 1 2δy

δxδz−

ρv+∂(ρv)

∂y 1 2δy

δxδz

+

ρw−∂(ρw)

∂z 1 2δz

δxδy−

ρw+ ∂(ρw)

∂w 1 2δz

δxδz (2)

The rate of increase of mass inside the element from equation (1) is now equated to the net rate of flow of into the element across its faces from equation (2). All terms of the resulting mass balance are arranged on the left hand side of the equals sign and the expression is divided by the element volumeδxδyδz [2]. This yields

∂ρ

∂t +∂(ρu)

∂x +∂(ρv)

∂y +∂(ρw)

∂z = 0 (3)

Or in vector notation

∂ρ

∂t + div(ρu) = 0 (4)

Equation (4) is the unsteady, three-dimensional mass conservation or continuity equation at a point in a compressible fluid [1].

For an incompressible fluid the densityρ is constant and equation (4) becomes

divu= 0 (5)

2.1.2 Momentum equation

Newton’s second law states that [1]:

Rate of increase of momentum of fluid particle = Sum of forces on fluid particle

(18)

The rates of increase of x-, y- and z-momentum per unit volume of a fluid particle are given by

ρDu

Dt ρDv

Dt ρDw

Dt (6)

Now, the state of stress of a fluid element is defined in terms of the pressure and the nine viscous stress components which is shown in Figure 2a.

Figure 2: (a) Stress components on three faces of fluid element (b) stress components in the x-direction

The pressure, a normal stress, is denoted by p. Viscous stresses are denoted by τ. The suffix notationτ_ij is applied to indicate the direction of the viscous stresses. The sufficesi andj inτ_ij indicate that the stress component acts in thej- direction on a surface normal to the i-direction. First we consider the x-components of the forces due to pressurep and stress componentsτxx,τyx andτzxwhich is shown in Figure 2b. Forces aligned with the direction of a co-ordinate axis get a positive sign and those in the opposite direction a negative sign. The net force in the x-direction is the sum of the force components acting in that direction on the fluid element.

On the pair of faces (E, W) we have

p− ∂p

∂x 1 2δx

−

τ_xx−∂τ_xx

∂x 1 2δx

δyδz

+

−

p+ ∂p

∂x 1 2δx

+

τxx+∂τxx

∂x 1 2δx

δyδz=

−∂p

∂x+ ∂τxx

∂x

δxδyδz (7) The net force in the x-direction on the pair of faces (N, S) is

(19)

−

τ_yx−∂τyx

∂y 1 2δy

δxδz+

τ_yx+∂τyx

∂y 1 2δy

δxδz= ∂τyx

∂y δxδyδz (8) And the net force in thex-direction on faces T and B is given by

−

τ_zx−∂τzx

∂z 1 2δz

δxδy+

τ_zx+∂τzx

∂z 1 2δz

δxδy= ∂τzx

∂z δxδyδz (9) The total force per unit volume on the fluid due to these surface stresses is equal to the sum of equations (7), (8) and (9) divided by the volumeδxδyδz.

∂(−p+τ_xx)

∂x +∂τ_yx

∂y +∂τ_zx

∂z (10)

In addition, the body forces are not consider in the above explanation. In further detail their overall effect can be included by defining a source S_{M x} of x-momentum per unit volume per unit time. Thex-component of the momentum equation is found by setting the rate of change ofx-momentum of the fluid particle in equation(6) equal to the total force in the x-direction on the element due to surface stresses in equation(10) plus the rate of increase ofx-momentum due to sources

ρDu

Dt = ∂(−p+τxx)

∂x +∂τ_yx

∂y + ∂τzx

∂z +SM x (11)

Similarly we can verify they-component of the momentum equation is given by

ρDv

Dt = ∂τ_xy

∂x +∂(−p+τ_yy)

∂y +∂τ_zy

∂z +S_{M y} (12)

And thez-component of the momentum equation by

ρDw

Dt = ∂τ_xz

∂x +∂τ_yz

∂y +∂(−p+τ_zz)

∂z +S_{M z} (13)

The source terms S_{M x}, S_{M y} and S_{M z} in above equations include contributions due to body forces only.

2.1.3 Navier-Stokes equations for a Newtonian fluid

The most useful forms of the conservation equation for fluid flows are obtained by introducing a suitable model for the viscous stresses τ_ij. In many fluid flows the viscous

(20)

stresses can be expressed as functions of the local deformation rate or strain rate [1]. In three dimensional flows, the local rate of deformation is composed of the linear deformation rate and the volumetric deformation rate. The rate of linear deformation of a fluid element has nine components in three dimensions, six of which are independent in isotropic fluid. they are denoted by the symbolsij.

In a Newtonian fluid the viscous stresses are proportional to the rates of deformation [2].

The three dimensional form of Newton’s law of viscosity for compressible flows involves two constants of proportionality: 1) dynamic viscosity, µ, to relate stresses to linear deformation, and 2) viscosity, λ, to relate stresses to the volumetric deformation. The viscous stress components, of which six are independent, are

τ_xx= 2µ∂u

∂x +λ divu τ_yy= 2µ∂v

∂y +λ divu τ_zz = 2µ∂w

∂z +λ divu

τxy =τyx=µ ∂u

∂y +∂v

∂x

τxz =τzx=µ ∂u

∂z +∂w

∂x

τyz =τzy =µ ∂v

∂z +∂w

∂y

(14) The second viscosity λ is not known much because of its effect is small in practice.

Substitution of the above shear stresses equation(14) into equations(11), (12) and (13) yields the so-called Navier-Stokes equations.

ρDu

Dt =−∂p

∂x+ ∂

∂x

2µ∂u

∂x +λ divu

+ ∂

∂y

µ ∂u

∂y +∂v

∂x

+ ∂

∂z

µ ∂u

∂z +∂w

∂x

+S_{M x} (15)

ρDv

Dt =−∂p

∂y + ∂

∂x

µ ∂u

∂y +∂v

∂x

+ ∂

∂y

2µ∂v

∂y +λ divu

+ ∂

∂z

µ ∂v

∂z +∂w

∂y

+S_{M y} (16)

ρDw

Dt =−∂p

∂z + ∂

∂x

µ ∂u

∂z +∂w

∂x

+ ∂

∂y

µ ∂v

∂z +∂w

∂y

+ ∂

∂z

2µ∂w

∂z +λ divu

+S_{M z} (17)

(21)

Now, rearrange the viscous stress terms as follows:

∂

∂x

2µ∂u

∂x+λ divu

+ ∂

∂y

µ ∂u

∂y +∂v

∂x

+ ∂

∂z

µ ∂u

∂z +∂w

∂x

= div(µgradu) + [sM x]

the viscous stresses in they- andz-component equations can be rearrange in a similar manner. And defining a new source terms by

SN =SM + [sM] (18)

Finally the Navier-Stokes equations can be written in the most useful form is ρDu

Dt =−∂p

∂x+ div(µgradu) +SN x (19) ρDv

Dt =−∂p

∂y + div(µgradv) +S_{N y} (20) ρDw

Dt =−∂p

∂z + div(µgradw) +SN z (21) Here the source termsSN x,SN y and SN z in above equations include contributions due to body forces. By solving these equations, the pressure and velocity of the fluid can be predicted throughout the flow.

2.2 Overview of Computational Fluid Dynamics(CFD)

Fluid dynamics is the science of fluid motion. The study of the fluid flow can be possible in three various ways as 1) Experimental 2) Theoretical and 3) Numerically. The numerical approach is called Computational fluid dynamics. CFD uses numerical methods and algorithms to solve and analyze problems that involve fluid flows by using computers [14]. The working principle of CFD based on three elements as the pre-processor, the solver and the post-processor.

• Pre-processor: Pre-processor consists of the input of the flow problem to a CFD program by means of an operator friendly interface and the subsequent transfor- mation of this input into a form suitable for use by the solver. The region of fluid to be analyzed is called the computational domain and it is made up of a number of discrete elements called the mesh (or grid). After the mesh generation, to define the properties of fluid and to specify appropriate boundary conditions [1].

• Solver: Solver calculates the solution of the CFD problem by solving the governing equations. The equations governing the fluid motion are Partial Differential Equations(PDE), made up of combinations of the flow variables (e.g. velocity and

(22)

pressure) and the derivatives of these variables. Computers cannot directly produce a solution of it. Hence the PDEs must be transformed into algebraic equations [4]. This process is known as numerical discretisation. There are four methods for it as 1) Finite difference method 2) Finite element method and 3) Finite volume method and 4) Spectral method. The finite difference method and the finite volume method both produce solutions to the numerical equations at a given point based on the values of neighboring points, whereas the finite element method produces equations for each element independently of all other elements. In the present work we have used ANSYS FLUENT 12.1 which is based on finite volume method.

• Post-processor: It used to visualize and quantitatively process the results from the solver part [1]. In a CFD package, the analyzed flow phenomena can be presented in vector plots or contour plots to display the trends of velocity, pressure, kinetic energy and other properties of the flow.

The following figure shows the schematic view of the CFD.

Figure 3: Overview of the CFD

(23)

2.3 What is turbulence ?

Turbulence is a phenomenon of fluid flow that occurs when momentum effects dominate viscous effects (high Reynolds number) [20]. It is usually triggered by some kind of dis- turbance, like flow around an object. Turbulence is characterized by random fluctuating motion of the fluid masses in three dimensions and is characterized by randomly fluctuating velocity fields at many distinct length and time scales. The fluctuating velocity fields manifest themselves as eddies (or regions of swirling motion)[4]. The free surface flows occurring in nature is almost always turbulent.

Turbulent flow is irregular, random and chaotic. The flow consists of a spectrum of different scales (eddy sizes) where largest eddies are of the order of the flow geometry [16].

At the other end of the spectra we have the smallest eddies which are by viscous forces dissipated into internal energy. Turbulent flow is dissipative, which means that kinetic energy in the small (dissipative) eddies are transformed into internal energy. The small eddies receive the kinetic energy from slightly larger eddies. The slightly larger eddies receive their energy from even larger eddies and so on. The largest eddies extract their energy from the mean flow. This process of transferred energy from the largest turbulent scales (eddies) to the smallest is called cascade process[9].

In laminar flows, viscous effects dominate momentum effects. The fluid can be thought of as flowing in layers, all of which are parallel to each other. Simple laminar flows often times permit analytic solutions to the Navier-Stokes equations. The presence of turbulence introduces many difficulties in obtaining a solution because of its inherently wide range of length and time scales.

2.4 History of turbulence

The historical overview of the study of turbulence, beginning with Leonardo da Vinci in the fifteenth Century. The first turbulence modelling may be traced back to his drawings.

But there seems to have been no substantial progress in understanding until the late19^th Century, beginning with Boussinesq in the year 1877 [6]. He introduced the idea of an eddy viscosity in addition to molecular viscosity. His hypothesis that ’turbulent stresses are linearly proportional to mean strain rates’ is still the cornerstone of most turbulence models. In 1894 Osborne Reynolds’ experiments, briefly described above and his semi- nal paper of 1894 are among the most influential results over produced on the subject of turbulence. In addition, it is interesting to note that at approximately the same time as Reynolds was proposing a random description of turbulent flow, Poincar´ewas finding

(24)

that relatively simple nonlinear dynamical systems were capable of exhibiting chaotic random-in-appearance behavior that was. in fact, completely deterministic [5].

Following Reynolds’ introduction of the random view of turbulence and proposed use of statistics to describe turbulent flows, essentially all analysis were along these lines. The first major result was obtained by Prandtl in 1925 in the form of a prediction of the eddy viscosity (introduced by Boussinesq) and the idea of a mixing length for determining the eddy viscosity. The next major steps in the analysis of turbulence were taken by G. I.

Taylor during the 1930. The literature says that he was the first researcher to utilize a more advanced level of mathematical rigor, and he introduced formal statistical methods involving correlations, Fourier transforms and power spectra into the turbulence is a random phenomenon and then proceeds to introduce statistical tools for the analysis of homogeneous isotropic turbulence[5].

In 1941 the Russian statistician A. N. Kolmogorov published three papers that provide some of the most important and most often quoted results published by Kolmogorov in a series of papers in 1941. The K41 theory provides two specific, testable results: the 2/3 law which leads directly to the prediction of aK^−5/3decay rate in the inertial range of the energy spectrum, and the 4/5 law that is the only exact results for N.-S. turbulence at high Re. Kolmogorov scale is another name for dissipation scales [4]. These scales were predicted on the basis of dimensional analysis as part of the K41 theory. In addition, in 1942 Kolmogorov developed thek−ω concept which provides the turbulent length scale,k^1/2/ω where1/ω is the turbulent time scale. In 1945 Prandtl theorized an eddy viscosity which is dependent on turbulent kinetic energy.

The first full-length books on turbulence theory began to appear in the 1950s. The best known of these are due to Batchelor, Townsend and Hinze. All of these treat only the statistical theory and heavily rely on earlier ideas of Prandtl, Taylor, Von K´arm´an. but often intermixed with the somewhat different views Kolmogorov, Obukhov and Landau.

A number of new techniques were introduced beginning in the late 1950s with the work of Kraichnan who utilized mathematical methods from quantum field theory in the analysis of turbulence. In 1963 the MIT meteorologist E. Lorenz published a paper, based mainly on machine computation that would eventually lead to a different way to view turbulence. In particular, this work presented a deterministic solution to a simple model of the Navier-Stokes equations [5].

Two other aspects of turbulence experimentation in the 70s and 80s are significant. The

(25)

first of these was detailed testing of the Kolmogorov ideas, the outcome of which was general confirmation, but not in complete detail. The second aspect of experimentation during this period involved increasingly more studies of flows exhibiting complex behav- iors beyond the isotropic turbulence [5]. By the beginning of the 1970s attention began to focus on more practical flows such as wall-bounded shear flows (especially boundary- layer transition), flow over and behind cylinders and spheres, jets, plumes, etc. During this period results such as those of Blackwelder and Kovasznay, Antonia et al., Reynolds and Hussain and the work of Bradshaw and coworkers are well known.

From the standpoint of present-day turbulence investigations probably the most important advances of the 1970s and 80s were the computational techniques. The first of these was large-eddy simulation (LES) as proposed by Deardorff in 1970. This was rapidly followed by the first direct numerical simulation (DNS) by Orszag and Patterson in 1972, and introduction of a wide range of Reynolds-averaged Navier–Stokes (RANS) approaches also beginning around 1972 (see e.g., Launder and Spalding and Launder et al.). It was immediately clear that DNS was not feasible for practical engineering problems (and probably will not be for at least another 10 to 20 years beyond the present), and in the 70s and 80s this was true as well for LES [5]. The reviews by Ferziger and Reynolds emphasize this. Thus, great emphasis was placed on the RANS approaches despite their many obvious shortcomings that we will note in the sequel. But by the beginning of the 1990s computing power was reaching a level to allow consideration of using LES for some practical problems if they involved sufficiently simple geometry, and since then a tremendous amount of research has been devoted to this technique.

Indeed, many new approaches are being explored, especially for construction of the required subgrid-scale models. These include the dynamic models of Germano et al. and Piomelli [5]. By far the most extensive work on two-equation models has been done by Launder and Spalding (1972). Launder’s k−εmodel is as well known as the mixing- length model and is the most widely used two-equation models. In 1974, Launder and Sharma was improve thek−εmodel and so called standardk−εmodel. In 1970 Saffman formulated ak−ω model without any prior knowledge of Kolmogorov’s work and that enjoys advantages over thek−ωmodel, especially for integrating through the viscous sub layer and for predicting effects of adverse pressure gradient. Wilcox and Alber (1972), Saffman and Wilcox (1974), Wilcox and Traci (1976), Wilcox and Rubesin (1980) and Wilcox (1988a) have pursed further development and application of k−ω models [4].

(26)

3 Turbulence models

Nowadays turbulent flows may be computed using several different approaches. Either by solving the Reynolds-averaged Navier-Stokes equations with suitable models for turbulent quantities or by computing them directly. The main approaches are summarized below.

3.1 Classification of turbulence models

Turbulent flows are characterized by velocity fields which fluctuate rapidly both in space and time. Since these fluctuations occur over several orders of magnitude it is compu- tationally very expensive to construct a grid which directly simulates both the small scale and high frequency fluctuations for problems of practical engineering significance.

Two methods can be used to eliminate the need to resolve these small scales and high frequencies: Filtering and Time averaging [10].

Figure 4: Turbulence models classification

Above Figure 4 presents the overview of turbulence models commonly available in CFD.

Generally, simulations of flow can be done by filtering or averaging the Navier-Stokes equations.

Filtering

The main idea behind this approach is to filter the time-dependent Navier-Stokes equation in either Fourier space or configuration space. A simulation using this approach is

(27)

known as aLarge Eddy Simulation(LES).

The filtering process creates additional unknown terms which must be modeled in order to provide closure to the set of equations. These terms are the sub-grid scale stresses and several models for these stresses. The simplest of these is the model originally proposed by Smagorinsky in which the sub-grid scale stresses (SGS) are computed using an isotropic eddy viscosity approach. The eddy viscosity is then calculated from an algebraic expression involving the product of a model constantCS, the modulus of the rate of strain tensor, and an expression involving the filter width. The problem with this approach is that there is no single value of the constant C_S which is universally applicable to a wide range of flows [10]. In addition, in the Dynamic Smagorinsky Model (DSM), the C_S is dynamically computed during the simulation using the information provided by the smaller scales of the resolved fields. C_S determined in this way varies with time and space and this allows the Smagorinsky model to cope with transitional flows and to include near-wall damping effects in a natural manner. In the next section we will get more details about Large Eddy Simulation.

Time averaging

In the Time averaging or Reynolds averaging approach all flow variables are divided into a mean component and a rapidly fluctuating component and then all equations are time averaged to remove the rapidly fluctuating components. In the Navier-Stokes equation the time averaging introduces new terms which involve mean values of products of rapidly varying quantities. These new terms are known as the Reynolds Stresses, and solution of the equations initially involves the construction of suitable models to represent these Reynolds Stresses [4]. There are two sub categories for time averaging approach:

Eddy-viscosity models (EVM) and Reynolds stress models.

Eddy-viscosity models

One assumes that the turbulent stress is proportional to the mean rate of strain. Further more eddy viscosity is derived from turbulent transport equations (usually k + one other quantity).

• Zero equation model:- The mixing length model is a zero equation models based on Reynolds averaged Navier-Stokes equations. It is one of the oldest turbulence model which was developed in the beginning of the this century. we assume the kinematic turbulent viscosityν_t, which can be expressed as a product of a turbulent velocity scaleϑand a turbulent length scale `[18].

(28)

ν_t=Cϑ` (22) Where C is a dimensionless constant of proportionality. And the dynamics turbulent viscosity is given by

µ_t=Cρϑ`

The kinetic energy of turbulence is contained in the largest eddies and turbulence length scale`. For such flows it is correct to state that, if the eddy length scale is

`,

ϑ=c `

∂U

∂y

(23) Where c is a dimensionless constant and ^∂U_∂y is the mean velocity gradient. Com- bining equations (22) and (23) and absorbing the two constants C and c into a new length scale`m we obtain

ν_t=`²_m

∂U

∂y

(24) This is Prandtl’s mixing length model. This model easy to implement and cheap in terms of computing resources. And also it is good to predict thin shear layers like jets, mixing layers, wakes and boundary layers. The mixing length model is completely incapable of describing flows with separation and recirculation. it is only calculates mean flow properties and turbulent shear stress.

• One equation models:- The Spalart-Allmaras model is one equation turbulence models because its solve a single transport equation that determines the turbulent viscosity. This is in contrast to many of the early one-equation models that solve an equation for the transport of turbulent kinetic energy and required an algebraic prescription of a length scale. The Spalart-Allmaras model also allows for reasonably accurate predictions of turbulent flows with adverse pressure gradients.

Furthermore, it is capable of smooth transition from laminar to turbulent flow at user specified locations. The Spalart-Allmaras model is an empirical equation that models production, transport, diffusion and destruction of the turbulent viscosity [3]. The Spalart-Allmaras model is suitable for aerospace applications involving wall-bounded flows and in the turbomachinery applications. In complex geome- tries it is difficult to define the length scale, so the model is unsuitable for more general internal flows.

• Two equation models :- Two equation turbulence models are one of the most common type of turbulence models. Models like theκ−model and theκ−ωmodel have become industry standard models and are commonly used for most types of engineering problems. By definition, two equation models include two extra

(29)

transport equations to represent the turbulent properties of the flow. One of the transported variables is the turbulent kinetic energy,κ, and the second transport variable varies depending on what type of two-equation model it is. Common choices are the turbulent dissipation,, or the specific dissipation, ω [4]. We will discuss in more detail later.

Reynolds stress models

The Reynolds stress model (RSM) is the most elaborate type of turbulence model. The RSM closes the Reynolds-averaged Navier-Stokes equations by solving transport equations for the Reynolds stresses, together with an equation for the dissipation rate. This means that five additional transport equations are required in 2D flows, in comparison to seven additional transport equations solved in 3D. Since the RSM accounts for the effects of streamline curvature, swirl, rotation, and rapid changes in strain rate in a more rigorous manner than one-equation and two-equation models, it has greater potential to give accurate predictions for complex flows [8].

Detached Eddy Simulation(DES)

Another approach is known as Detached Eddy Simulation (DES). This was first proposed by Spalart, in an attempt to combine the most favourable aspects of RANS and LES.

DES reduces to a RANS calculation near solid boundaries and a LES calculation away from the wall. ANSYS Fluent 12.1 offers a RANS/LES hybrid model based on the Spalart-Allmaras turbulence model near the wall and a one-equation SGS turbulence model away from the wall which reduces to an algebraic turbulent viscosity model for the SGS turbulence far from the wall [10].

Extend of modelling for certain CFD approaches for turbulence are illustrated in the Figure 5. It is clearly seen that the DNS and the LES models are computing fluctuation quantities resolve shorter length scales than models solving RANS equations. Hence they have the ability to provide better results. However they have a demand of much greater computer power than those models applying RANS methods [21].

(30)

Figure 5: Extend of modelling for certain types of turbulent models [21]

3.2 Large eddy simulation (LES)

Large eddy simulation (LES) is classified as a space filtering method in CFD. LES directly computes the large-scale turbulent structures which are responsible for the transfer of energy and momentum in a flow while modelling the smaller scale of dissipative and more isotropic structures. In order to distinguish between the large scales and small scales, a filter function is used in LES. A filter function dictates which eddies are large by introducing a length scale, usually denoted as∆in LES, the characteristic filter cutoff width of the simulation [10]. All eddies larger than ∆are resolved directly, while those smaller than∆are approximated.

3.2.1 Filtering of Navier-Stokes equations

In LES, the flow parametersφis separated into a filtered, resolved partφ¯and a sub-filter, unresolved part,φ⁰,

φ= ¯φ+φ⁰ (25)

Define a spatial filtering operation by means of a filter functionG(X, X⁰,∆)as follows:

φ(X, t)¯ ≡ Z ∞

−∞

Z ∞

−∞

Z ∞

−∞

G(X, X⁰,∆)φ(X⁰, t)dx⁰₁dx⁰₂dx⁰₃ (26) Whereφ(X, t)¯ = filtered function,φ(X, t)= original (unfiltered) function and∆= filter

(31)

cutoff width

Here, the overbar indicates spatial filtering, not time-averaging. Equation (26) shows that filtering is an integration,just like time-averaging in the development of the RANS equations, only in the LES the integration is not carried out in time but in three- dimensional space [1].

As mentioned, the filter function dictates the large and small eddies in the flow. This is done by the localized function G(X, X⁰,∆). This function determines the size of the small scales

G(X, X⁰,∆) =

(1/∆³ kX−X⁰k ≤∆/2

0 kX−X⁰k>∆/2 (27)

Various filtering methods exist, the top-hat filter is common in LES. The function represents Eq. (27). The top-hat filter is used in finite volume implementations of LES. The cutoff width is intended as an indicative measure of the size of eddies that are retained in the computations and the eddies that are rejected. In CFD computations with the finite volume method it is pointless to select a cutoff width that is smaller than the grid size. The most common selection is to take the cutoff width to be of the same order as the grid size. In three-dimensional computations with grid cells of different length∆x, width∆y and height∆z the cutoff width is often taken to be the cube root of the grid cell volume,

∆ = p³

∆x∆y∆z (28)

Now, the unsteady Navier-Stokes equations for a fluid with constant viscosity µ are as follows:

∂ρ

∂t + div(ρu) = 0 (29)

∂(ρu)

∂t + div(ρuu) =−∂p

∂x+µdiv(grad(u)) +Su (30)

∂(ρv)

∂t + div(ρvu) =−∂p

∂y +µdiv(grad(v)) +S_v (31)

∂(ρw)

∂t + div(ρwu) =−∂p

∂z +µdiv(grad(w)) +S_w (32) If the flow is also incompressible we havediv(u) = 0, and hence the viscous momentum source termsS_u,S_v andS_w are zero.Filtering of above equations,

∂ρ

∂t + div(ρ¯u) = 0 (33)

(32)

∂(ρu)¯

∂t + div(ρuu) =−∂p¯

∂x+µdiv(grad(¯u)) (34)

∂(ρ¯v)

∂t + div(ρvu) =−∂p¯

∂y +µdiv(grad(¯v)) (35)

∂(ρw)¯

∂t + div(ρwu) =−∂p¯

∂z +µdiv(grad( ¯w)) (36) This equation set should be solved to yield the filtered velocity field u,¯ ¯v and w¯ and filtered pressure field p. We need to compute convective terms of the form¯ div(ρφu) on the left hand side, but we only have available the filtered velocity field u,¯ v,¯ w¯ and pressure fieldp¯[1]. To make some progress we write,

div(ρφu) = div( ¯φ¯u) + (div(ρφu)−div( ¯φ¯u))

The first term on the right hand side can be calculated from the filtered φ¯and u¯ fields and the second term is replaced by a model. Substitution into above equation and some rearrangement yields the LES momentum equations:

∂(ρ¯u)

∂t + div(ρ¯u¯u) =−∂p¯

∂x +µdiv(grad(¯u))−(div(ρuu)−div(¯uu))¯ (37)

∂(ρ¯v)

∂t + div(ρ¯v¯u) =−∂p¯

∂y +µdiv(grad(¯v))−(div(ρvu)−div(¯vu))¯ (38)

∂(ρw)¯

∂t + div(ρw¯u) =¯ −∂p¯

∂z +µdiv(grad( ¯w))−(div(ρwu)−div( ¯w¯u)) (39) In these equations, the first two terms on the left hand side denote the rate of change and convective fluxes of filtered x−,y−and z−momentum. And third and forth terms on the right hand side denote the gradients in thex−,y−andz−directions and diffusive fluxes of filtered x−, y−and z−momentum. The last terms are caused by the filtered operation. They can be considered as a divergence of a set of stresses τij. In suffix notation thei−component of these terms can be written as follows:

div(ρuiu−divρu¯iu) =¯ ∂(ρuiu−ρu¯iu)¯

∂x +∂(ρuiv−ρu¯iv)¯

∂y +∂(ρuiw−ρu¯iw)¯

∂z = ∂τ_ij

∂xj

(40)

(33)

W here τ_ij =ρu_iu−ρu¯_iu¯ =ρu_iu_j−ρu¯_iu¯_j (41) The termτij is known as the subgrid scale (SGS) Reynolds Stress. Physically, the right hand side of Eq. (41) represents the large scale momentum flux due to turbulence motion.

The nature of these contributions can be determined with the aid of a decomposition of a flow variableφ(x, t) as the sum of (i) the filtered functionφ(x, t)¯ and (ii)φ⁰(x, t).

φ(x, t) =φ(x, t) +¯ φ⁰(x, t) (42) Using this decomposition in Eq. (41) we can write the SGS stresses as follows:

τ_ij =ρu_iu_j−ρu¯_iu¯_j = (ρu¯_iu¯_j−ρu¯_iu¯_j) +ρu¯_iu⁰_j+ρu⁰_iu¯_j+ρu⁰_iu⁰_j (43)

thus, we find that the SGS stresses contain three groups of contributions:

• Leonard stresses L_ij: L_ij =ρu¯_iu¯_j−ρu¯_iu¯_j which represent the interaction between two resolved scale eddies to produce small scale turbulence.

• cross-stressesCij: Cij =ρu¯iu⁰_j+ρu⁰_iu¯jare the cross-stress terms that describe the interaction between resolved eddies and small-scale eddies.

• LES Reynolds stresses Rij: Rij = ρu⁰_iu⁰_j is the subgrid scale stress that represents the interactions between two small scale eddies

3.2.2 Smagorinksy-Lilly SGS model

To approximate the SGS Reynolds stress, a SGS model can be employed. The most commonly used SGS models in LES is the Smagorinsky-Lilly model. In a flow, it is the shear stress and the viscosity of the flow that cause the chaotic and random nature of the fluid motion. Thus, in the Smagorinsky-Lilly model, the effects of turbulence are represented by the eddy viscosity based on the well known Boussinesq hypothesis.

The Boussinesq hypothesis relates the Reynolds stress to the velocity gradients and the turbulent viscosity of the flow [8]. It is therefore assumed that the SGS Reynolds stress R_ij is proportional to the modulus of the strain rate tensor of the resolved flow S¯ij = ¹₂(∂u¯i/∂xj+∂u¯j/∂xi)

Rij =−2µ_SGSS¯ij +1

3Riiδij =−µ_SGS ∂u¯i

∂x_j +∂u¯j

∂x_i

+ 1

3Riiδij (44) where µ_SGS is the SGS eddy viscosity. Leonard stresses and cross-stresses are lumped together with the LES reynolds stresses in the current versions of the finite volume

(34)

method. The whole stress τ_ij is modeled as a single entity by means of a single SGS turbulence model.

τij =−2µ_SGSS¯ij +1

3τiiδij =−µ_SGS ∂u¯i

∂x_j +∂u¯j

∂x_i

+1

3Riiδij (45) The Smagorinksy-Lilly SGS model builds on Prandtl’s mixing length model and assumes that we can define a kinematic SGS viscosityν_SGS, which can be described in terms of the one length scale and one velocity scale and is related SGS viscosity by ν_SGS =µ_SGS/ρ.

Since the size of the SGS eddies is determined by the details of the filtering function, the obvious choice for the length scale is the filter cutoff width∆. The velocity scale is expressed as the product of the length scale∆and the average strain rate of the resolved flow ∆× kSk, where¯ kSk¯ =

q

2 ¯SijS¯ij. Thus, the SGS viscosity is evaluated as follows:

µSGS =ρ(CSGS∆)²kSk¯ =ρ(CSGS∆)² q

2 ¯SijS¯ij (46)

WhereCSGS = constants and S¯ij = ¹₂ ∂u¯i

∂xj +^∂_∂x^u^¯^j

i

where CSGS is the Smagorinsky constant that changes depending on the type of flow.

For isotropic turbulent flow, theC_SGS value is usually around 0.17 to 0.21. Basically, the Smagorinsky SGS model simulates the energy transfer between the large and the subgrid- scale eddies. Energy is transferred from the large to the small scales but backscatter sometimes occurs where flow becomes highly anisotropic, usually near to the wall [4]. To account for backscattering, the length scale of the flow can be modified using Van Driest damping, and suggested that C_SGS = 0.1 is most appropriate for this type of internal flow calculation.

The Smagorinksy model has been successfully applied to various flows as it is relatively stable and demands less computational resources among the SGS models. But some disadvantages of the model have been reported,

• Too dissipative in laminar regions.

• Requires special near wall treatment and laminar turbulent transition.

• CSGS is not uniquely defined.

• Backscatter of flow is not properly modelled

Germano and Lilly conceived a procedure in which the Smagorinsky model constant, C_SGS, is dynamically computed based on the information provided by the resolved scales of motion. Hence, the dynamic SGS model has been introduced. This model employs a

(35)

similar concept as the Smagorinsky model, with the Smagorinsky constantC_SGSreplaced by the dynamic parameter C_dym [8]. The parameter C_dym is computed locally as a function of time and space, which automatically eliminates the problem of using constant CSGS. In the dynamic SGS model, another filter is introduced which takes into account of the energy transfer in the dissipation range. Performing the double filtering allows the subgrid coefficient to be calculated locally based on the energy drain in the smallest scales. Generally, the dynamic model predicted better agreement with experimental work in region of transition flow and the near wall region.

Some advantages of the dynamic model over the Smagorinsky-Lilly models are,

• Dynamic SGS automatically uses a smaller model parameter in isotropic flows.

• Near the wall, the model parameters need to be reduced; the dynamic SGS model adapts these parameters accordingly.

3.3 Reynolds-averaged Navier-Stokes equations

The basic tool required for the derivation of the Reynolds-averaged Navier-stokes (RANS) equations from the instantaneous Navier–Stokes equations is the Reynolds decomposition. Reynolds decomposition refers to separation of the flow variable into the mean component and the fluctuating component[1]. The following rules will be useful while deriving the RANS equation. We begin by summarising the rules which govern time averages of fluctuating properties φ = Φ + ´φ and ψ = Ψ + ´ψ and their summation, derivatives and integrals:

φ⁰ =ψ⁰ = 0 Φ =φ ∂φ

∂s = ∂Φ

∂s Z

φds= Z

Φds

φ+ψ= Φ + Ψ φψ = ΦΨ +φ⁰ψ⁰ φΨ = ΦΨ φ⁰Ψ = 0 (47) In addition, div and grad are both differentiations, the above rules can be extended to a fluctuating vector quantity a=A+ ´a and its combinations with a fluctuating scalar φ= Φ + ´φ:

diva= divA; div(φa) = div(φa) = div(ΦA) + div(φ⁰a⁰);

div gradφ= div grad Φ (48)

Now, we consider the instantaneous continuity and Navier-Stokes equations in a Carte- sian co-ordinate system so that the velocity vector u has x-component u, y-component v and z-component w:

(36)

divu= 0 (49)

∂u

∂t + div(uu) =−1 ρ

∂p

∂x +µdiv(grad(u)) (50)

∂v

∂t + div(vu) =−1 ρ

∂p

∂y +µdiv(grad(v)) (51)

∂w

∂t + div(wu) =−1 ρ

∂p

∂z +µdiv(grad(w)) (52) This system of equations governs every turbulent flow, but we investigate the effects of fluctuations on the mean flow using the Reynolds decomposition in equations (49), (50), (51) and (52) and replace the flow variables u and p by the sum of a mean and fluctuating component. Thus

u=U+u⁰ u=U +u⁰ v=V +v⁰ w=W +w⁰ p=P+p⁰

Then the time average is taken, applying the rules stated in equations (47) and (48).

Considering the continuity equation (49), first we note thatdivu= divU. This yields the continuity equation for the mean flow:

divU= 0 (53)

A similar process is applied on the x-momentum equation (50). The time averages of the individual terms in this equation can be written as follows:

∂u

∂t = ∂U

∂t div(uu) = div(UU) + div(u⁰v⁰)

−1 ρ

∂p

∂x =−1 ρ

∂p

∂x νdiv(grad(u)) =ν div(grad(U)) Substitution of these results gives the time-average x-momentum equation

∂U

∂t + div(UU) + div(u⁰u⁰) =−1 ρ

∂P

∂x +νdiv(grad(U)) (54) Repetition of this process on equations (51) and (52) yield the time-average y- and z-momentum equations:

(37)

∂V

∂t + div(VU) + div(v⁰u⁰) =−1 ρ

∂P

∂y +νdiv(grad(V)) (55)

∂W

∂t + div(WU) + div(w⁰u⁰) =−1 ρ

∂P

∂z +νdiv(grad(W)) (56) Note that the terms (I), (II), (IV) and (V) in equations (54), (55) and (56) also appear in the instantaneous equations (50), (51) and (52), but the process of time averaging has introduced new terms (III) in the resulting time-average momentum equations [1]. This terms is product of fluctuating velocities and are associated with convective momentum transfer due to turbulent eddies. And put these terms on the right hand side of equations (54), (55) and (56) to reflect their role as additional turbulent stresses on the mean velocity components U, V and W:

∂U

∂t + div(UU) =−1 ρ

∂P

∂x +νdiv (grad(U)) +1

ρ

"

∂(−ρu⁰²)

∂x +∂(−ρu⁰v⁰)

∂y + ∂(−ρu⁰w⁰)

∂z

#

(57)

∂V

∂t + div(VU) =−1 ρ

∂P

∂y +νdiv (grad(V)) +1

ρ

"

∂(−ρu⁰v⁰)

∂x +∂(−ρv⁰²)

∂y + ∂(−ρv⁰w⁰)

∂z

#

(58)

∂W

∂t + div(WU) =−1 ρ

∂P

∂z +νdiv (grad(W)) +1

ρ

"

∂(−ρu⁰w⁰)

∂x +∂(−ρv⁰w⁰)

∂y +∂(−ρw⁰²)

∂z

#

(59)

The extra stress terms have been written out as follows. They result from six additional stresses among of three normal stresses

τ_xx=−ρu⁰² τ_yy=−ρv⁰² τ_zz =−ρw⁰² (60) and three shear stresses

τxy =τyx=−ρu⁰v⁰ τxz =τzx=−ρu⁰w⁰ τyz =τzy =−ρv⁰w⁰ (61) These extra turbulent stresses are called the Reynolds stresses. The normal stresses involve the respective variances of the x-, y- and z-velocity fluctuations.They are always

Numerical Investigation of Flow Past a Circular Cylinder and in a Staggered Tube Bundle Using Various Turbulence Models