Simulation of particle trajectories in turbulent flow over a backward-facing step

 

 

Iain Barton

Lecturer in Fluid Mechanics, Department of Mechanical Engineering, Brunei University, Uxbridge, Middlesex, UB8 3PH, United Kingdom

 

 

---

ABSTRACT

The stochastic Monte Carlo scheme is an excellent method to investigate turbulent particle dispersion albeit with a high computational cost. It is a popular method used in conjunction with the k-ϵ turbulence model. The traditional stochastic technigue tends to fail for simulations of heavy particles probably because the crossing trajectories effect, the inertia effect, and the continuity effect are ignored. Models attempting to incorporate two of these effects have been included into the present stochastic technique. The majority of the studies using the stochastic Monte Carlo scheme have considered axial-symmetric jets. The present results, however, investigate a simulation for dilute particulate turbulent flow over a backward-facing step. Satisfactory agreement is achieved with experimental data. The integration of the particle trajectories was achieved using a recently developed predictor-corrector Lagrangian tracking scheme.

---

 

 

Nomenclature

Ccoefficient term

Eturbulent roughness parameter

Ffraction constant

FrFroude number

Grate of production of turbulent kinetic energy

Lnon-dimensionalised length scale

ReReynolds number

Ssource term

StkStokes number

Tnon-dimensionalised time

Utime-averaged velocity component

ddiameter

gacceleration of gravity

kturbulent kinetic energy

l actual length scale

ppressure

u'velocity function

 

Greek

∆x, ∆ycell lengths in x and y co-ordinates

∆ttime-step

Γ diffusivity

αconstant term for eddy's length

ßconstant term for eddy's time

εdissipation rate of kinetic energy

κ von Karman constant

μviscosity

pdensity

σnon-dimensionalised diffusivity constant

τshear stress

 

Superscript

+ non-dimensionalised turbulent velocity and length scales

ntime level n

ppredicted level

 

Subscript

Bboundary point

aaerodynamic relaxation time

e eddy

f fluid

i, jtensor flexing

tturbulent

pparticle

rctraverse time

wwall

w, ww, wwwwest, west-west, 3 × west

 

Introduction

Two-phase flows occur in the environment, industry, and the laboratory. Examples range from particulate pollution, silting in rivers, to important industrial processes such as coal and oil-fired furnace combustion. Cyclone separators, spray-systems, and filtration units are other applications where the behaviour of two-phase flow is important and the behaviour of the dispersed phase is crucial for efficient operation. In almost all cases the flow field is turbulent and how the particles behave in a turbulent flow is therefore of great interest.

Fundamentally, there are two approaches to deal with the particle phase. They are the Eulerian and Lagrangian approaches. The Eulerian approach treats the particle phase as a continuum.^1,2 The Lagrangian approach predicts individual 'representative' particle trajectories.^3,4

It is probably a fair comment to say that the Lagrangian approach is favoured by the engineering community for industrial applications. The advantages and disadvantage of the two approaches have been discussed previously.⁵ A brief review of the Lagrangian approach is given below because it is the approach adopted in the present study. An early industrial example was investigated by Lockwood et al.⁴ Lockwood's study predicted the flow for a hypothetical axial-symmetric coal-fired cylindrical furnace incorporating turbulence modelling and thermal effects. This work was developed to make comparisons with an actual furnace unit.⁶ Subsequent industrial examples in principle are very similar to these early examples. Gosman and Ioannides⁷ predicted a similar furnace design using the stochastic technique which is discussed below. A full three-dimensional front wall coal-fired furnace has been predicted⁸ with only limited agreement with the experimental data which had a large scope of error. Generally, however, numerical predictions investigating combustion design do not attempt to model all the physical features of the flow and only consider, for instance, isothermal particle-laden flow in order to get a rough guide of the fluid-flow behaviour and particle transport.^9,10 Other numerical industrial studies using the Lagrangian approach include applications such as cyclone separation,¹¹ and cooling tower simulation.^3,12

Besides industrial applications the Lagrangian approach can be an excellent way to investigate particle dispersion in a turbulent flow using a Monte Carlo stochastic scheme, albeit with a high computational cost. The stochastic scheme simulates the interaction between phases basically by generating estimated values of the eddy's lifetime, size, and velocity perturbation. The particle experiences the velocity perturbation until the eddy 'dies' or the particle crosses the eddy, after which a new eddy lifetime and velocity perturbation are generated and the procedure is repeated. An early example¹³ applied the basic principles of the standard stochastic scheme to time-dependent spray simulations using some significant simplifications. A popular method, established by the work of Gosman and Ioannides,⁷ is to use the stochastic technique with the two-equation turbulence model of Launder and Spalding¹⁴ where the life-time and size of the eddy as well as the velocity perturbation is determined from the turbulent flow properties.^15-19 Recently, Lu et al.²⁰ have reviewed various improvements to the standard stochastle technique and recommended further alternations.

The flow problem considered in this study is planar turbulent flow over a backward-facing step. The experimental study used for the simulation has a very dilute doping of particles. Therefore we can make an additional simplification that the flow pattern is essentially uncoupled. It is one of the simplest forms of turbulent flow separation because the separation is caused by the change in geometry and the geometry is simple. In addition, this type of flow is used in a variety of ways in various applications, usually to create a recirculation region or a sudden change in pressure. However, our understanding of this flow is still incomplete. The flow configuration has been used considerably to develop and test turbulence models, in particular the k-ϵ model.^21-30 It has not, however, been used to any great extent in the study of turbulent dispersion.

The present study develops the model of Gosman and Ioannides⁷ for predicting the behaviour of particles in a turbulent flow. The model is applied to the application of dilute particulate turbulent flow over a backward-facing step.

 

Formulation of the problem

Governing equations for turbulent flow

The governing equations for turbulent flow can be modelled assuming that the flow properties can be represented by a mean and fluctuating component and in planar form. We make a further assumption that, because there is a dilute presence of the particles, this means there is only an insignificant amount of momentum exchange between phases. This gives the following governing equations using Reynolds averaging: Continuity:

The above equations use tensor-indexing, where ρand μ, are the density and viscosity of the fluid. The time-averaged velocity components are given by U_i, fluctuation terms are expressed as u'. These fluctuation terms lead to extra stress terms, They are called the Reynolds stress terms which have to be modelled and solved. The model is based on the mixing length theory of Boussinesq and Prandtl, where the stress terms are modelled by

The model includes the turbulence energy, k, multiplied by the delta function, δ_ij. The term, μ_t, is the 'turbulent viscosity' which is defined as

where the term, C_μ, is an empirical constant, and the turbulence dissipation rate is ϵ. The governing differential equations for k and ϵcan be expressed as:

where k and εare substituted for φand μ_t/σ_k and μ_t/σ_εare substituted for Γ_φ, respectively. The terms, σ_k and σ_ε, are empirical constants. The source terms for k is as:

where G is the rate of production of turbulent kinetic energy it is expressed:

the terms C₁and C₂are further empirical constants. The values of the empirical coefficients used in the present research were taken from Launder and Spalding¹⁴ and are given as C_μ = 0.09, C₁ - 1.44, C₂= 1.92, σ_k = 1.0 and σ_ε= 1.3.

The turbulent wall boundary conditions are more complex than laminar boundary conditions partly because velocity gradients are very steep. If the Reynolds number is low the steep gradient can be resolved using a low Reynolds number turbulence model, in conjunction with a fine grid close to the wall. The present research uses the generally preferred method of wall functions.

The shear stress along the wall, τ_w, is expressed as

where the boundary point B is adjacent to the wall, therefore u_b is the velocity at this point, y_b is the distance from the wall to this point, and so on. Therefore the new terms introduced are the parameter, E, is the 'roughness parameter' and the term, κ, is the von Karman constant. The non-dimensionalised distance, y_b⁺, is found from:

The significance of = 11.6 is that it is the intersection of the linear wall-function and the logarithmic law which is applicable in the main boundary of the turbulent flow. Therefore if < 11.6 the point is assumed to be in the viscous sub-layer.

The dissipation at the interior boundary node is modelled by ignoring convection and diffusion terms and is set to:

The generation term in the turbulence kinetic energy equation is assumed to satisfy

where ∂u/∂y is the gradient normal to the wall boundary.

These terms are incorporated using the approach described previously.³¹ For instance, the dissipation and prediction terms in the turbulence kinetic energy equation are found by integrating across the cell, where ∆x and ∆y are the length scales for the control volume.

This avoids errors caused by the steep gradients close to the wall boundary as well as decoupling the k_B term from ε_B- The term is given by

The generation term, G_b, in the turbulence kinetic energy equation is simply modelled assuming that the shear stress term, τ_w, is found using equation (9) and the gradient term, (∂u/∂y) is equal to u_B/y_B.

Governing equation for particulate motion

In calculating the particle trajectories, the particles are assumed to be spherical and non-rotating. The motion of the particles is assumed to be dominated by drag and gravity. Also, the particles are assumed not to affect the fluid flow. The various forces and magnitudes have been reviewed elsewhere.^32-34 The present study modifies Stokes law (valid for low Re_p) with the coefficient term used³³ which attempts to describe the behaviour of the wake behind the particle.

The equation of motion of a small solid spherical particle is therefore given as, in non-dimensional form:

The u_p and u_fare the velocity vectors for the particle and instantaneous fluid velocity. The equation is in non-dimensionalised form and T is non-dimensionalised time. The non-dimensionalised terms C_df, Stk and Fr are introduced below. The coefficient term, C_df, in the above equation is:

The C_df coefficient is used to modify the drag term for ultra-Stokesian drag. The drag coefficient uses an important particle aerodynamic parameter, Re_p, the particle Reynolds number which describes the velocity-slip between the particle and fluid.

where d_p is the particle's diameter. The equation of motion for a particle has two important dimensionless parameters, namely the Stokes number, Stk, and the Froude number, Fr. They are defined, respectively, as:

where U and L represent the fluid's velocity and length scales. Note that the Stk number uses the term, ρ_p, which is the particle's density. The Froude number is given by

where g is the acceleration due to gravity. At present, the validity of the equation of motion for a solid spherical particle is doubtful partly because some significant terms have been found empirically through experiment. Nevertheless, it should describe the motion of the particle for Re_p < 100 fairly accurately.

Particulate dispersion model

The present particulate dispersion model is based on the Lagrangian Monte Carlo stochastic method of Gosman and Ioannides⁷ and Shuen et al.^15,35 The fundamental principle of the model is that a particle passes through an eddy which is characterised by a velocity perturbation, a lifetime and a size. Nearly all the models assume that the velocity perturbation is randomly generated forming a Gaussian distribution where the standard deviation of a velocity fluctuation is given by

The majority of the models estimate the 'size' of the eddy, l_e, in the following form:

where α is a constant. Kallio and Stock³⁶ state that αvaries from 0.15 to 2.0 by making comparisons with the experimental data. The majority of other studies are in agreement for this range of αand are summarised in Table 1. The table includes the results of Adeniji-Fashola and Chen,³⁷ Mostafa and Mongia,³⁸ Chen and Crowe,³⁹ Lu et al.,⁴⁰ and Modarress et al.⁴¹

 

 

Generally, the 'lifetime' of the eddy, t_e, is calculated from

where ß₁is a constant. The above term is discussed in Shuen et al.^15,35 Alternatively, the lifetime is calculated from the expression

where β₂ is again a constant. This is detailed in Gosman and Ioannides.⁷ In the second formulation of the eddy lifetime uses the randomly generated velocity perturbation |u'| to calculate the eddy lifetime. Assuming isotropic turbulent conditions it follows that:

If this is substituted into equation (23) a value of β₁can be estimated. Various values of β₁are summarised in Table 2. In this study we assume the eddy lifetime is in the form of equation (22).

 

 

When the size and lifetime of the eddy have been estimated the velocity perturbation is applied until the eddy 'dies' or the particle traverses the eddy. The time taken to traverse the eddy, t_rc, can be estimated^38,39 using the expression:

where |u_r| is the slip velocity between the particle and fluid velocity. More accurately, Stokes' law can be integrated and the solution can be inverted. This approach is adopted,^71,15 where the traverse time is given as

where t_a is the 'aerodynamic' particle relaxation time which is given by

If l_e/ (t_a |u_r|) < 1 then the particle is assumed to be trapped inside the eddy until the eddy dies. This formulation for the traverse time is only appropriate if there are no significant body forces, such as gravity, present. If there was a constant body force present the equation of motion can be integrated but not inverted. However, the traverse time can be solved using the Regula Falsi method. Alternatively, the trajectory of the eddy can be calculated and thus the relative displacement of the particle from the centre of the eddy is known, This approach appears to have been adopted.¹⁹ Unfortunately, this approach places further constraints on the model because the time-step has to be small enough to allow integration of the trajectory within the eddy.

The present research assumes mid-range values of eddy size and lifetime in comparison with previous studies using α= 1.5 and ß = 1.5. Therefore for the present research t_e = l_e/k^1/2. The traverse-time is not estimated but is found directly by numerical integration which is appropriate if the body force is relatively simple in form and a predictor-corrector integration is used so the effects of ultra-Stokesian drag can be incorporated without introducing a significant error. The time-step value is estimated from:

where F is a small fraction set to about F = 0.2; u_p and V_p are the velocity components for the particle; the function 'min' selects the minimum value. This method for incorporating body forces into the model is one of the improvements made in the present research to the model of Shuen et al.^15,35 The failure of the model to incorporate body-force was recognised by Graham⁴² and describes the error as 'the inertia effect'. Also, Graham points out that the standard stochastic model has two additional errors which are described as 'the crossing trajectory effect' and 'the continuity effect'. The continuity effect is observed in the experiments where dispersion is stronger in the streamwise direction. Zhang et al.¹⁶ have also observed this effect and assumed non-isotropic turbulence conditions for the particles and set the streamwise velocity perturbation to:

It follows therefore that the standard deviation of the velocity perturbation in the streamwise direction is σ_u = k^l/2 and for other directions σ_v,w = (l/2k)^1/2. This simple model is adopted in the present research where the streamwise direction is taken to be the local fluid velocity vector.

 

Numerical methodology

Solution of the flow-field

The governing equations are solved in primitive form (u,v,p) using the Semi-Implicit Method for Pressure Linked-Equations (SIMPLE) methodology.⁴³ The numerical methodology used in the present study is similar to previous studies.^44,45 Discretization is achieved using the Second-Order Upwind Difference scheme instead of the hybrid differencing scheme⁴⁶ used.⁴⁴ The SOUD scheme was found to be more successful in predicting grid-independent results from previouc; investigations^5,45

Non-uniform grid distributions were used for simulations with clustering near-solid boundaries, especially around the step. Grid independent results are presented for the flow field. The grid dependency of computations was tested using similarly clustered grids with 40×40, 80×60 and 100×80 grid points.

The outlet condition was found by extrapolation. First-order and quadratic extrapolations recommended by Peric³¹ were found to give fairly poor solutions in the outlet region. The outflow boundary condition developed and discussed⁵¹ was applied. In short, the extrapolated velocities at the outlet are calculated using the following fit.

where a₁, a₂, and a₃are the coefficients for the fit. In the extrapolation the four velocity positions upstream of the outflow boundary are used. The velocity position the furthest away from the outflow boundary is used as the datum position for x, and the velocity values at the other three positions are used fcr extrapolation calculation. The ∆x term is the cell length adjacent to the outflow boundary. Therefore, if uniform cells are used near the exit region, the extrapolated velocity value is estimated by (using compass notation),

where u_exit is the extrapolated outlet velocity. The exit velocity is then corrected to ensure that overall flux is conserved. This formulation appears to reduce numerical errors near the outflow boundary, as detailed.⁵

The convergence criterion was assumed to have been met when the average absolute changes in the velocity and pressure fields had reduced five orders of magnitudes.

Integration of the particulate trajectories

The principle of calculating a particle's trajectory in a laminar flow field is well-established and straight-forward. There are examples where the equation of motion for a particle in an oscillating laminar flow field is integrated.^47,48 More recently, Ruetsch and Meiburg⁴⁹ integrated the equation of motion for a particle in idealistic analytical flow-fields and analysed the various forces acting on the particle. A similar study by Barton⁵⁰ integrated the equation of motion for particles in laminar flows over a backward-facing step. The integration of the motion of the particle has a complexity related to Stokes' law. In order to illustrate this point, assume that the motion of the particle is only governed by drag:

where u_f is the velocity of the instantaneous fluid velocity. The equation of motion is therefore a simple first-order differential equation which is difficult to solve because the fluid velocity is a function of position and/or time for all but simple cases.

For most cases equation (32) can be accurately solved using a fourth-order Runge-Kutta technique.⁵¹ The Runge-Kutta scheme has the draw-back that it is limited by the stability requirement that the time-step, ∆t, meets the condition ∆t < 2t_a which is similar to other explicit Eulerian schemes.⁵² Alternatively, a less accurate exponential Lagrangian tracking scheme can be used. The exponential scheme was developed partly to overcome the stability requirement associated with the Runge-Kutta scheme.⁵³ Initially, the study below considers integrating equation (32) using the well-established exponential scheme and develops the scheme to give improved accuracy, as detailed.⁵⁴ The original exponential scheme simply assumes that the fluid velocity remains constant until the next time-step. This allows equation (32) to be integrated to form the solution:

where At is the time-step interval, n + 1 is the final time-level and n is the initial time-level. Equation (33) can then be integrated again to find the position of the particle, x_p.

This approach has been adopted.⁵⁵ However, the displacement is usually found using the following simpler expression:⁵⁶

Equation (33) predicts that particles with a large t_a value tend to give the limit u_pⁿ⁺¹ = u_pⁿ and for small t_a values the equation tends to the limit u_pⁿ⁺¹ = u_fⁿ. The first limit is physically acceptable. The second limit is not acceptable because the particle's velocity lags behind the fluid's velocity by the numerical time-step. A physically acceptable limit would be u_pⁿ⁺¹ = u_fⁿ⁺¹ but the term is unknown until the particle steps forward in both time and displacement. This is the fundamental problem of the exponential scheme which the present study addresses by developing a new predictor-corrector exponential scheme. The scheme is outlined below, with further details.^5,54 The scheme has been successfully applied for the study of particulate coupled laminar flow over a backward-facing step.⁵⁷

The predictor-corrector scheme assumes the fluid velocity varies with time over the preceding and succeeding time-steps with the relationship:

The U_f^np value is set to the most recently predicted value of the fluid velocity. When the above u_f(t) equation is substituted into equation (32) then the general solution for the particle's velocity and position can be expressed as:

By setting A = 0, B = 0 the general solution reduces to original scheme. The limits for the new scheme for small t_a values are u_pⁿ⁺¹ = u_fⁿ⁺¹for the velocity and for the displacement.

So far the predictor-corrector exponential scheme has only been presented for the problem of integrating Stokes Law. Extending the scheme to incorporate additional forces such as gravity is discussed elswhere ⁵⁴ It is felt by the author that a predictor-corrector integration method has to be used when predicting trajectories for the present methodology, because when the particle exits an eddy and enters another eddy, it is possible that a large velocity slip will occur. If this velocity slip is used for the particle Reynolds number then artificially high values will be produced. The alternative option would be to set the particle Reynolds number to the previous time-step producing an artificially low value.

 

Computational results

Turbulent flow over a backward-facing step validation

The turbulence model was validated by considering a standard turbulent flow problem established at the 1980-1981 Stanford Conference.⁵⁸ The geometry has an inlet channel five step-heights long and an expansion number, E = 3/2. The channel is 30 step-heights long overall. The inlet mean velocity, turbulence intensity and turbulence dissipation profiles are specified five step-heights upstream of the sudden expansion and are matched to available experimental data. The inlet Reynolds number is Re = 132 000 where the velocity scale is based on the inlet centre-line velocity and the outlet channel height. Various predictions of the problem using the standard k - εmodel can be found elsewhere,^21-30 and the problem was investigated experimentally by Kim et al.⁵⁹

Using the standard k - εmodel,¹⁴ the flow was predicted to reattach at x/h ≈ 6.2 in the present research which is in agreement with the standard k - εmodel results.²⁹ Further details can be found.⁵ However, this value is an under-prediction relative to the mean experimental reattachment point at x/h ε 7.1.⁶⁰ The underestimate of the recirculation length in comparison with the experimental data is a well-known property of the standard k - εturbulence model. Thangam and Speziale²⁸ have assessed some improvements for the model using the present problem as their main test case.

Turbulent particulate dispersion results

The models discussed in the previous sections are now applied to the problem described previously.^60-61 The problem is illustrated in Figure 1. The flow configuration is a turbulent flow over a backward-facing step with an expansion number of two. Results are presented for the flow Reynolds numbers of Re = 15 000 and Re = 64 000. (The Reynolds number is based on the step-height, h, and the mean inlet centre-line velocity, U_o. The predictions model the experiments by assuming similar inlet conditions which are set 30h upstream of the sudden expansion. The main channel is 25h long. The computational grid uses 100x80 grid points with clustering near the solid boundaries. There is a greater grid density of lines in the main channel compared with the irilet channel.

 

 

The predictions integrate 1 000 particle trajectories per inlet cell, therefore a simulation tracks a total of 40 000 particle trajectories. Similar to an experimental study, the numerical study has to 'sample' or 'detect' a particle entering a volume in the flow, The 'sample' volume's diameter was set to 0.01h. This approach makes a numerical simulation costly because such a high number of particle trajectories have to be integrated to obtain reasonable results.

In the studies,^60,61 experimental data is given for 1 mm oil particles and starch particles which are approximately spherical with an average diameter of 70 μm. Despite the fairly significant differences in sizes and densities the mean U-velocity profile results were found to be fairly similar to the experimental results. In the numerical simulations the profiles were found to be virtually identical; see Figures 2 and 3. Figure 2 shows the predictions and experimental data of the U-velocity profiles for Re = 64 000 using the 1 μm oil particles and Figure 3 shows the U-velocity profiles for Re = 64 000 using the 70 μm starch particles. (Note: the length scale of 2 step-heights on the figures represents a velocity value of U_o). In Figure 4 the U-velocity profiles are shown for Re = 15 000 using the 70 μm starch particles.

Overall, the predictions of U-velocity profiles are in satisfactory agreement with the experimental data. However, considering that the U-velocity profiles essentially demonstrate the accuracy of the k - εmodel, altering the stochastic model should have only a minor effect on the U-velocity profiles. A significant problem with these comparisons is that the numerical predictions are demonstrating the behaviour of the k - εturbulence model as well as the stochastic model and it is not straight-forward to isolate one model from the other.

Unlike the U-velocity profile predictions shown in Figure 2 and 3, where different particles are used for the same Reynolds number, there are fairly significant differences in the predictions of the streamwise turbulence results, U'². (Note that: there is no experimental data for y-component velocities or stress terms). The profiles are shown in Figures 5 and 6 for the predictions and experimental data of the U'²/U_o'² profiles. Figure 5 shows profiles for the Re = 64 000 flow using 1 mm oil particles and Figure 6 shows the profiles for the same Reynolds number using the 70 μmstarch particles. (Note that: a length scale of l½ step-heights represents a U'²/U_o² value of 0.1 on the figures). The agreement with the experimental data is satisfactory, especially for the profiles near x/h = 1 and x/h = 9. There is a slightly better agreement for the 70 μm particle results between the simulations and the experimental data than for the 1 μmparticle results. The simulation using the 70 μm particles have smaller turbulence intensity values in comparison with the simulation using the 1 μmparticles. Both sets of predictions of streamwise turbulence tend to overestimate the experimental data. The streamwise turbulence results for Re = 15 000 using the 70 μm starch particles are shown in Figure 7. Again, the predictions and the experimental data are in fairly acceptable agreement with the predictions: the greatest differences occur for the profiles near x/h = 5 where the predictions tend to over-estimate the experimental data.

 

Concluding remarks

It is important to use a predictor-corrector integration scheme if particle trajectories are predicted moving from one eddy to another.

The stochastic model used in the present study incorporates two developments in comparison with the standard model. The first is incorporating external forces acting on the particle by integrating the movement of the eddy as well as the movement of the particle. The second is modelling the 'continuity effect' by assuming anisotropy, i.e. a stronger velocity perturbation in the local stream-wise direction than in the cross-stream direction.

The stochastic method for turbulent particle dispersion is applied to a simple application of flow over a backward-facing step. In comparison with experimental data, the predictions are similar; differences are probably due to inadequacies in the turbulence model. Nevertheless, the two improvements made to the standard model work, without causing a significant increase in computational cost or problems. It is not possible to argue from the present results whether the improved stochastic method gives better results than the standard model because of the limitations of the turbulence model. However, it is clear that the present model avoids some obvious modelling flaws.

 

References

1. FH Harlow & AA Amsden. Numerical calculation of multiphase fluid flow. J. Comput. Phys., 17, 1975, pp.19-52.         [ Links ]

2. DB Spalding. Numerical computation of multiphase fluid flow and heat transfer. Recent Advances in Numerical Methods in Fluids, Taylor C & Morgan K (eds), 1, 1980, pp.139-168.

3. C Benocci, JM Buchlin & P Weinacht. A prediction method for the air-droplets flow in the inlet section of a natural draught cooling tower. Technical Memorandum 40EAT 8602/CB - J MB - PW/NT, von Karman Institute, Belgium, 1986.

4. FC Lockwood, AP Salooja & SA Syed. A prediction method for coal-fired furnaces. Combust. Flame, 38, 1980, pp.1-15.         [ Links ]

5. IE Barton. Computation of Incompressible Particulate Flows. PhD Thesis, University of Manchester, Manchester, 1995.         [ Links ]

6. AS Abbas, SS Koussa & FC Lockwood. The prediction of the particle laden gas flows. Imperial College report, FS/80/1, 1980.

7. AD Gosman & E Ioannides. Aspects of computer simulation of liquid-fuelled combustors. J. Energy, 7, 1983, pp.482-490.         [ Links ]

8. RK Boyd & JH Kent. Three-dimensional furnaces computer modelling. Twenty-first Symposium (International) on Combustion, The Combustion Institute, 1986, pp.265-274.

9. I Celik. Isothermal prediction of particle and gas flow in a coal-fired reactor. Particulate Sci. Tech., 6, 1988, pp.53-68.         [ Links ]

10. JS Anagnostopoulous & GC Bergeles. Discrete-phase effects on the flow field of a droplet-laden swirling jet with recirculation: A numerical study. Int. J Heat Fluid Flow, 13, 1992, pp.141-151.         [ Links ]

11. F Boysan, WH Ayers & J Swithenbank. A fundamental mathematical modelling approach to cyclone design. Trans Inst. Chem. Engng., 60, 1982, pp.222-230.         [ Links ]

12. P Weinacht & JM Buchlin. A numerical model for gas-droplet flow application to liquid spray and cooling towers. Technical Note 144 D/1982/0238/249, von Karman Institute, Belgium, 1982.

13. JK Dukowicz. A particle-fluid numerical model for liquid sprays. J. Comput. Phys., 35, 1980, pp.225-253.         [ Links ]

14. BE Launder & DB Spalding. The numerical computation of turbulence flows. Comput. Meth. Appl. Mech. Engng., 3, 1974, pp.269-289.         [ Links ]

15. JS Shuen, ASP Solomon, QF Zhang & GM Faeth. Structure of particle-laden jets: Measurements and predictions. AIAA J, 23, 1985, pp.396-404.         [ Links ]

16. QF Zhang, JS Shuen ASP Solomon & GM Faeth. Structure of ducted particle-laden turbulent jets. AAIA J, 23, 1985, pp.1123-1125.         [ Links ]

17. AA Mostafa & HC Mongia. On the interaction of particles and turbulent fluid flow. Int. J Heat Mass Transfer, 31, 1988, pp.2063-2075.         [ Links ]

18. E Mastorakos, J McGuirk & AMKP Taylor. The origin of turbulence acquired by heavy particles in a round, turbulent jet. Part. Part. Syst. Charact., 7, 1990, pp.203-208.         [ Links ]

19. CP Chen, HM Shang &Y Jiang. An efficient pressure-velocity procedure for gas-droplet two-phase flow calculations. Int. J Numer. Meth. Fluids, 15, 1992, pp.233-245.         [ Links ]

20. QQ Lu, JR Fontaine & G. Aubertin. Particle motion in two-dimensional confined turbulent flows. Aerosol Sci. Tech., 17, 1992, pp. 169-185.         [ Links ]

21. C Taylor, CE Thomas & K Morgan. Modelling flow over a backward-facing step using the F.E.M. and the two-equation model of turbulence. Int. J Numer. Meth. Fluids, 1, 1981, pp.295-304.         [ Links ]

22. LP Hackman, GD Raithby & AB Strong. Numerical predictions of flows over backward-facing steps. Int. J Numer. Meth. Fluids, 4, 1984, pp.711-724.         [ Links ]

23. RS Amano & P Goel. Computations of turbulent flow beyond backward-facing steps using Reynolds-stress closure. AIAA J, 23, 1985, pp.1356-1361.         [ Links ]

24. CG Speziale. On non-linear k - l and k - e models of turbulence. J Fluid Mech., 178, 1987, pp.459-475.         [ Links ]

25. CG Speziale h T Ngo. Numerical solution of turbulent flow past a backward-facing step using a nonlinear k - e model. Int. J Engng. Sci., 26, 1988, pp.1099-1112.         [ Links ]

26. JL Sohn. Evaluation of FIDAP on some classical laminar and turbulent benchmarks. Int. J Numer. Meth. Fluids, 8, 1988, pp.1469-1490.         [ Links ]

27. S Thangam & N Hur. A highly-resolved numerical study of turbulent separated flow past a backward-facing step. Int. J Engng. Sci., 29, 1991, pp.607-615.         [ Links ]

28. S Thangam & CG Speziale. Turbulent flow past a backward-facing step: A critical evaluation of two-equation models. AIAA J., 30, 1992, pp.1314-1320.         [ Links ]

29. CG Speziale & S Thangam. Analysis of an RNG based turbulence model for separated flows. Int. J. Engng. Sci., 30, 1992, pp.1379-1388.         [ Links ]

30. A Pentaris, K Nikolados & S Tsangaris. Development of projection and artificial compressibility methodologies using the approximate factorization technique. Int. J Numer. Meth. Fluids, 19, 1994, pp.1013-1038.         [ Links ]

31. M Peric. A finite volume method for the prediction of three-dimensional fluid flow in complex ducts. PhD Thesis, Imperial College, London, 1985.         [ Links ]

32. GB Wallis. One-dimensional Two-phase Flow. New York, McGraw-Hill, 1969.

33. R Clift, JR Grace & ME Weber. Bubbles, Drops and Particles. Academic, New York, 1978.

34. SL Soo. Multiphase Fluid Dynamics. Science Press, Gower Technical, 1982.

35. JS Shuen, LD Chen & GM Faeth. Evaluation of a stochastic model of particle dispersion in a turbulent round jet. AIChE J., 29, 1983, pp.167-170.         [ Links ]

36. GA Kallio & DE Stock. Turbulent particle dispersion: a comparison between Lagrangian and Eulerian modelling approaches. Gas-Solid Flows. JT Jurewicz (ed.), ASME, 1986, pp.23-29.

37. A Adeniji-Fashola & CP Chen. Modelling of confined turbulent fluid-particle flows using Eulerian and Lagrangian schemes. Int. J Heat Mass Transfer, 33, 1990, pp.691-701.         [ Links ]

38. AA Mostafa & JHC Mongia. On the modelling of turbulent evaporating sprays: Eulerian versus Lagrangian approach. Int. J. Heat Mass Transfer, 30, 1987, pp.2583-2593.         [ Links ]

39. PP Chen & C.T. Crowe. On the Monte Carlo method for modelling particle dispersion in turbulence. Gas-Solid Flows. JT Jurewicz (ed.), ASME, 1984, pp.37-41.

40. QQ Lu, JR Fontaine & G Aubertin. Numerical study of the solid particle motion in grid-generated turbulent flows. Int. J Heat Mass Transfer, 36, 1993, pp.79-87.         [ Links ]

41. D Modarress, H Tan & S. Elghobashi. Two-component LDA measurement in a two-phase turbulent jet. AIAA J., 22, 1984, pp.624-630.         [ Links ]

42. DI Graham. Improved eddy interaction models with random length and time scales. Int. J. Multiphase Flow, 24, 1998, pp.335-345.         [ Links ]

43. SV Patankar. Numerical Heat Transfer and Fluid Flow. Hemisphere. Washington DC, 1980.

44. IE Barton. Laminar flow past an enclosed and open backward-facing step. Phys. Fluids A., 6, 1994, pp.4054-4056.         [ Links ]

45. IE Barton. A numerical study of flow over a confined backward facing step. Int. J. Numer. Methods Fluids, 21, 1995, pp.653-665.         [ Links ]

46. DB Spalding. A novel finite difference formulation for differential expressions involving both first and second derivatives. Int. J Numer. Methods Engrng., 4, 1972, pp.551-559.         [ Links ]

47. PR Schöneborn. The interaction between a single particle and an oscillating fluid. Int. J Multiphase Flow, 2, 1975, pp.307-317.         [ Links ]

48. RA Herringe. A study of particle motion induced by two-dimensional liquid oscillations. Int. J Multiphase Flow, 3, 1977, pp.243-253.         [ Links ]

49. GR Ruetsch & E. Meilburg. On the motion of small spherical bubbles in two-dimensional vortical flows. Phys. Fluids. A, 5, 1993, pp.2326-2341.         [ Links ]

50. IE Barton. Computation of particle tracks over a backward-facing step. J. Aerosol Sci., 26, 1995, pp.887-901.         [ Links ]

51. CY Chow. An Introduction to Computational Fluid Mechanics. John Wiley &Sons, 1979.

52. DA Anderson, JC Tannehill & RH Pletcher. Computational Fluid Mechanics and Heat Transfer. Hemisphere, Washington, 1984.

53. JP Boris & ES Oran. Numerical simulation of reactive flow. Elsevier, 1987.

54. IE Barton. Exponential Lagrangian tracking schemes applied to Stokes law. Trans. ASME, J Fluids Eng., 118, 1996, pp.85-89.         [ Links ]

55. RJ Litchford & SM Jeng. Efficient statistical transport model for turbulent particle dispersion in sprays. AIAA J, 29, 1991, pp.1443-1451.         [ Links ]

56. CT Crowe, MF Sharma & DE Stock. The Particle-Source-In Cell (PSI-Cell) model for gas-droplet flows. Trans. ASME, J Fluids Engng, 99, 1977, pp.325-332.         [ Links ]

57. IE Barton. Computation of dilute particulate laminar flow over a backward facing step. Int. J Numer. Methods Fluids, 22, 1996, pp.211-221.         [ Links ]

58. SJ Kline, BJ Cantwell & GM Lilley (eds). Backward Facing Step Flows. Proc. 1980-1981 AFOSR-HTTM Stanford Conference on Complex Turbulent Flows. Stanford University Press, 1981.

59. J Kim, SJ Kline & JP Johnston. Investigation of a reattaching turbulent shear layer: Flow over a baclward-facing step. Trans. ASME. J Fluids Engng., 102, 1980, pp.302-308.         [ Links ]

60. B Ruck & B Makiola. Particle dispersion in a single-sided backward-facing step flow. Int. I. Multiphase Flow, 14, 1988, pp.787-800.         [ Links ]

61. B Ruck & B Makiola. Dynamics of micrometer-sized particles in a separated step flow. Physics of Separated Flows - Numerical, Experimental and Theoretical Aspects. Gersten K (ed.) Vieweg, 1993, pp.56-65.

 

 

Received January 1999
Final version September 1999