Module: Research Methods (MACE61078)

Date of Submission: 19th April 2013.

Submitted to

Prof Hector Iacovides

Submitted By

Jahedul Islam Chowdhury

Student ID# 8572073

MSc Thermal Power and Fluid Engineering

School of Mechanical, Aerospace and Civil Engineering

The University of Manchester

Contents

1

Introduction

1

2

Analytical Wall Function Treatments

2

2.1

UMIST-A Scheme

3

2.2

UMIST N Scheme

5

3

Wall Function approaches using two equation model

7

4

Wall Function for oscillating flow

10

4.1

Periodic Channel flow- UMIST-N Approach

10

4.2

Oscillating Turbulent flow- UMIST-A Approach

12

5

Direct Numerical Simulation

13

6

Conclusions

15

Bibliography

16

Nomenclature

k

Turbulent kinetic energy

?

Dissipation rate of k

�"

specific dissipation rate of k

C?

turbulent model constant

Cl

log-law constant

Nu

Nusselt number

Re

Reynolds number

ReDh

Reynolds number based on hydraulic diameter

Reτ

Reynolds number based on frictional velocity

U+

Nondimensional velocity

y+

Dimensionless distance

?

Dynamic viscosity

ρ

Density of the fluid

?t

Eddy or turbulent viscosity

?t*

Dimensionless eddy or turbulent viscosity

y*

Dimensionless distance to the nearest wall

Bo

Buoyancy parameter

τw

Wall shear stress

Pk

Production of k

Lx

Length in x direction

Ly

Length in y direction

Lz

Length in z direction

uτ*

Instantaneous velocity fluctuation

T

Wall temperature

1. Introduction

The Computational Fluid Dynamics (CFD) is the computer solutions of a number of fundamental equations to predict pretty much any type of fluid motion. CFD may sometimes be used to investigate the physical behaviour of flow, but the common goal is to provide an input to engineering analysis and design. The CFD has developed very rapidly in last two decades because of the advancement of technology and increase of computing processor capacity of computer. Considering the development of technology, Researchers focuses toward more complex calculations. With the high-speed supercomputers, better solutions of the complicated flow can be achieved.

Now a days, CFD and experimental analysis are extensively used as a design tool in most of the industries. The position of CFD is currently behind the experimental analysis due to the fact that CFD does not produce absolute results. The reason for this is that the numerical methods, which govern the solutions in a CFD problem, rely on several modelling assumptions that may not have been validated to a satisfactory level. However, CFD presently offers itself as a powerful design tool and even more so in the future because of: firstly, dangerous or expensive trial and error experiments whereas CFD reduce the cost and time that usually require for the experiment. Secondly, the numerical schemes and physical models that are the building blocks of CFD are continually improving. Moreover, it can be able to produce reliable and accurate results on one particular design compare to the exact one. For this reason, most of the researchers prefer CFD rather than Experimental analysis.

The CFD analysis starts with dividing or discretizing the geometry to be modelled into usually a large number of small computational cells. Discretization is the method of approximating the differential governing equations by a system of algebraic equations for the variables at some set of discrete locations in space and time. The discrete locations are normally referred to as the "grid" or the "mesh". This analysis is iteratively solving equations for the conservation of mass, heat, and momentum, and the dissipation of motion into turbulence. The additional effects mentioned above are added into the solution procedure if needed. The results are often presented graphically to assist in visualizing what's actually happening in the flow domain

Most flows of relevance to engineering applications are clearly turbulent, and most CFD codes are designed to compute and predicting the behaviour of turbulent flow. Most commonly, this is accomplished by the adoption of a turbulence model that makes use of certain assumptions to approximate flow behaviour using more tractable equations. This precludes the possibility of predicting every feature of the flow, but the parameters providing the greatest engineering relevance may be estimated at a greatly reduced computational cost.

CFD approximations of turbulent flow near the walls generally employ one of two broad strategies to resolve the very influential, complex, but thin near-wall viscosity-affected sub-layer. One method uses a fine numerical mesh and a turbulence model which take viscous influences into account; the other use a “wall functions’’. However, because of the complexity of the three-dimensional time-dependent flows in engines, the use of low-Reynolds number turbulence models to resolve the near wall regions is computationally very expensive. Thus wall functions are applied in these regions, coupled to a high Reynolds number turbulence model for the bulk flow. The functions used apply the log-law for velocity variation across the near wall regions derived from simple wall shear flows. However, because of the oscillating nature of the flow in engines due to the piston motion, it is not clear that such wall functions are appropriate, as the near-wall flow will also vary in an unsteady fashion.

The aim of the literature survey of this project is to study the wall function to an oscillating turbulent flow imposed on top of a steady through flow and analyse earlier works on wall function near the wall. And hence, ultimately, gather the knowledge to develop an alternative analytical wall function that does properly take the effects of oscillation into account. The scope of the literature survey is to cover the study of different turbulence modelling approaches i.e. two equation models (k-ε, k-ω), Direct Numerical Simulation (DNS) etc. However, major focus of the research would be the near wall treatment of the turbulent flow for instances, wall function, sub-grid approaches (UMIST-A and UMIST-N). Moreover, previous works of the wall function for oscillating turbulent flow will also be covered by this survey.

2. Analytical Wall function Treatments

Two separate projects have done on developing more general wall-functions at UMIST. The pathways of schemes are different with each other. One is UMIST�"A which is an analytical scheme, (Unified Methodology for Integrated Sub-layer Transport�"Analytical), the other is UMIST�"N, based on numerical approach [1].

2.1 UMIST �"A Scheme:

An approach was developed at UMIST by Craft et al. [2] takes into account the pressure gradients, buoyancy forces, and changes in the thickness of viscous sub-layer among other parameters by the analytical integration of the momentum and energy equations. This approach is known as the Unified Modelling through Integrated Sub-layer Treatment �" an Analytical approach (UMIST-A). This method reduced the computing cost by one or two order of magnitudes compared to a standard low-Reynolds-number approach.

UMIST-A scheme was designed for the safety studies specially convection (forced, mixed or natural) flow on near vertical surfaces, strong variations of molecular transport properties across the VSL, and laminarization, i.e. a marked thickening of the VSL in buoyancy-aided mixed convection.

Figure1 Viscosity distribution assumed over near wall cell.

Turbulent viscosity (b) Molecular viscosity [1]

The prescribed viscosity distribution was taken as the form for .The coefficients and Cl are 0.09 and 2.55 as conventional 1-equation turbulence model and .It was proved that evaluating k in the at the sub-layer interface, yv much less stable than adopting the nodal value and also lead to greater dependence on the size of the near-wall cell because to extrapolate values to yv need to use the information from the wall than yP. However, the variation of molecular properties across the sub-layer due to temperature changes also taken into consideration. The paper showed another area where it was felt appropriate to improve current practice was in the prescription of the kinetic energy dissipation rate next to the wall, Chieng and Launder [3] had approximated the exact result of Jones and Launder [4]

To make the treatment sensitive to laminarization, Researchers made two choices, namely appropriate parameters to use as detector and operand. After extensive testing, it was claimed that the ratio of the shear stress between the wall and the edge of the sub-layer, λ was the best detector. In case operand, a direct choice was adopted: the mean level of dissipation rate over the near-wall control volume was adjusted by a weighting function F that in turn was a function of λ.

UMIST-A was tested for two fully developed pipe flow with different near wall cell size. However, at a lower Reynolds number of 6753 the experimental data of Kudva and Sesonske [5] lie above the supposedly ‘universal’ log-law, as do predictions with the low-Reynolds-number (LRN) k�"ε model of Launder and Sharma [6]. Crucially, current wall function results also agree with the data. Moreover, adding to most such schemes, it could be noted that the UMIST�"A scheme results show scarcely any sensitivity to the size of the near-wall cell. The paper shows the Variation of Nusselt number in mixed-convection up-flow in a vertical tube along with UMIST�"A treatment, UMIST�"A without F(λ), experimental data of Li [7], Dittus�"Boelter, LRN calculation in a graph. Result shows, that Wall function did agree with the experimental data. One run was done by UMIST-A without F(λ) results 20% increase in Nusselt number.

This scheme was tested for the Nusselt number variation along a heated annulus wall at a low buoyancy parameter Bo = 0.22 and 2.89 .and compared with experiments of Jackson et al. [8]; LRN: low-Re model calculations; standard wall function predictions; Analytical Wall function. It was found that at a low level of buoyancy the full low-Reynolds-number model of Launder and Sharma [6] and both wall function treatments broadly capture the small increase in Nusselt number. Nevertheless, at higher values of the buoyancy parameter, both the UMIST�"A scheme and the full low-Reynolds-number do capture the now substantial increase in Nusselt number found in the experiments, whereas the standard wall functions do not.

The overall features of the UMIST- A are:

Convective transport (conventional wall treatment ignore it in many cases) expressed by a simplified form.

In a case where buoyancy is an important fact, buoyant force of the near wall cell can be obtained by integrating a fit to the analytical temperature profile over the cell rather than the force purely on the temperature at the near wall cell itself.

If the viscous sub-layer thickness exceeds the cell thickness, yn (as it may do in limited regions if a structured grid is adopted), although a further reformulation is essentially needed, but the identical principle can still be used.

In conclusion to this paper, UMIST-A is based on the analytical solution of simplified near-wall momentum and temperature equations, accounting for pressure gradients and other force fields such as buoyancy. This scheme approaches have been applied to a range of flows in which standard log-law based wall functions are known to perform badly. The overall methodology of UMIST�"A, is closer to that employed in standard wall-functions, thus possibly making it the easier scheme to incorporate into existing flow solvers.

2.2 UMIST-N Scheme:

The analytical approach is not enough equipped in case where the velocity profile parallel with the wall undergoes strong skewing across the sub-layer, for instances, non-orthogonal impingement of flow on a bank of heat-exchanger tubes, for flow with strong streamline curvature where the strain stress relation provided by eddy-viscosity model does not depict the turbulence generation process etc. Considering above difficulties, a numerical scheme, UMIST�"N, has been developed by Gant [9] and Craft et al. [10]. The scheme is much more similar to the low-Reynolds-number models in that the wall-function cell is itself sub-divided into, typically, 30 thin slices as shown in figure 2. With suitable simplifications, the mean flow and turbulence differential equations are solved effectively as a one-dimensional problem across this fine grid. Which ultimately generate data required as of wall function quantities such as wall shear stress, averaged source term etc. and supply appropriate boundary conditions for the whole-field solution carried out on the primary grid.

Figure 2. Treatment of near wall control volume in UMIST- N [1]

There were several simplifications made by the Gant [9] in order to reduction in computer time that one looks form the wall function.

The pressure gradient parallel to the wall is considered uniform across all the sub-grids, equal to the pressure gradient across the near-wall cell of the primary grid.

The velocity component normal to the wall can be obtained by continuity rather than by solving the momentum equation normal to the wall.

The paper shows that two alternative turbulence models have been employed within the above numerical treatment. One is the LRN k�"ε model of Launder and Sharma [6] and other is cubic non-linear eddy viscosity model of Craft et al. [11]. However, result shows that the cubic terms in the latter model make it far more sensitive to streamline curvature than a linear EVM.

Figure 3. Nusselt number variation on a flat plate beneath an axisymmetric impinging jet. (a) Linear k�"ε EVM. (b) Cubic non-linear EVM (Craft et al [11]). Symbols: experiments of Baughn et al [12]; heavy line: full LRN treatment; other lines: UMIST�"N, with different near-wall cell sizes [1]

The paper shows the first test case is for a turbulent jet impinging orthogonally onto a uniformly heated flat plate. The jet discharges from a long smooth pipe whose exit is four diameters above the plate. Figure 3 Shows the resultant variation of Nusselt number over the plate from the stagnation point (r=0) outwards. As is seen in Figure 3(a), the computed Nusselt number at the stagnation point is more than twice the measured value. However, the computation results of wall function and full LRN are close enough. Again same test was done for the non-linear EVM (Craft et al [11]) instead of Linear k�"ε EVM. Now the results from computation and experiment are much closer. The only major difference between the UMIST�"N results and those of the complete LRN treatment is in the computer time required. It is suggested from the paper that the wall-function result with the same grid density takes less than one eighth of the time required for the complete low-Reynolds-number model.

Another test example considered heat transfer from a mildly heated disc spinning about its own axis. The disc’s rotation induces a radially outward motion that peaks outside the VSL, the tangential velocity increases rapidly.

Figure 4. Radial velocity profile for spinning disc in wall-layer coordinates. Solid line: LRN calculation; broken line with symbols: UMIST�"N; chain line: log-law; other lines: standard wall-function treatments [1].

Result shows that the induced radial velocity predicted with UMIST-N (using the linear EVM of Launder and Sharma [6]) agrees very closely with the results of the corresponding LRN computation. However, plot of integral Nusselt numbers shows the negligible difference among the treatments except the computation time. Eventually, the complete low-Reynolds-number model required thirteen times more computation time than UMIST�"N.

In closing to this paper, UMIST�"N is based on a local one-dimensional numerical solution of the governing equations. It is clearly perform badly for the standard log-law based wall functions. This scheme is, in principle, the more general approach, but is inevitably more computationally expensive terms of CPU time and storage requirements.

3. Wall function approaches using two equation model

The paper titled ‘‘The analysis of wall function approaches using two equations turbulence model’’ was done by the Pe´rez-Segarra et al [13].The objective of this paper was to obtain tools which were able to approximate both fluid flow and heat transfer with less time and storage consumption, reduction of grid-sensitivity and a relatively good accuracy. However, considering the scope of literature review, fluid flow sections are going to consider only.

Authors focused on WF approaches, using for inner nodes different high Reynolds number two-equation models, depending on the dissipative variable (k�"ε or k�"ω treatments). Several assumptions and numerical corrections were taken into consideration for near wall treatments. The governing equations of the time-averaged continuity, momentum and energy of the fluid flow were considered. A low- Reynolds �"number model was tested and compared to the Ince�"Launder [14] k�"ε LRN model and with Wilcox [15] k�"ω LRN model in order to find out the limitation and improvements of the numerical and computational phenomena of the wall functions. However, the governing equations were solved using finite volume techniques over a staggered discretization. A structured grid, fully implicit time integration and SIMPLE method were employed to the numerical approach. Diffusion and convection terms were calculated with the central different and UPWIND and SMART scheme for the channel flow and UPWIND and Power law for backward flow respectively. All computations were carried out using a pseudo-transient iterative algorithm, applying the biggest time-step which ensures convergence. Post processing was done using generalized Richardson extrapolation for h-refinement studies and the grid convergence index (GCI) method due to verify the numerical results provided by the continuous (LRN) models analysed in his work.

Different cases were tested shown in this paper. For instances, A channel flow at moderate Reynolds number (Reτ =590, ReDh ≈ 43,000), A backward facing step flow (BFS) at Reynolds number of ReH= 37,500. The inlet conditions used for test cases were identical. For channel flow, , in addition to a turbulence intensity of 7% (kin= 0.072 ) and a specific dissipation rate and a fixed inlet temperature of 300K.

8H

(a) (b)

Figure 5. (a) Channel and (b) backward facing step test cases [13]

Channel Flow case:

The results were plotted as U+ against y+ for different mashes and eventually compared with the DNS data by Moser et al. [16]. The computation results of wall function showed that standard approach (SWF1: wall function treatment proposed bay Launder [17] and k�" ϵ modification (SWF2: modification of SWF1) were failed when near-wall node was placed at viscous regions. Though, improvement achieved when near-wall node was place at buffer layer and sublayer, above two prediction was not appropriate below y+<10. However, results from k-ω models were acceptable manner. Meanwhile, plot of k+ against y+ showed that turbulent kinetic energy predictions was comparatively good than the DNS data, especially with the WWF2 treatment. It was surprisingly found that turbulent suddenly increased at near-wall nodes placed between 5 < y+ < 20 for the second wall cell for WWF1 (ω wall function 1), computations. Overall, All treatments agreed with DNS data by Moser et al [16] for near-wall cells placed at logarithmic region since they are using pure log-law for y+ > 30. On the other hand, LRN computation was done with the finest mesh. Due to finer grid approach, it was found the CPU time consumption was much higher than the WF approach, though the accuracy was more or less the same for all treatments.

Backward Facing step (BFS) Flow Case:

This type of flow require complex and related to different wall mesh grid. However, for WF treatments, seven different meshes were used to test grid-sensitivity of models whereas for LRN, it was five. It was found that using coarsest meshes, WWF1 and WWF2 treatments present a better behaviour than SWF1 and SWF2 approaches. In case of Wall function computation, If meshes refined, k�"ω approaches tend to the asymptotic results provided by the high Reynolds number k�"ω model by Wilcox et al.[18] whereas k�"ϵ predictions were becoming less accurate. However, LRN model provides different result. In order to achieve asymptotic results, near-wall cells placed at y+= 0.5 in both walls (development channel and downstream of the step). This is the main reason to explain the extremely dense mesh used and the high CPU time needed for both simulations.

In conclusion of this paper it can be said that, researcher proposed different strategies for the improvement of the grid sensitivity for standard SWF1 treatment [17] which were based on the enhancing ϵ profile using a two-layer integration (SWF2 approach); changing dissipative variable from ϵ to ω and also from k�" ϵ (Launder, D. Spalding [19]) to k�" ω (Wilcox [18]) for HRN model would be used for inner nodes (WWF1 approach) and, finally, designing a blending function to reduce grid-sensitivity on an ω platform (WWF2 approach).

4. Wall Functions for Oscillating Flow

Many Researchers were seeking for resolving near wall region by wall function for oscillation flow. Among of which UMIST-N and UMIST-A wall functions for periodic flow application will be considered only.

Periodic flow �" the UMIST-N Approach

UMIST-N near wall treatment applied to periodic channel flow was done by Richards [20]. Richards used a parabolic coded PASSABLE (PArabolic Solution Scheme Applied to Boundary Layer Equations) for the simulation. This code is first developed by Addad [21] in his PhD thesis at UMIST and later modified by Richards for periodic flow.

Channel flow with walls that were a distance 2δ apart was used for the investigation. Flow was driven by pressure gradient in the x direction. The problem was solved from y=0 at wall to y=δ at the halfway of the channel. Flow was considered independent of z direction and x direction as well. RANS equations for continuity and conservation of momentum were used for the channel flow. However, governing equations of the flow were discretised by the using of finite volume method. The actual grid used within his work (figure 6) contained 98, 18 and 50 nodes for the for low- Reynolds- number, high-Reynolds-Number, and sub-grid approach respectively.

(b) (c)

Figure 6. Grid for (a) Low-Reynolds-number (b) High-Reynolds-number and (c) sub-grid mesh by Richards [20]

However, Richards used separate under-relaxation factor for main and sub-grid as follows:

(b)

Figure 7. Under-relaxation factor for (a) main grid and (b) sub-grid [20]

This paper shows that, he used four configurations of turbulence models for his simulation work. Firstly, named low-Re k-ε where a low-Reynolds-number standard k-ε model was employed without using sub-grid. Second case was called k-ε with log law which was similar to the first configuration but now a high-Reynolds-number standard k-ε model was used with the log-law boundary condition at near the wall. However, Richards third approach was to use sub-grid k-ε. Now, a high-Reynolds-number standard k-ε model was used in the main grid with the near-wall boundary condition derived from the sub-grid solution and the sub-grid employed the low-Reynolds-number standard k-ε model. Finally, sub-grid k-ω configuration was used for the periodic flow. He employed the 1988 k-ω model in the main grid with the near-wall boundary condition that derived from the sub-grid solution. Same model also used for the sub-gird.

The tests were run over major two cases. These were steady channel flow and periodic channel flow. Furthermore, periodic channel flow was driven from periodically variable pressure gradient and periodically variable bulk flow rate. Results shows in the paper that the periodic cases agreed well with the pressure gradient fluctuations and bulk flow variation of the DNS study of Kawamura & Homma [22]. It was clearly shown that the UMIST-N was successfully applied for the computation of periodic flow.

The results for the log-law of wall, k-ε sub-grid and low-Reynolds-number k-ε configurations were similar for the steady case. However, in the unsteady case the sub-grid approach was given better results than the log-law and had the similarities to the low-Reynolds number approach. The main advantage of the UM IST-N approach was its ability to achieve to lower refinement the near-wall grid compared to the low-Reynolds number grid without hindering the main grid. The k-ω on the other hand had acceptable results overall but the performance was significantly influenced due to the questionable behaviour between the main and sub-grid. This was due to ω tending to infinity at the wall and the averaged production of ω at the sub-grid could not be obtained. Thus, Richards [24] was come into conclusion not to apply UMIST-N approach to a k- ω model.

In conclusion of this paper, Results obtained using the UMIST-N sub grid model are considerably more accurate in terms of the near-wall variations in mean velocity and turbulence parameters than those given by the low-low. However, this approach requires two separate grids and there are some uncertainties associated with the interface between them.

4.2 Oscillating Flow- UMIST-A Approach

UMIST-A analytical wall function approach was investigated by Bauly [23]. The aim of the project was to see whether this would give the results of equal accuracy to those obtained with UMIST-N. Flow considered in this research was 1-D channel flow. The flow conditions were exactly the same as those tested by Richards [20]. The simulation was run by the code called PASSABLE (Parabolic Solution Scheme Applied to Boundary Layer Equations). This originally modified for the UMIST-N periodic flow. In addition to this code, UMIST-A wall function was inserted to it.

Bauly used the same grids and wall boundaries on k�"ε and k-ω as Richards [20] did use in UMIST-N periodic flow approach. In addition to the UMIST-A, he added source terms to the model. For instances, - was added to the U momentum equation, was added to the K transport equation. However, Bauly used the same under relaxation factor for the main and sub-grid as Richards except for ε=1.01 instead of 1.0 in the main grid.

Paper shows that, the cases which tested by Richards [20] were re-run in order to make sure the code works properly. Bauly [23] only considered sub-grid models. . The code was run for both steady and unsteady flow situation. For steady case, the values of U+ and k+ were plotted against the non-dimensional distance y/δ and compared with the DNS results by Kim et al [24]. The k�"ε, 1988 k-ω, and high-Reynolds �"number main grid with log-law wall function produced same results as Richards but results from low-Reynolds-number configuration was unacceptable. Though, UMIST-A wall function underpredicts the values of U+, the results were close to the DNS and thus acceptable. Moreover, the values of k+ also underpredicts for UMIST-A at near the wall because of the higher dissipation rate of k. Result of U+ and τw+ were plotted against period in the case of periodic flow and found underpredicted than the DNS results by Kawamura et al [22]. This was because of the pressure gradient of U momentum equation of the wall shear stress not only drives the flow but also acts to overcome the wall shear stress. It was shown in the paper that, results of U+ for UMIST-A are less accurate away from the wall because of the faulty calculation of turbulent shear stress. Results from the unsteady case shows that UMIST-A wall function produced significant difference in production of k+ than the other models.

In summary to this paper, the application of UMIST-A wall function to the oscillating flow produced poor results especially for the periodic flow. The reason was primarily noted by the Researcher that the UMIST-A was originally developed for the heat transfer type of flow whereas current investigation carried out on the isothermal flow driven by the pressure gradient. Moreover, turbulent kinetic energy results were extremely poor for the periodic flow. It was assumed that reversing of velocity gradient influenced by the backward flow in an oscillating flow that gave the wrong calculation of the turbulent kinetic energy. However, it can be noted that the additional terms added to the standard wall functions did not adequately express the effects of oscillations.

5. Direct Numerical Simulation (DNS)

Direct Numerical Solution (DNS) solves the Navier-Stokes equation numerically without the use of any turbulence models. Thus, DNS solves for the whole range of spatial and temporal scales of turbulence from the smallest to the largest scale. To do this a very fine mesh is required. The number of nodes required for a DNS simulation is proportional to Re9/4, where Re is the Reynolds number. Hence, DNS simulations are very computationally expensive and are usually used to simulate flows of simple geometries and at low Reynolds number. DNS is highly accurate and its data is used to validate turbulence models.

The paper considered here is DNS of an unsteady turbulent channel flow driven by a temporary sinusoidal pressure gradient by Kawamura et al [22]. Most of the DNS analysis assumed the flow is steady but here, a temporally sinusoidal pressure gradient was imposed on the turbulent channel flow. The computational domain specified as Lx=25.6δ, Ly= 2δ and Lz= 3.2δ. the basic continuity and momentum equation of the channel was used where sinusoidal pressure gradient was defined as: , where variables were normalized by channel half width δ and the friction velocity uτs for the steady state. The amplitude A was set to 6.0. The time of one cycle was selected so that an unsteady effect appears while flow reversal would absent, hence, T+ was 6.0. Reynolds number, Reτs =180 was taken for the situation at A=0. However, a fully developed in stream-wise direction flow was assumed and periodic boundary conditions were applied to the x and z directions.

Results from the simulation were plotted in order to get the mean velocity profiles, turbulent kinetic energy, visualization of the turbulence structure and two point correlation. The mean velocity was averaged over the channel section and turbulence quantities were within each phase over 15 cycles. It was shown that, mean velocity increased when pressure gradient positive and decreased for the negative gradient. The results of mean velocity which normalized by the instantaneous uτ* showed that profile were flattened at the centre of the channel during acceleration phase and were more peaked in the deceleration part. Mean velocity was re-plotted in the semi-logarithmic coordinates and was found mostly larger than then well-known steady stale logarithmic distribution and a kind of laminarization observed. However, friction coefficients was plotted against the instantaneous Reynolds number and fond larger than the steady state correlation in phase 1-2 whereas it was smaller in most of the period. The dependence of wall shear stress on Reynolds number was more declarable than the steady case.

The normalized turbulence kinetic energy distribution was plotted across the channel. The results were pretty different from the steady state. It was larger than the steady state when pressure gradient decreasing and smaller when increasing. Moreover, peak of the turbulent energy for both acceleration and deceleration periods were increased away from the wall with compare to the steady state. Authors also normalized the kinetic energy by instantaneous friction velocity and showed that it was larger than the steady state for deceleration period whereas smaller for acceleration. On the other hand, total shear stress normalized by instantaneous friction velocity in order to depict the difference between acceleration and deceleration period. Total values of shear stress for acceleration were decreased while increased in deceleration period. The budget of the turbulent kinetic energy was plotted for the acceleration period (phase 6 and 8) and deceleration period (phase 12 and 14). In the acceleration section, the peak of production was quite smaller than the steady state value of 0.25 because of Reynolds stress decrease at this point. But for the deceleration period, the pressure gradient term becomes negative and thus, shear stress increased significantly. As a result of this effect, the turbulent kinetic energy and production also increased.

The plot of high and low speed streaks in acceleration and deceleration period were shown in paper. In case of acceleration, high and low speed strakes long enough so that they were connected from upstream to the downstream boundaries. However for the deceleration case, the structures were broken into several complex segments and elongated in the span-wise direction. Moreover, in the central region, many ramped structures were appeared. In the streamwise vorticity fluctuations, it was elongated spanwise for deceleration and there was no elongated streaky structure observed for the acceleration period. Two point correlations in x direction was used by the Authors to examine the adequacy of the computational domain. They come into conclusion that, streamwise domain was not large enough to apply the periodic boundary conditions for some of the phases.

6. Conclusions

A critical study on various turbulence models have been carried out for gathering the basic concepts of the computation. In CFD it is more challenge to treat where the viscous effects become important. It is found from the studies that using wall function to resolve near wall region is relatively less expensive rather than the low Reynolds number turbulence number. The famous k and ? wall function treatments are widely used in industrial CFD simulation. Study shows that this model gives less accurate result in near wall region. That is why; k-�" model has been introduced, which was reported to perform better than ? based schemes in the near wall adverse pressure gradient flow. However, recently, more advanced schemes (UMIST-A and UMIST-N), which remove the assumption of log-law have been developed. It’s been studied that analytic integration of mean momentum equation for UMIST-A can be used to estimate the wall shear stress and the production of k. However, the second approach (UMIST-N: numerical wall function) involves putting a 1-D ‘sub-grid’ across each main near wall cell. This scheme basically solves the transport equations for each sub-grid and quantities such as wall shear stress, production and dissipation were computed from these local solutions. On the other hand, solving near wall region in case of periodic flow is more challenging than simple shear flow. Study shows that UMIST-N and UMIST-A schemes were adopted and developed for the case of oscillating flow. The former one has successfully implemented to the period flow but latter one not to be the up to mark. Finally, DNS of unsteady flow has also been studied. The results of this method are more accurate and most of the papers studied here compare the results from turbulence models with that of DNS one.

Date of Submission: 19th April 2013.

Submitted to

Prof Hector Iacovides

Submitted By

Jahedul Islam Chowdhury

Student ID# 8572073

MSc Thermal Power and Fluid Engineering

School of Mechanical, Aerospace and Civil Engineering

The University of Manchester

Contents

1

Introduction

1

2

Analytical Wall Function Treatments

2

2.1

UMIST-A Scheme

3

2.2

UMIST N Scheme

5

3

Wall Function approaches using two equation model

7

4

Wall Function for oscillating flow

10

4.1

Periodic Channel flow- UMIST-N Approach

10

4.2

Oscillating Turbulent flow- UMIST-A Approach

12

5

Direct Numerical Simulation

13

6

Conclusions

15

Bibliography

16

Nomenclature

k

Turbulent kinetic energy

?

Dissipation rate of k

�"

specific dissipation rate of k

C?

turbulent model constant

Cl

log-law constant

Nu

Nusselt number

Re

Reynolds number

ReDh

Reynolds number based on hydraulic diameter

Reτ

Reynolds number based on frictional velocity

U+

Nondimensional velocity

y+

Dimensionless distance

?

Dynamic viscosity

ρ

Density of the fluid

?t

Eddy or turbulent viscosity

?t*

Dimensionless eddy or turbulent viscosity

y*

Dimensionless distance to the nearest wall

Bo

Buoyancy parameter

τw

Wall shear stress

Pk

Production of k

Lx

Length in x direction

Ly

Length in y direction

Lz

Length in z direction

uτ*

Instantaneous velocity fluctuation

T

Wall temperature

1. Introduction

The Computational Fluid Dynamics (CFD) is the computer solutions of a number of fundamental equations to predict pretty much any type of fluid motion. CFD may sometimes be used to investigate the physical behaviour of flow, but the common goal is to provide an input to engineering analysis and design. The CFD has developed very rapidly in last two decades because of the advancement of technology and increase of computing processor capacity of computer. Considering the development of technology, Researchers focuses toward more complex calculations. With the high-speed supercomputers, better solutions of the complicated flow can be achieved.

Now a days, CFD and experimental analysis are extensively used as a design tool in most of the industries. The position of CFD is currently behind the experimental analysis due to the fact that CFD does not produce absolute results. The reason for this is that the numerical methods, which govern the solutions in a CFD problem, rely on several modelling assumptions that may not have been validated to a satisfactory level. However, CFD presently offers itself as a powerful design tool and even more so in the future because of: firstly, dangerous or expensive trial and error experiments whereas CFD reduce the cost and time that usually require for the experiment. Secondly, the numerical schemes and physical models that are the building blocks of CFD are continually improving. Moreover, it can be able to produce reliable and accurate results on one particular design compare to the exact one. For this reason, most of the researchers prefer CFD rather than Experimental analysis.

The CFD analysis starts with dividing or discretizing the geometry to be modelled into usually a large number of small computational cells. Discretization is the method of approximating the differential governing equations by a system of algebraic equations for the variables at some set of discrete locations in space and time. The discrete locations are normally referred to as the "grid" or the "mesh". This analysis is iteratively solving equations for the conservation of mass, heat, and momentum, and the dissipation of motion into turbulence. The additional effects mentioned above are added into the solution procedure if needed. The results are often presented graphically to assist in visualizing what's actually happening in the flow domain

Most flows of relevance to engineering applications are clearly turbulent, and most CFD codes are designed to compute and predicting the behaviour of turbulent flow. Most commonly, this is accomplished by the adoption of a turbulence model that makes use of certain assumptions to approximate flow behaviour using more tractable equations. This precludes the possibility of predicting every feature of the flow, but the parameters providing the greatest engineering relevance may be estimated at a greatly reduced computational cost.

CFD approximations of turbulent flow near the walls generally employ one of two broad strategies to resolve the very influential, complex, but thin near-wall viscosity-affected sub-layer. One method uses a fine numerical mesh and a turbulence model which take viscous influences into account; the other use a “wall functions’’. However, because of the complexity of the three-dimensional time-dependent flows in engines, the use of low-Reynolds number turbulence models to resolve the near wall regions is computationally very expensive. Thus wall functions are applied in these regions, coupled to a high Reynolds number turbulence model for the bulk flow. The functions used apply the log-law for velocity variation across the near wall regions derived from simple wall shear flows. However, because of the oscillating nature of the flow in engines due to the piston motion, it is not clear that such wall functions are appropriate, as the near-wall flow will also vary in an unsteady fashion.

The aim of the literature survey of this project is to study the wall function to an oscillating turbulent flow imposed on top of a steady through flow and analyse earlier works on wall function near the wall. And hence, ultimately, gather the knowledge to develop an alternative analytical wall function that does properly take the effects of oscillation into account. The scope of the literature survey is to cover the study of different turbulence modelling approaches i.e. two equation models (k-ε, k-ω), Direct Numerical Simulation (DNS) etc. However, major focus of the research would be the near wall treatment of the turbulent flow for instances, wall function, sub-grid approaches (UMIST-A and UMIST-N). Moreover, previous works of the wall function for oscillating turbulent flow will also be covered by this survey.

2. Analytical Wall function Treatments

Two separate projects have done on developing more general wall-functions at UMIST. The pathways of schemes are different with each other. One is UMIST�"A which is an analytical scheme, (Unified Methodology for Integrated Sub-layer Transport�"Analytical), the other is UMIST�"N, based on numerical approach [1].

2.1 UMIST �"A Scheme:

An approach was developed at UMIST by Craft et al. [2] takes into account the pressure gradients, buoyancy forces, and changes in the thickness of viscous sub-layer among other parameters by the analytical integration of the momentum and energy equations. This approach is known as the Unified Modelling through Integrated Sub-layer Treatment �" an Analytical approach (UMIST-A). This method reduced the computing cost by one or two order of magnitudes compared to a standard low-Reynolds-number approach.

UMIST-A scheme was designed for the safety studies specially convection (forced, mixed or natural) flow on near vertical surfaces, strong variations of molecular transport properties across the VSL, and laminarization, i.e. a marked thickening of the VSL in buoyancy-aided mixed convection.

Figure1 Viscosity distribution assumed over near wall cell.

Turbulent viscosity (b) Molecular viscosity [1]

The prescribed viscosity distribution was taken as the form for .The coefficients and Cl are 0.09 and 2.55 as conventional 1-equation turbulence model and .It was proved that evaluating k in the at the sub-layer interface, yv much less stable than adopting the nodal value and also lead to greater dependence on the size of the near-wall cell because to extrapolate values to yv need to use the information from the wall than yP. However, the variation of molecular properties across the sub-layer due to temperature changes also taken into consideration. The paper showed another area where it was felt appropriate to improve current practice was in the prescription of the kinetic energy dissipation rate next to the wall, Chieng and Launder [3] had approximated the exact result of Jones and Launder [4]

To make the treatment sensitive to laminarization, Researchers made two choices, namely appropriate parameters to use as detector and operand. After extensive testing, it was claimed that the ratio of the shear stress between the wall and the edge of the sub-layer, λ was the best detector. In case operand, a direct choice was adopted: the mean level of dissipation rate over the near-wall control volume was adjusted by a weighting function F that in turn was a function of λ.

UMIST-A was tested for two fully developed pipe flow with different near wall cell size. However, at a lower Reynolds number of 6753 the experimental data of Kudva and Sesonske [5] lie above the supposedly ‘universal’ log-law, as do predictions with the low-Reynolds-number (LRN) k�"ε model of Launder and Sharma [6]. Crucially, current wall function results also agree with the data. Moreover, adding to most such schemes, it could be noted that the UMIST�"A scheme results show scarcely any sensitivity to the size of the near-wall cell. The paper shows the Variation of Nusselt number in mixed-convection up-flow in a vertical tube along with UMIST�"A treatment, UMIST�"A without F(λ), experimental data of Li [7], Dittus�"Boelter, LRN calculation in a graph. Result shows, that Wall function did agree with the experimental data. One run was done by UMIST-A without F(λ) results 20% increase in Nusselt number.

This scheme was tested for the Nusselt number variation along a heated annulus wall at a low buoyancy parameter Bo = 0.22 and 2.89 .and compared with experiments of Jackson et al. [8]; LRN: low-Re model calculations; standard wall function predictions; Analytical Wall function. It was found that at a low level of buoyancy the full low-Reynolds-number model of Launder and Sharma [6] and both wall function treatments broadly capture the small increase in Nusselt number. Nevertheless, at higher values of the buoyancy parameter, both the UMIST�"A scheme and the full low-Reynolds-number do capture the now substantial increase in Nusselt number found in the experiments, whereas the standard wall functions do not.

The overall features of the UMIST- A are:

Convective transport (conventional wall treatment ignore it in many cases) expressed by a simplified form.

In a case where buoyancy is an important fact, buoyant force of the near wall cell can be obtained by integrating a fit to the analytical temperature profile over the cell rather than the force purely on the temperature at the near wall cell itself.

If the viscous sub-layer thickness exceeds the cell thickness, yn (as it may do in limited regions if a structured grid is adopted), although a further reformulation is essentially needed, but the identical principle can still be used.

In conclusion to this paper, UMIST-A is based on the analytical solution of simplified near-wall momentum and temperature equations, accounting for pressure gradients and other force fields such as buoyancy. This scheme approaches have been applied to a range of flows in which standard log-law based wall functions are known to perform badly. The overall methodology of UMIST�"A, is closer to that employed in standard wall-functions, thus possibly making it the easier scheme to incorporate into existing flow solvers.

2.2 UMIST-N Scheme:

The analytical approach is not enough equipped in case where the velocity profile parallel with the wall undergoes strong skewing across the sub-layer, for instances, non-orthogonal impingement of flow on a bank of heat-exchanger tubes, for flow with strong streamline curvature where the strain stress relation provided by eddy-viscosity model does not depict the turbulence generation process etc. Considering above difficulties, a numerical scheme, UMIST�"N, has been developed by Gant [9] and Craft et al. [10]. The scheme is much more similar to the low-Reynolds-number models in that the wall-function cell is itself sub-divided into, typically, 30 thin slices as shown in figure 2. With suitable simplifications, the mean flow and turbulence differential equations are solved effectively as a one-dimensional problem across this fine grid. Which ultimately generate data required as of wall function quantities such as wall shear stress, averaged source term etc. and supply appropriate boundary conditions for the whole-field solution carried out on the primary grid.

Figure 2. Treatment of near wall control volume in UMIST- N [1]

There were several simplifications made by the Gant [9] in order to reduction in computer time that one looks form the wall function.

The pressure gradient parallel to the wall is considered uniform across all the sub-grids, equal to the pressure gradient across the near-wall cell of the primary grid.

The velocity component normal to the wall can be obtained by continuity rather than by solving the momentum equation normal to the wall.

The paper shows that two alternative turbulence models have been employed within the above numerical treatment. One is the LRN k�"ε model of Launder and Sharma [6] and other is cubic non-linear eddy viscosity model of Craft et al. [11]. However, result shows that the cubic terms in the latter model make it far more sensitive to streamline curvature than a linear EVM.

Figure 3. Nusselt number variation on a flat plate beneath an axisymmetric impinging jet. (a) Linear k�"ε EVM. (b) Cubic non-linear EVM (Craft et al [11]). Symbols: experiments of Baughn et al [12]; heavy line: full LRN treatment; other lines: UMIST�"N, with different near-wall cell sizes [1]

The paper shows the first test case is for a turbulent jet impinging orthogonally onto a uniformly heated flat plate. The jet discharges from a long smooth pipe whose exit is four diameters above the plate. Figure 3 Shows the resultant variation of Nusselt number over the plate from the stagnation point (r=0) outwards. As is seen in Figure 3(a), the computed Nusselt number at the stagnation point is more than twice the measured value. However, the computation results of wall function and full LRN are close enough. Again same test was done for the non-linear EVM (Craft et al [11]) instead of Linear k�"ε EVM. Now the results from computation and experiment are much closer. The only major difference between the UMIST�"N results and those of the complete LRN treatment is in the computer time required. It is suggested from the paper that the wall-function result with the same grid density takes less than one eighth of the time required for the complete low-Reynolds-number model.

Another test example considered heat transfer from a mildly heated disc spinning about its own axis. The disc’s rotation induces a radially outward motion that peaks outside the VSL, the tangential velocity increases rapidly.

Figure 4. Radial velocity profile for spinning disc in wall-layer coordinates. Solid line: LRN calculation; broken line with symbols: UMIST�"N; chain line: log-law; other lines: standard wall-function treatments [1].

Result shows that the induced radial velocity predicted with UMIST-N (using the linear EVM of Launder and Sharma [6]) agrees very closely with the results of the corresponding LRN computation. However, plot of integral Nusselt numbers shows the negligible difference among the treatments except the computation time. Eventually, the complete low-Reynolds-number model required thirteen times more computation time than UMIST�"N.

In closing to this paper, UMIST�"N is based on a local one-dimensional numerical solution of the governing equations. It is clearly perform badly for the standard log-law based wall functions. This scheme is, in principle, the more general approach, but is inevitably more computationally expensive terms of CPU time and storage requirements.

3. Wall function approaches using two equation model

The paper titled ‘‘The analysis of wall function approaches using two equations turbulence model’’ was done by the Pe´rez-Segarra et al [13].The objective of this paper was to obtain tools which were able to approximate both fluid flow and heat transfer with less time and storage consumption, reduction of grid-sensitivity and a relatively good accuracy. However, considering the scope of literature review, fluid flow sections are going to consider only.

Authors focused on WF approaches, using for inner nodes different high Reynolds number two-equation models, depending on the dissipative variable (k�"ε or k�"ω treatments). Several assumptions and numerical corrections were taken into consideration for near wall treatments. The governing equations of the time-averaged continuity, momentum and energy of the fluid flow were considered. A low- Reynolds �"number model was tested and compared to the Ince�"Launder [14] k�"ε LRN model and with Wilcox [15] k�"ω LRN model in order to find out the limitation and improvements of the numerical and computational phenomena of the wall functions. However, the governing equations were solved using finite volume techniques over a staggered discretization. A structured grid, fully implicit time integration and SIMPLE method were employed to the numerical approach. Diffusion and convection terms were calculated with the central different and UPWIND and SMART scheme for the channel flow and UPWIND and Power law for backward flow respectively. All computations were carried out using a pseudo-transient iterative algorithm, applying the biggest time-step which ensures convergence. Post processing was done using generalized Richardson extrapolation for h-refinement studies and the grid convergence index (GCI) method due to verify the numerical results provided by the continuous (LRN) models analysed in his work.

Different cases were tested shown in this paper. For instances, A channel flow at moderate Reynolds number (Reτ =590, ReDh ≈ 43,000), A backward facing step flow (BFS) at Reynolds number of ReH= 37,500. The inlet conditions used for test cases were identical. For channel flow, , in addition to a turbulence intensity of 7% (kin= 0.072 ) and a specific dissipation rate and a fixed inlet temperature of 300K.

8H

(a) (b)

Figure 5. (a) Channel and (b) backward facing step test cases [13]

Channel Flow case:

The results were plotted as U+ against y+ for different mashes and eventually compared with the DNS data by Moser et al. [16]. The computation results of wall function showed that standard approach (SWF1: wall function treatment proposed bay Launder [17] and k�" ϵ modification (SWF2: modification of SWF1) were failed when near-wall node was placed at viscous regions. Though, improvement achieved when near-wall node was place at buffer layer and sublayer, above two prediction was not appropriate below y+<10. However, results from k-ω models were acceptable manner. Meanwhile, plot of k+ against y+ showed that turbulent kinetic energy predictions was comparatively good than the DNS data, especially with the WWF2 treatment. It was surprisingly found that turbulent suddenly increased at near-wall nodes placed between 5 < y+ < 20 for the second wall cell for WWF1 (ω wall function 1), computations. Overall, All treatments agreed with DNS data by Moser et al [16] for near-wall cells placed at logarithmic region since they are using pure log-law for y+ > 30. On the other hand, LRN computation was done with the finest mesh. Due to finer grid approach, it was found the CPU time consumption was much higher than the WF approach, though the accuracy was more or less the same for all treatments.

Backward Facing step (BFS) Flow Case:

This type of flow require complex and related to different wall mesh grid. However, for WF treatments, seven different meshes were used to test grid-sensitivity of models whereas for LRN, it was five. It was found that using coarsest meshes, WWF1 and WWF2 treatments present a better behaviour than SWF1 and SWF2 approaches. In case of Wall function computation, If meshes refined, k�"ω approaches tend to the asymptotic results provided by the high Reynolds number k�"ω model by Wilcox et al.[18] whereas k�"ϵ predictions were becoming less accurate. However, LRN model provides different result. In order to achieve asymptotic results, near-wall cells placed at y+= 0.5 in both walls (development channel and downstream of the step). This is the main reason to explain the extremely dense mesh used and the high CPU time needed for both simulations.

In conclusion of this paper it can be said that, researcher proposed different strategies for the improvement of the grid sensitivity for standard SWF1 treatment [17] which were based on the enhancing ϵ profile using a two-layer integration (SWF2 approach); changing dissipative variable from ϵ to ω and also from k�" ϵ (Launder, D. Spalding [19]) to k�" ω (Wilcox [18]) for HRN model would be used for inner nodes (WWF1 approach) and, finally, designing a blending function to reduce grid-sensitivity on an ω platform (WWF2 approach).

4. Wall Functions for Oscillating Flow

Many Researchers were seeking for resolving near wall region by wall function for oscillation flow. Among of which UMIST-N and UMIST-A wall functions for periodic flow application will be considered only.

Periodic flow �" the UMIST-N Approach

UMIST-N near wall treatment applied to periodic channel flow was done by Richards [20]. Richards used a parabolic coded PASSABLE (PArabolic Solution Scheme Applied to Boundary Layer Equations) for the simulation. This code is first developed by Addad [21] in his PhD thesis at UMIST and later modified by Richards for periodic flow.

Channel flow with walls that were a distance 2δ apart was used for the investigation. Flow was driven by pressure gradient in the x direction. The problem was solved from y=0 at wall to y=δ at the halfway of the channel. Flow was considered independent of z direction and x direction as well. RANS equations for continuity and conservation of momentum were used for the channel flow. However, governing equations of the flow were discretised by the using of finite volume method. The actual grid used within his work (figure 6) contained 98, 18 and 50 nodes for the for low- Reynolds- number, high-Reynolds-Number, and sub-grid approach respectively.

(b) (c)

Figure 6. Grid for (a) Low-Reynolds-number (b) High-Reynolds-number and (c) sub-grid mesh by Richards [20]

However, Richards used separate under-relaxation factor for main and sub-grid as follows:

(b)

Figure 7. Under-relaxation factor for (a) main grid and (b) sub-grid [20]

This paper shows that, he used four configurations of turbulence models for his simulation work. Firstly, named low-Re k-ε where a low-Reynolds-number standard k-ε model was employed without using sub-grid. Second case was called k-ε with log law which was similar to the first configuration but now a high-Reynolds-number standard k-ε model was used with the log-law boundary condition at near the wall. However, Richards third approach was to use sub-grid k-ε. Now, a high-Reynolds-number standard k-ε model was used in the main grid with the near-wall boundary condition derived from the sub-grid solution and the sub-grid employed the low-Reynolds-number standard k-ε model. Finally, sub-grid k-ω configuration was used for the periodic flow. He employed the 1988 k-ω model in the main grid with the near-wall boundary condition that derived from the sub-grid solution. Same model also used for the sub-gird.

The tests were run over major two cases. These were steady channel flow and periodic channel flow. Furthermore, periodic channel flow was driven from periodically variable pressure gradient and periodically variable bulk flow rate. Results shows in the paper that the periodic cases agreed well with the pressure gradient fluctuations and bulk flow variation of the DNS study of Kawamura & Homma [22]. It was clearly shown that the UMIST-N was successfully applied for the computation of periodic flow.

The results for the log-law of wall, k-ε sub-grid and low-Reynolds-number k-ε configurations were similar for the steady case. However, in the unsteady case the sub-grid approach was given better results than the log-law and had the similarities to the low-Reynolds number approach. The main advantage of the UM IST-N approach was its ability to achieve to lower refinement the near-wall grid compared to the low-Reynolds number grid without hindering the main grid. The k-ω on the other hand had acceptable results overall but the performance was significantly influenced due to the questionable behaviour between the main and sub-grid. This was due to ω tending to infinity at the wall and the averaged production of ω at the sub-grid could not be obtained. Thus, Richards [24] was come into conclusion not to apply UMIST-N approach to a k- ω model.

In conclusion of this paper, Results obtained using the UMIST-N sub grid model are considerably more accurate in terms of the near-wall variations in mean velocity and turbulence parameters than those given by the low-low. However, this approach requires two separate grids and there are some uncertainties associated with the interface between them.

4.2 Oscillating Flow- UMIST-A Approach

UMIST-A analytical wall function approach was investigated by Bauly [23]. The aim of the project was to see whether this would give the results of equal accuracy to those obtained with UMIST-N. Flow considered in this research was 1-D channel flow. The flow conditions were exactly the same as those tested by Richards [20]. The simulation was run by the code called PASSABLE (Parabolic Solution Scheme Applied to Boundary Layer Equations). This originally modified for the UMIST-N periodic flow. In addition to this code, UMIST-A wall function was inserted to it.

Bauly used the same grids and wall boundaries on k�"ε and k-ω as Richards [20] did use in UMIST-N periodic flow approach. In addition to the UMIST-A, he added source terms to the model. For instances, - was added to the U momentum equation, was added to the K transport equation. However, Bauly used the same under relaxation factor for the main and sub-grid as Richards except for ε=1.01 instead of 1.0 in the main grid.

Paper shows that, the cases which tested by Richards [20] were re-run in order to make sure the code works properly. Bauly [23] only considered sub-grid models. . The code was run for both steady and unsteady flow situation. For steady case, the values of U+ and k+ were plotted against the non-dimensional distance y/δ and compared with the DNS results by Kim et al [24]. The k�"ε, 1988 k-ω, and high-Reynolds �"number main grid with log-law wall function produced same results as Richards but results from low-Reynolds-number configuration was unacceptable. Though, UMIST-A wall function underpredicts the values of U+, the results were close to the DNS and thus acceptable. Moreover, the values of k+ also underpredicts for UMIST-A at near the wall because of the higher dissipation rate of k. Result of U+ and τw+ were plotted against period in the case of periodic flow and found underpredicted than the DNS results by Kawamura et al [22]. This was because of the pressure gradient of U momentum equation of the wall shear stress not only drives the flow but also acts to overcome the wall shear stress. It was shown in the paper that, results of U+ for UMIST-A are less accurate away from the wall because of the faulty calculation of turbulent shear stress. Results from the unsteady case shows that UMIST-A wall function produced significant difference in production of k+ than the other models.

In summary to this paper, the application of UMIST-A wall function to the oscillating flow produced poor results especially for the periodic flow. The reason was primarily noted by the Researcher that the UMIST-A was originally developed for the heat transfer type of flow whereas current investigation carried out on the isothermal flow driven by the pressure gradient. Moreover, turbulent kinetic energy results were extremely poor for the periodic flow. It was assumed that reversing of velocity gradient influenced by the backward flow in an oscillating flow that gave the wrong calculation of the turbulent kinetic energy. However, it can be noted that the additional terms added to the standard wall functions did not adequately express the effects of oscillations.

5. Direct Numerical Simulation (DNS)

Direct Numerical Solution (DNS) solves the Navier-Stokes equation numerically without the use of any turbulence models. Thus, DNS solves for the whole range of spatial and temporal scales of turbulence from the smallest to the largest scale. To do this a very fine mesh is required. The number of nodes required for a DNS simulation is proportional to Re9/4, where Re is the Reynolds number. Hence, DNS simulations are very computationally expensive and are usually used to simulate flows of simple geometries and at low Reynolds number. DNS is highly accurate and its data is used to validate turbulence models.

The paper considered here is DNS of an unsteady turbulent channel flow driven by a temporary sinusoidal pressure gradient by Kawamura et al [22]. Most of the DNS analysis assumed the flow is steady but here, a temporally sinusoidal pressure gradient was imposed on the turbulent channel flow. The computational domain specified as Lx=25.6δ, Ly= 2δ and Lz= 3.2δ. the basic continuity and momentum equation of the channel was used where sinusoidal pressure gradient was defined as: , where variables were normalized by channel half width δ and the friction velocity uτs for the steady state. The amplitude A was set to 6.0. The time of one cycle was selected so that an unsteady effect appears while flow reversal would absent, hence, T+ was 6.0. Reynolds number, Reτs =180 was taken for the situation at A=0. However, a fully developed in stream-wise direction flow was assumed and periodic boundary conditions were applied to the x and z directions.

Results from the simulation were plotted in order to get the mean velocity profiles, turbulent kinetic energy, visualization of the turbulence structure and two point correlation. The mean velocity was averaged over the channel section and turbulence quantities were within each phase over 15 cycles. It was shown that, mean velocity increased when pressure gradient positive and decreased for the negative gradient. The results of mean velocity which normalized by the instantaneous uτ* showed that profile were flattened at the centre of the channel during acceleration phase and were more peaked in the deceleration part. Mean velocity was re-plotted in the semi-logarithmic coordinates and was found mostly larger than then well-known steady stale logarithmic distribution and a kind of laminarization observed. However, friction coefficients was plotted against the instantaneous Reynolds number and fond larger than the steady state correlation in phase 1-2 whereas it was smaller in most of the period. The dependence of wall shear stress on Reynolds number was more declarable than the steady case.

The normalized turbulence kinetic energy distribution was plotted across the channel. The results were pretty different from the steady state. It was larger than the steady state when pressure gradient decreasing and smaller when increasing. Moreover, peak of the turbulent energy for both acceleration and deceleration periods were increased away from the wall with compare to the steady state. Authors also normalized the kinetic energy by instantaneous friction velocity and showed that it was larger than the steady state for deceleration period whereas smaller for acceleration. On the other hand, total shear stress normalized by instantaneous friction velocity in order to depict the difference between acceleration and deceleration period. Total values of shear stress for acceleration were decreased while increased in deceleration period. The budget of the turbulent kinetic energy was plotted for the acceleration period (phase 6 and 8) and deceleration period (phase 12 and 14). In the acceleration section, the peak of production was quite smaller than the steady state value of 0.25 because of Reynolds stress decrease at this point. But for the deceleration period, the pressure gradient term becomes negative and thus, shear stress increased significantly. As a result of this effect, the turbulent kinetic energy and production also increased.

The plot of high and low speed streaks in acceleration and deceleration period were shown in paper. In case of acceleration, high and low speed strakes long enough so that they were connected from upstream to the downstream boundaries. However for the deceleration case, the structures were broken into several complex segments and elongated in the span-wise direction. Moreover, in the central region, many ramped structures were appeared. In the streamwise vorticity fluctuations, it was elongated spanwise for deceleration and there was no elongated streaky structure observed for the acceleration period. Two point correlations in x direction was used by the Authors to examine the adequacy of the computational domain. They come into conclusion that, streamwise domain was not large enough to apply the periodic boundary conditions for some of the phases.

6. Conclusions

A critical study on various turbulence models have been carried out for gathering the basic concepts of the computation. In CFD it is more challenge to treat where the viscous effects become important. It is found from the studies that using wall function to resolve near wall region is relatively less expensive rather than the low Reynolds number turbulence number. The famous k and ? wall function treatments are widely used in industrial CFD simulation. Study shows that this model gives less accurate result in near wall region. That is why; k-�" model has been introduced, which was reported to perform better than ? based schemes in the near wall adverse pressure gradient flow. However, recently, more advanced schemes (UMIST-A and UMIST-N), which remove the assumption of log-law have been developed. It’s been studied that analytic integration of mean momentum equation for UMIST-A can be used to estimate the wall shear stress and the production of k. However, the second approach (UMIST-N: numerical wall function) involves putting a 1-D ‘sub-grid’ across each main near wall cell. This scheme basically solves the transport equations for each sub-grid and quantities such as wall shear stress, production and dissipation were computed from these local solutions. On the other hand, solving near wall region in case of periodic flow is more challenging than simple shear flow. Study shows that UMIST-N and UMIST-A schemes were adopted and developed for the case of oscillating flow. The former one has successfully implemented to the period flow but latter one not to be the up to mark. Finally, DNS of unsteady flow has also been studied. The results of this method are more accurate and most of the papers studied here compare the results from turbulence models with that of DNS one.