Modeling Turbulent Flows with Free-Surface in Unbaffled Agitated

Mar 10, 2006 - Nexia Solutions Limited, Sellafield, Seascale CA20 1PG, U.K.. In the present study, numerical simulations of turbulent flow with free-s...
0 downloads 0 Views 726KB Size
Ind. Eng. Chem. Res. 2006, 45, 2881-2891

2881

Modeling Turbulent Flows with Free-Surface in Unbaffled Agitated Vessels Jennifer N. Haque, Tariq Mahmud,* and Kevin J. Roberts Institute of Particle Science and Engineering, School of Process, EnVironmental and Materials Engineering, The UniVersity of Leeds, Leeds LS2 9JT, U.K.

Dominic Rhodes Nexia Solutions Limited, Sellafield, Seascale CA20 1PG, U.K.

In the present study, numerical simulations of turbulent flow with free-surface vortex in unbaffled vessels agitated by a paddle impeller and a Rushton turbine, which were investigated experimentally by Nagata (John Wiley & Sons: New York, 1974) and Ciofalo et al. (Chem. Eng. Sci. 1996, 51, 3557-3573), respectively, have been carried out. A homogeneous multiphase flow model coupled with a volume-of-fluid (VOF) method for interface capturing has been applied to determine the shapes of the gas-liquid interface and to compute the turbulent flow fields in unbaffled vessels. Turbulence is modeled using the k-/k-ω based shear-stress transport model (Menter, F. R. AIAA J. 1994, 32, 1598-1605) and a second-moment differential Reynoldsstress transport model. Calculations are carried out using the ANSYS CFX-5.7 CFD code (ANSYS: Canonsburg, PA, 2004). Validation of the predictions is effected against the measured free-surface profiles (Nagata, 1974; Ciofalo et al., 1996), and the mean velocity distributions (Nagata, 1974). The predicted liquid surface profiles using the VOF method in conjunction with both turbulence models are generally in good agreement with measurements. As for the mean velocity components, the Reynolds-stress transport model predictions are superior than those obtained using the shear-stress transport model. 1. Introduction Agitated vessels used in the chemical and process industries commonly contain baffle(s) of various designs, which include, for example, four flat vertical plates fixed to the vessel wall with or without a gap and a single beavertail or finger baffle away from the wall, to prevent vortex formation at the surface of the liquid and to enhance mixing. However, vessels without baffles are preferred in a number of applications in the food and dairy, nuclear fuel processing, pharmaceutical, and speciality chemicals industries. Unbaffled agitated vessels are used for crystallization to reduce particle attrition,1 blending of viscous liquids to avoid formation of dead zones,2 solid-liquid mass transfer to achieve high rates,3 fermentation,4 suspension of solid particles, and disengagement of gas bubbles from liquid to reduce foaming.5 The hydrodynamic and mixing characteristics of agitated vessels are significantly influenced by the configuration of the vessel, the agitator type, and the operating conditions, which, in turn, affect the product properties. The increasing importance of product quality in these industry sectors demands improved understanding of the influence of hydrodynamic and mixing characteristics of the vessel. A vast amount of experimental and computational fluid dynamic (CFD) modeling work on the hydrodynamics and mixing in baffled vessels agitated mostly by Rushton turbines, but also with other agitator types such as paddle and pitched blade impellers, has been reported in the literature (see reviews in references 6-8). In contrast, a very limited number of experimental and computational studies have been carried out in unbaffled vessels. In most experimental studies,9-13 the liquid surface was covered with a lid to prevent vortex formation and air entrainment. In these studies, measurements of mean and fluctuating velocities were carried out using LDV and PIV * Author for correspondence. Tel.: +44 (0)113 3432431. Fax: +44 (0)113 3432405. E-mail: [email protected].

systems. There is, however, a dearth of experimental data in unbaffled vessels with a free liquid surface in the literature. Nagata14 reported measurements of mean velocity distributions using a Pitot-static probe and surface profile in such a system agitated by a paddle impeller. Ciofalo et al.15 conducted only measurements of the free-surface profiles at different agitation speeds in an unbaffled vessel with a Rushton turbine. In agitated vessels, strong mean axial and radial velocity components are useful because they provide the flow structure necessary for effective mixing. However, in an unbaffled vessel, the impeller generates strong swirling motion, resulting in relatively weak axial and radial flows which lead to poor mixing. Previous velocity measurements have revealed that the mean tangential velocity is much larger than the mean axial and radial velocities, and the radial distributions of tangential velocity at different heights conform to that of a combined vortex profile consisting of an inner forced-vortex and an outer free-vortex. As for the turbulence structure, high levels of anisotropy exist throughout the vessel in the absence of baffles,12,13 whereas in baffled vessels, the turbulence is anisotropic in the impeller stream16,17 and largely isotropic in the bulk-flow region.12,18,19 In the absence of a lid in contact with the liquid, a vortex is formed on the free liquid surface, and the vortex depth increases with increasing impeller speed. The shape of the free-surface is controlled by the hydrodynamics of the vessel established at a given impeller speed.15 Modeling turbulent flow in an unbaffled vessel is highly challenging because of the complexity of the flow, which is characterized by strong turbulence anisotropy, streamline curvature, and rotation, and the problem is exacerbated in the presence of a free liquid surface. The computational studies2,9,11,20-22 reported in the literature have focused on the vessels fitted with a lid in contact with the liquid surface, thus preventing vortex formation. Among these authors, Alcamo et al.,2 Armenante and co-workers,9,11 and Dong et al.20 compared

10.1021/ie051021a CCC: $33.50 © 2006 American Chemical Society Published on Web 03/10/2006

2882

Ind. Eng. Chem. Res., Vol. 45, No. 8, 2006

the predicted mean velocity components and turbulent kinetic energy with in-house experimental data,9-11,13 while Shekhar and Jayanti21 used the measurements obtained by Dong et al.10 Mununga et al.22 reported only the predictions of gross flow quantities such as power number, flow number, and pumping efficiency, as well as some qualitative comparisons with experimental data.14 In these studies, the turbulence was represented by the isotropic eddy-viscosity based standard highReynolds-number k- model, with the exception of those by Armenante and co-workers9,11 and Shekhar and Jayanti,21 who also used an algebraic Reynolds-stress model (ASM) and a lowReynolds-number k- model, respectively. Alcamo et al.2 adopted only the large-eddy simulation (LES) method. The predicted11,20,21 distributions of mean tangential velocity obtained using the high-Re k- turbulence model were generally in good agreement with measurements away from the impeller shaft in the free-vortex region, but in the region close to the shaft, i.e., in the forced-vortex region, the distributions were significantly lower than those measured. There was a tendency in the predicted tangential velocity distributions to approach forcedvortex (solid-body rotation) profiles across the whole tank, whereas the ASM and the low-Re k- turbulence model captured the correct shape of the tangential velocity profile. Varying levels of agreement between the high-Re k- model predicted and measured mean axial and radial velocities were reported by these authors. The magnitude of the axial velocity is usually significantly smaller than the other two components, which led Armenante et al.11 and Alcamo et al.2 to doubt the reliability of measurements for this quantity. However, Armenante et al.11 reported that the magnitudes of the axial and radial velocities were generally well-predicted by the ASM in the region above the impeller, while the predicted velocities were only in qualitative agreement with measurements below the impeller. On the other hand, Alcamo et al.2 refrained from comparing the LES predicted axial velocity with measurements. The LES predictions were in excellent agreement with the measured mean tangential velocity but were less satisfactory for the radial velocity. Numerical simulations of turbulent flow with free-surface vortex in unbaffled vessels are rare in the literature. The only study of this type of flow in unbaffled vessels agitated by a Rushton turbine and a paddle impeller was carried out by Ciofalo et al.15 using the standard k- and a differential Reynolds-stress turbulence model. The authors developed an iterative solution technique to determine the shape of the vortex by imposing a tentative free-surface profile (usually a flat surface) as an initial guess. The initial surface profile was used to generate a body-fitted computational mesh for the calculation of the flow field and surface pressures, which were then used to compute the new shape of the liquid surface. The procedure was repeated until convergence of the free-surface profile was obtained. The development of this methodology must be regarded as a significant step forward toward modeling complex free-surface flow in unbaffled agitated vessels. However, the approach is not expedient from a computational point of view because it requires generation of new body-fitted meshes as the free-surface profile evolves with successive iterations. Ciofalo et al.15 validated the predicted surface profiles against their own measurements for a Rushton turbine and that obtained from Nagata’s theory for a paddle impeller,14 showing reasonably good agreement. The predicted mean velocity components were compared with the experimental data of Nagata.14 The predicted tangential velocity profiles using the standard k- turbulence model erroneously approached solid-body rotation across the

whole tank. On the other hand, a satisfactory agreement between the predicted and measured flow fields was obtained with the Reynolds-stress model. In this study, numerical simulations of turbulent flows with free-surface in two unbaffled vessels, one agitated by a sixbladed Rushton turbine and the other agitated by a paddle impeller, which were studied experimentally by Ciofalo et al.15 and Nagata,14 respectively, have been carried out using the general purpose proprietary CFD code ANSYS CFX-5.7.23 A homogeneous multiphase and free-surface flow model based on explicit modeling of the gas-liquid interface has been applied to determine the profiles of depressed liquid surface and the flow fields at different agitator speeds. To the best of the authors’ knowledge, such a modeling approach has only been applied to simulate flow with a flat, but wavy, free-surface in a conventional baffled agitated vessel by Serra et al.24 However, no flow field predictions were reported by the authors. In the present study, calculations are performed using a homogeneous flow model coupled with the eddy-viscosity based shear-stress transport model of Menter25 and a differential Reynolds-stress model for turbulence for the entire flow field in the vessel consisting of air and water. The computed freesurface profiles of liquid are validated against the existing experimental data,14,15 and the three components of mean velocity are compared with the measurements reported by Nagata.14 2. Mathematical Models 2.1. Prediction Procedure for Free-Surface Flow. The procedure adopted in the present study for the numerical simulation of liquid flow with a free-surface in contact with gas requires capture of the gas-liquid interface and calculation of flow fields in both the phases. A homogeneous multiphase and free-surface flow model, which is a limiting case of the general Eulerian-Eulerian multiphase flow modeling approach, has been applied to determine the shapes of the liquid surface and the flow fields in unbaffled agitated vessels. In this approach,23,26 the two fluids, air and water, in the agitated vessels share the same velocity and turbulence fields within the whole domain, which can be determined by solving a single set of governing transport equations with the volume-weighted mixture density and viscosity. This formulation relies on the fact that the air-water interface is distinct and well-defined everywhere and that there is no entrainment of one fluid into the other. Thus, in computational cells away from the air-water interface, the flow variables and fluid properties are representative of either air or water, and in the cells containing the interface, they are representative of mixtures of the phases depending on their volume fractions. The turbulent flow calculation procedure is based on numerical solutions of the three-dimensional, steady-state, Reynoldsaveraged continuity and Navier-Stokes (RANS) equations. Supplementary equations are solved to determine the interface between the gas and liquid phases and the turbulent momentum fluxes (Reynolds stresses). The governing equations for homogeneous multiphase flow of Newtonian fluids are expressed in the next section in concise forms using Cartesian tensor notation. 2.2. Governing Equations for Multiphase Flow. The Reynolds-averaged continuity and momentum equations are

∂ (F u ) ) 0 ∂xi m i and

(1)

Ind. Eng. Chem. Res., Vol. 45, No. 8, 2006 2883

([

])

∂ui ∂uj ∂ ∂ (Fmuiuj) ) µm + ∂xi ∂xi ∂xj ∂xi ∂P ∂ + (-Fmu′iu′j) + Fmgj + Fj (2) ∂xj ∂xi where ui and ui′ are the Reynolds-averaged (mean) and fluctuating velocity components, respectively, in the xi direction; P is the mean pressure; gj is the acceleration due to gravity; Fj is the body forces arising from the centrifugal and coriolis terms; and Fm and µm are the volume-weighted mixture density and viscosity, respectively. The Reynolds stresses, Fmu′iu′j, are obtained using two different closure models as described below. 2.3. Free-Surface Model. To determine the air-water interface, an interface-capturing scheme based on the volumeof-fluid (VOF) method23,26 is used. In this approach, the computation domain is discretized using a fixed mesh which spans the interface, and the volume fraction of each fluid is determined throughout the flow field by solving a continuity equation for one of the phases. This equation in the absence of any source terms can be expressed as

∂ (F φ u ) ) 0 ∂xi k k i

(3)

where φk is the volume fraction of phase k. For the present homogeneous multiphase system consisting of air and water, the above equation is solved for the volume fraction of the liquid phase (φl). The volume fraction of the gas phase (φv) can be obtained from the following:

φl + φv ) 1

(4)

In a given computational cell, the volume fraction of unity represents pure water and zero represents pure air. The airwater interface is determined by identifying the cells where the volume fraction is 0 < φl < 1. 2.4. Turbulence Modeling. As mentioned before, modeling turbulent flow in an unbaffled agitated vessel is highly challenging because of the complexity associated with strong turbulence anisotropy, streamline curvature, and rotation. This is further compounded by the likelihood of flow being laminar in regions away from the impeller even at high impeller rotational speeds corresponding to Re > 104.21 In previous numerical simulations of flows in both baffled and unbaffled agitated vessels, the turbulent momentum fluxes were generally represented by the standard high-Re k- turbulence model, which is based on the Boussinesq isotropic eddy-viscosity concept with a linear stress-strain relationship. The limitations of this type of turbulence model for the prediction of complex flows with streamline curvature, strong swirl and vortices, body forces, and in the low-Re near-wall region are well-known. There is, therefore, a requirement for more accurate representation of the turbulent transport of momentum in order to achieve an improved prediction of the flow field in agitated vessels. This goal may be realized through the use of advanced turbulence closures such as a second-moment differential Reynolds-stress turbulence (RST) model and the LES method. Although the latter approach requires large computing resources for calculations of complex three-dimensional flows, there has been a limited application to simulate flows in baffled (see, for example, refs 27,28) and in unbaffled2 vessels. The use of the LES method has not yet been prevalent for the computation of engineering flows. In the present study, calculations were

performed using an RST model, based on the modeled differential equations for the transport of the individual Reynolds stresses, and an improved eddy-viscosity based turbulence model, namely, the shear-stress transport (SST) model.25 The performances of these two turbulence closures were evaluated against the reported experimental data.14,15 2.4.a. SST Model. The SST model25,29 is a combination of the eddy-viscosity based k-ω turbulence model of Wilcox30 applied to the near-wall region and the standard high-Re k- model applied away from the wall. In previous studies,25,30,31 the k-ω model has demonstrated a superior performance in the near-wall regions versus the k- model. This led Menter25 to develop a hybrid k-/k-ω model, referred to as the baseline (BSL) model, which is based on a common formulation (see, for details, ref 25) using the k- and ω-equation with the new model constants obtained from blending the original constants of the underlying turbulence models. The values of the turbulent kinetic energy, k, and the specific dissipation rate (also known32 as turbulence frequency), ω ()/ k), can be obtained by solving their modeled transport equations:

[( ) ] [( ) ]

µt ∂k ∂ ∂ (Fuik) ) µ+ + Pk - β*Fkω ∂xi ∂xi σk ∂xi

(5)

µt ∂ω ∂ ∂ (Fuiω) ) µ+ + ∂xi ∂xi σω ∂xi ω 1 ∂k ∂ω + R Pk - βFω2 (6) 2(1 - F1)Fσω,k- ω ∂xi ∂xi k It should be noted that subscript m is omitted for simplicity. F1 is a blending function, which can be calculated using the expression given by Menter.25 Away from the solid surface, F1 ) 0 (which activates the k- model), and inside the boundary layer, F1 ) 1 (activating the k-ω model). The constants of the new model can be obtained as a linear combination of the corresponding constants of the original k-ω and the transformed k- models. The SST model25 adopts the hybrid formulation of the BSL model, and in addition, the eddy-viscosity in the boundary layer is obtained by taking into account the effect of the transport of the principal turbulent shear stress. The kinematic eddy-viscosity is obtained based on Bradshaw’s assumption that the shear stress is linearly related to the turbulent kinetic energy in the boundary layer as

Vt )

a1k max (a1ω,SF2)

(7)

where S is the invariant measure of the strain rate, F2 is a second blending function defined by Menter25 with the same limiting values as those of F1, and the constant a1 is assigned a value of 0.31. 2.4.b. RST Model. The RST model employed in this study is based on the modeled differential transport equations for the individual Reynolds stresses and is expressed in general form as

[(

)

]

∂ 2 k2 ∂Fu′iu′j ∂ (Fuku′iu′j) ) µ + csF + ∂xk ∂xk 3  ∂xk Pij + Φij - Fij + Gij (8) where the various terms from left to right represent, respectively, convection, diffusion, production, pressure-strain redistribution,

2884

Ind. Eng. Chem. Res., Vol. 45, No. 8, 2006

viscous dissipation rate, and an additional source term arising from Coriolis forces. One of the advantages of the RST model compared with the eddy-viscosity based models is its ability to simulate the additional stress anisotropy arising from Coriolis forces.33 The production term, Pij, can be expressed as

(

Pij ) - Fu′iu′k

∂uj ∂ui + Fu′ju′k ∂xk ∂xk

)

(9)

The diffusion term is modeled via a gradient-diffusion approximation. The stress dissipation process is assumed to be isotropic and is modeled in terms of the rate of dissipation of turbulent kinetic energy (). The pressure-strain redistribution term, Φij, typically consists of three parts, namely, the “slow”, “rapid”, and “wall-reflection”. A number of approaches for modeling the above terms have been proposed in the literature (see, for details, the review in ref 32). In the present study, the slow and rapid parts are represented by the nonlinear pressure-strain model of Speziale, Sarker, and Gatski,34 referred to as the SSG model, which is quadratic in the Reynolds stresses and contains no wall-reflection term. The SSG model was found34,35 to capture correct flow behavior in the near-wall equilibrium layer without recourse to the wall-reflection term. Speziale et al.34 have reported that this model performs better than the linear pressure-strain models of Launder, Reece, and Rodi36 (known as the LLR models) in predicting rotating homogeneous shear flows. The Reynolds stress transport equations are closed by computing the turbulent energy dissipation rate, , via solving a transport equation compatible with the RST model, and the turbulent kinetic energy, k, is obtained directly from the normal stresses as: ui′ui′/2. 2.5. Numerical Details. The computation domain is discretized using an unstructured mesh. The RANS equations for the mean velocity components, the volume fraction equation, and the transport equations for turbulence model quantities are discretized using an element-based finite-volume method.37 In this method, the governing differential equations are integrated over control volumes covering the whole domain. A highresolution scheme38 is used to discretize the convection terms in order to reduce the numerical diffusion errors. This scheme is a blend between the first-order upwind and a second-order upwind-biased scheme, with the blending factor varying between zero and unity. The former differencing scheme is activated if the blending factor is zero, and the latter is activated if it is unity. The blending factor is calculated locally based on the solution field of previous iteration and is close to unity. The high-resolution scheme is always bounded and, therefore, does not cause overshoots or undershoots in the solutions. It should be noted that previous computational studies6,8,15 using structured meshes demonstrated that, for a fine mesh, the calculated velocities with the third-order accurate QUICK differencing scheme39 were similar to those obtained using the upwind scheme. For the solution of the volume fraction equation in the free-surface model, a compressive differencing scheme and a surface-sharpening algorithm23,26 are used to provide greater resolution of the free-surface profile by significantly reducing numerical diffusion. The discretized RANS, pressure correction, and volume fraction equations, together with the equations for turbulence model quantities, are solved iteratively using a coupled algebraic multigrid technique.40 The Rhie-Chow algorithm41 is employed to avoid checkerboard oscillations resulting from the use of nonstaggered mesh. The calculated

Figure 1. (a) Configuration of an unbaffled agitated vessel with (b) an eight-bladed paddle impeller and (c) a six-bladed Rushton turbine impeller. Table 1. Dimensions of Agitated Vessels and Impellers vessel geometry

Rushton turbine15

paddle impeller14

vessel diameter (T) impeller diameter (D) impeller clearance (C) blade height (b) blade width (a) initial liquid height (H)

0.19 m 0.095 m (T/2) 0.063 m (T/3) 0.019 m 0.023 m 0.19 m (1T)

0.585 m 0.30 m (T/2) 0.2925 m (T/2) 0.06 m 0.14 m 0.585 m (1T)

Cartesian velocities are transformed into corresponding axial (uz), radial (ur), and tangential (uθ) components in cylindrical coordinates (r, θ, z). Finally, the mean velocities are averaged in the circumferential direction for all 360° of blade positions for the purpose of comparison with experimental data. As mentioned before, computations were performed using the ANSYS CFX-5.7 CFD code. 3. Experimental Cases Two experimental14,15 cases reported in the literature were used in this study to simulate hydrodynamics in unbaffled vessels. The first vessel used was that employed in the experimental investigation carried out by Ciofalo et al.15 and was agitated by a six-bladed Rushton turbine. The vessel geometry is shown in Figure 1. The cylindrical vessel was flatbottomed with a diameter of T ) 0.19 m and an initial liquid height of H ) 1T. The impeller diameter was D ) T/2. A clearance C ) T/3 was used between the bottom of the vessel and the midsection of the impeller disk. The detailed dimensions of the vessel are given in Table 1. The impeller rotation was in the clockwise direction, and three rotational speeds (N) were used: 139, 194, and 240 rpm. The corresponding impeller Reynolds numbers, Re ()FND2/µ), are 20 908, 29 181, and 36 100; and the tip velocities, Utip ()πDN), are 0.69, 0.97, and

Ind. Eng. Chem. Res., Vol. 45, No. 8, 2006 2885

1.19 m/s, respectively. The free-surface profiles of the liquid were measured at three rotational speeds using a vertically adjustable finger traversing along the vessel radius. The results obtained provided free-surface heights at different radial locations of the vortex. No velocity was measured for these cases. The second vessel simulated was the one used in the experimental investigation conducted by Nagata.14 The vessel was agitated by an eight-flat-bladed paddle impeller as shown in Figure 1. The cylindrical vessel was flat-bottomed with a diameter of T ) 0.585 m with an initial liquid height of 1T. The impeller diameter was T/2. A clearance of T/2 was used between the bottom of the vessel and the midsection of the blade height. The details of the vessel are given in Table 1. The impeller rotation was in the clockwise direction, and the rotational speed was 72 rpm with a corresponding Re ) 108 000 and Utip ) 1.13 m/s. Three components of mean velocity were measured using a Pitot-static tube probe. In addition, experimental data for free-surface heights at different radial locations were reported. In both experimental cases, the working fluid used was water. 4. Computational Details In the present study, the unbaffled agitated vessels are modeled using a rotating reference frame attached to the impeller without recourse to any experimental flow data as the impeller boundary conditions. In this approach, the impeller and the shaft are kept stationary while the vessel wall and the bottom are assigned an angular velocity which is equal and opposite to the impeller rotational speed. The no-slip velocity condition together with appropriate wall functions is applied to all solid surfaces. In conjunction with the SST model, an automatic wall function23,42 is employed, which blends a low-Re k-ω near-wall formulation and the standard log-law based wall-function depending on the values of y+, the dimensionless distance from the wall. For the RST model, scalable wall functions23 are used which ensure that y+ does not fall below 11.06. At the top of the computation domain in the gas phase, the zero-shear boundary condition is used. As mentioned before, computations were carried out on unstructured meshes with inflated boundaries consisting of prismatic elements to resolve the boundary layer along the solid surfaces and tetrahedral elements covering the rest of the domain. Although the initial level of water in the vessel was 1T, the height of the computation domain was extended up to 1.5T in order to capture the air-water interface. A uniform mesh was used throughout the domain except in the proximity of the solid surfaces and in the region around the air-water interface. A finer mesh was required to track the interface accurately. The vessels agitated by the Rushton turbine and paddle impeller were represented using meshes consisting of 9 × 105 and 16 × 105 cells, respectively. Figure 2 shows a typical computational mesh used in the simulations. Mesh independence tests were carried out to establish the effect of mesh size on the calculated results. For example, in the case of the vessel with a Rushton turbine, calculations were performed with three mesh sizes consisting of 6 × 105, 9 × 105, and 12 × 105 cells to examine the effect on the solutions. The velocities obtained from calculations using 12 × 105 cells were nearly identical to those obtained using 9 × 105 cells. This mesh was used in all further simulations. A similar procedure was used to determine the most suitable mesh size for the simulations of the other vessel. Achieving converged solutions of the governing equations for flow with free-surface vortex in an unbaffled agitated vessel

Figure 2. Illustration of a computational mesh.

proved to be more difficult than those for flows with a flat surface in a baffled vessel and in an unbaffled vessel with a lid in contact with the liquid. To obtain converged results, the target value of the normalized residual for each variable was set to 10-5, as generally recommended in the CFX-5.7 User Manual.23 However, it was found that satisfying the target residual was not sufficient to ensure a good convergence of the calculated flow field. The solution continued to change slowly, particularly in and around the impeller jet, with further iterations even after the residuals had reached the target value. In the present study, during the final stage of calculation, the velocity distributions at a number of locations, including the impeller region and near the top surface of the liquid, and the shape of the free-surface were monitored with the progress of iteration. The results were accepted as fully converged when the field values were almost identical over the last 500 iterations. Such a stringent convergence test ensured that the calculated results were free from numerical errors resulting from the level of convergence. Typically, 1500-2000 iterations, depending on the turbulence model used, and a lower level of residuals, of the order of 10-6, were required to obtained good convergence of free-surface flows. The simulations were run on a Sun UltraSPARC III processor with 48Gb RAM and a clock speed of 900 MHz. A typical run time was 40 h using the SST model and 73 h using the RST model for the smaller size vessel agitated by a Rushton turbine15 with a computational mesh consisting of 9 × 105 cells.

2886

Ind. Eng. Chem. Res., Vol. 45, No. 8, 2006

Figure 3. Predicted flow patterns in the vessel with a Rushton turbine15 at 240 rpm using (a) RST model and (b) SST model.

Figure 4. Predicted contours of the free-surface of water at different agitation speeds in the vessel with a Rushton turbine.15

5. Results and Discussion 5.1. Vessel Agitated with a Rushton Turbine. 5.1.a. Predicted Flow Field. The predicted mean flow patterns in terms of velocity vectors in a vertical plane midway between two impeller blades in the vessel agitated by a Rushton turbine15 are shown in Figure 3. Predictions were obtained using both the SST and RST models of turbulence for the highest impeller rotational speed of 240 rpm with corresponding Re ) 36 100. As can be seen from these figures, the overall flow structures obtained using both turbulence models are similar. The impeller produces a radial jet with a downward deflection, which expands and entrains surrounding liquid while flowing toward the side wall of the vessel. The impeller stream impinges on the wall, producing a flow structure similar to that of a wall jet. After impingement, the liquid flows along the wall vertically toward the top surface and downward toward the bottom of the vessel, and then both streams return to the impeller, thus forming expected double recirculation zones, one above and the other below the impeller stream. The predicted flow patterns obtained using both turbulence models at lower impeller speeds, 139 and 194 rpm, reveal the same general trend, except in the case of the SST model, the downward inclination of the impeller jet reduces with decreasing speeds. The absence of velocity measurements for these experimental cases precluded quantitative comparison of the performances of the SST and RST models.

The predicted flow patterns for the same experimental cases obtained by Ciofalo et al.15 using an RST model together with an iterative method for the calculation of the liquid surface profile revealed similar features, except the impeller jets were horizontal. The downward inclination of the impeller jets shown in the predictions presented in Figure 3 was also observed experimentally in an unbaffled vessel with a lid13 and in the LES simulation2 of the same flow. It is interesting to note that, as found in previous experimental and computational studies, the impeller stream in a baffled vessel agitated by a Rushton turbine has a slight upward inclination. 5.1.b. Comparison Between Predicted and Measured FreeSurface Profiles. The predicted contours of the air-water interface for the three impeller speeds, 139, 194, and 240 rpm, obtained using the VOF method together with the RST model for turbulence, are shown in Figure 4. As can be seen, at the lowest speed of 139 rpm, the free-surface of water is fairly flat in the outer region of the vortex at r > 0.06 m with a little depression in the inner region around the impeller shaft. However, the depression of the free-surface becomes larger at increased speeds, since the centrifugal forces are larger. The shape of the free-surface is controlled by the flow field in the vessel established at a given impeller speed. Figure 5 shows quantitative comparisons between the freesurface profiles obtained from the predictions of the VOF method coupled with the SST and RST models, calculated from

Ind. Eng. Chem. Res., Vol. 45, No. 8, 2006 2887

Figure 5. Comparison between the predicted and measured15 free-surface profiles (∆ experimental, s SST model, s RST model, - - - theory14).

Figure 6. Predicted flow patterns in the vessel with a paddle impeller14 at 72 rpm using (a) RST model and (b) SST model.

the theory of Nagata14 and measurements15 at the same impeller speeds. Nagata14 developed equations for free-surface profile in an unbaffled vessel based on a simplified potential flow theory of the vortex geometry (details may also be found in ref 15). In general, the CFD predicted and the theoretical free-surface profiles are in good agreement with those determined experimentally, particularly at 139 rpm. It is, however, difficult to draw a definite conclusion on the performance of the SST and RST models with regard to the quality of the surface-profile predictions at high impeller speeds. As can be seen in Figure 5, there is better agreement between the measurements and the SST predictions in the inner region of the vortex, whereas in the outer region, the RST model provides a better match with the data. The free-surface model overpredicts the liquid height, and thus underpredicts the vortex depths, in the proximity of the impeller shaft. The elevation of the liquid surface at the vessel wall is balanced by the depression of the surface in the inner region. It should be noted that a small error in the prediction results in a large difference in liquid height in this region because of the weighting of the volume on radius. This discrepancy between the CFD predictions and the experimental data increases with increasing impeller speed. A similar observation was made by Ciofalo et al.15 with regard to the overprediction of the height of the liquid surface near the impeller shaft. It appears from Figure 5 that the liquid volumes under the predicted surface profiles are larger than the volume under the measured profile. The liquid volumes calculated using

the CFD predicted free-surface profiles are the same as the initial volume of the liquid with a flat surface (i.e., liquid height ) 1T), which confirmed that the volume was conserved in the CFD calculations. However, the liquid volumes calculated based on the measured profiles are somewhat smaller than the initial volume. This could be due to the errors associated with measurements of the liquid free-surface height. 5.2. Vessel Agitated with a Paddle Impeller. 5.2.a. Predicted Flow Field. Figure 6 shows the predicted overall mean flow patterns using the SST and RST models in the vessel agitated with a paddle impeller, which was studied experimentally by Nagata.14 As shown in Figure 6a for the case of the RST model, there is a strong radial jet, compared with that predicted using the SST model, issuing from the impeller which impinges on the vessel wall and then splitting into two streams moving along the wall in the upward and downward directions. Both streams then return to the impeller, thus forming the familiar double recirculation zones, as shown in Figure 3 for the case of the Rushton turbine. Nagata14 reported streamlines determined from the measured velocities in the region above the impeller. The measured tangential velocity profiles reveal that the radius of the forced-vortex, rc (which corresponds to the radial location of the maximum tangential velocity), for this flow occurs at rc ≈ 0.2T, which is close to that determined by Ciofalo et al.15 using the correlation of Yamamoto (as reported by Nagata14), 0.19T. The predicted flow pattern in the freevortex, r g rc, above the impeller is qualitatively similar to that

2888

Ind. Eng. Chem. Res., Vol. 45, No. 8, 2006

observed experimentally. However, in the inner forced-vortex region, r e rc, the predicted mean axial and radial velocities are extremely small in the regions above and below the impeller, as shown in details below, resulting in virtual dead zones in the absence of secondary flow. This feature was also observed by Ciofalo et al.15 in their numerical simulation of Nagata’s experiment. In contrast, the experimentally determined streamlines in the same region above the impeller revealed a weak counter-rotating recirculation zone, which was consistent with the upward direction of the measured axial velocities. Figure 6b shows the predicted mean flow pattern obtained using the SST model. It can be seen that the impeller jet is less intense than that predicted by the RST model and is deflected upward. The predicted radial velocities in the impeller stream using the SST model are smaller than those obtained with the RST model. The flow structure in the outer region is less welldefined, revealing the formation of small secondary vortices near the tip of the impeller blade and in the proximity of the vessel wall. However, the SST-model predicted flow structure in the inner region of the vessel above and below the impeller is similar to that shown in Figure 6a. 5.2.b. Comparison Between Predicted and Measured Velocities. Quantitative validation of the predicted flow fields is effected by comparing the calculated radial profiles of the tangential, axial, and radial components of mean velocity with the experimental data at different axial locations from the impeller in Figures 7, 8, and 9, respectively. All the stations are above the impeller, where experimentally determined velocities were reported by Nagata.14 Predictions are obtained using both the SST and RST models. The predicted and measured tangential velocity profiles at four heights, z ) 0.39, 0.44, 0.49, and 0.54 m, from the impeller are shown in Figure 7. At all four stations, the predicted tangential velocity, in agreement with measurements, increases linearly up to a certain radial distance and then decreases slowly toward the vessel wall. The shape of the velocity profiles conforms to that of a combined vortex consisting of an inner region, r e rc, of forced-vortex motion (or solid-body rotation) and an outer region, r g rc, of free-vortex motion and is invariant in the axial direction along the height of the liquid. As can be seen in Figure 7, in the inner core of the flow within the forcedvortex region at all four heights, the magnitudes of the tangential velocity predicted by both the SST and RST turbulence models are very similar, which is in very good agreement with the measured profiles. However, in the outer free-vortex region, the RST model predicted velocities are smaller compared with those obtained from the SST model and are in better agreement with measurements. The peaks of the measured and RST-model predicted tangential velocity profiles occur at the same radial location, at rc ≈ 0.2T. In contrast, the SST model overestimates the width of the forced-vortex region and the magnitudes of the peak tangential velocity resulting in a significant overprediction of the velocity in the free-vortex region. A similar performance of this turbulence model was found in our previous simulation43 of the flow generated by a Rushton turbine in the unbaffled vessel of Ciofalo et al.15 The predicted and measured radial profiles of the mean axial velocity are shown in Figure 8 at four stations above the impeller at z ) 0.39, 0.44, 0.49, and 0.54 m. As can be seen, in the outer free-vortex region (r g 0.2T), the predictions generated using the RST model are in good qualitative and, in some areas, good quantitative agreement with the experimental data. The width of the upward-directed flow along the vessel wall is overpredicted by this turbulence model, whereas the SST model

Figure 7. Comparison between the predicted and measured14 tangential velocities at different heights from the center of the paddle impeller (∆ experimental,14 - - - SST model, s RST model).

fails to capture the measured trends in the free-vortex region. In contrast with the measurements and RST model predictions, the axial velocity profiles obtained using the SST model indicate the formation of a counter-rotating recirculation zone, as revealed more lucidly in the predicted overall flow structure shown in Figure 6b. In the forced-vortex region (r e 0.2T), the axial velocities predicted by both turbulence models decrease gradually toward the impeller shaft, and in the inner core enveloping the shaft, the velocities are almost zero. The upwardand downward-directed axial velocities measured by Nagata14 at different heights in this region are not reproduced by these turbulence models. The numerical results obtained using an RST model for this experimental case by Ciofalo et al.15 also revealed a similar trend. It is worth noting that the magnitudes of the axial velocity are much smaller than the tangential velocity over the bulk of the flow field, and accurate measurements of this quantity using a Pitot-static probe are difficult. Alcamo et al.2 reported that they were unable to obtain reliable experimental data for the mean axial velocity using PIV in an unbaffled tank covered with a lid. Figure 9 shows comparisons between the predicted and measured radial velocity profiles in the impeller stream, at z )

Ind. Eng. Chem. Res., Vol. 45, No. 8, 2006 2889

Figure 10. Predicted contour of free-surface of water in the vessel with a paddle impeller.14

Figure 8. Comparison between the predicted and measured14 axial velocities at different heights from the center of the paddle impeller (∆ experimental,14 - - - SST model, s RST model).

Figure 11. Comparison between the predicted and measured14 free-surface profiles at 72 rpm (∆ experimental, s SST model, s RST model, - - theory14).

Figure 9. Comparison between the predicted and measured14 radial velocities at different heights from the center of the paddle impeller (∆ experimental,14 - - - SST model, s RST model).

0.29 m, and above the impeller blades, at z ) 0.39 m. Because of the way in which Nagata14 reported the measured radial velocity data, it was difficult to extract velocity profiles at other heights accurately for comparison purposes. The RST-model predicted radial velocities along the impeller stream agree very well with the measurements, in comparison with the SST model predictions. The latter model significantly underpredicts the radial velocity, resulting in a weaker impeller jet, as shown in Figure 6b. In the free-vortex region above the impeller at z ) 0.39 m, again the RST model predictions are in better agreement with the measurements compared with those generated by the SST model. In the forced-vortex region, in contrast with measurements, both turbulence models return virtually zero radial velocity.

Comparisons between the predicted and measured mean velocity components shown in Figures 7-9 clearly demonstrate that a second-moment closure for turbulence provides predictions of free-surface flow in an unbaffled vessel which are in better agreement with measurements than those of an isotropic eddy-viscosity based model. The advantages of using the former model compared with the latter type stem from its ability to simulate the anisotropy of the Reynolds stresses and the contribution of Coriolis forces via the term Gij in the stress transport equations. It is worth noting that previous numerical simulations11,15,20,21 of flows in unbaffled vessels covered with a lid showed that the standard high-Re k- turbulence model tended to predict erroneously the tangential velocity distributions approaching the forced-vortex profile, with the velocity increasing monotonically from the impeller shaft to the vessel wall.15 The present numerical results, however, reveal a superior performance of the SST model in predicting the vortex profiles compared with that of the standard k- model. 5.2.c. Predicted Free-Surface Profile. The predicted airwater interface obtained using the RST model is shown in Figure 10. The water surface remains fairly flat in the outer region of the vessel with a shallow depression in the central region around the impeller shaft. The characteristics of the free-surface profile

2890

Ind. Eng. Chem. Res., Vol. 45, No. 8, 2006

produced by the paddle impeller at 72 rpm (Re ) 108 000) are similar to those produced by the Rushton turbine at 139 rpm (Re ) 20 908). Comparisons between the CFD predicted and theoretical14 free-surface profiles with the measured one are shown in Figure 11. Both the SST and RST turbulence models and the theory overpredict the free-surface height near the impeller shaft, which is balanced by underprediction of the freesurface height in the central region. However, the elevation of liquid at the vessel wall predicted by both the CFD and theory is in good agreement with the measurement. 6. Concluding Remarks Numerical simulations of turbulent flow with free liquid surface in unbaffled vessels agitated by a Rushton turbine and a paddle impeller have been carried out using a homogeneous multiphase and free-surface flow model to determine the shape of the vortex and the flow field. In this approach, a homogeneous flow model together with the VOF method for interface capturing are used to compute the entire flow field in both phases in the vessel. The turbulence is represented via the SST and a differential RST model. The flow fields in the vessel agitated with the Rushton turbine15 at three rotational speeds predicted using both the SST and RST models consist of usual double recirculation loops. There is a downward inclination of the predicted impeller jet, which was also observed experimentally13 and in the LES simulation2 in an unbaffled tank with a lid. However, the downward deflection of the impeller jet predicted by the SST model reduces with decreasing impeller speed. The RST-model predicted flow field in the vessel agitated with the paddle impeller14 reveals a strong horizontal radial jet that impinges on the wall, and a flow structure with double recirculation loops is established. This is in good qualitative agreement with the measurement in the outer free-vortex region. However, in the inner forced-vortex region, the experimentally observed counterrotating recirculation zone above the impeller is not captured in the prediction. With the SST model, the predicted flow structure is not well-defined, revealing a number of secondary vortices near the impeller blades and the vessel wall. In general, the predicted vortex profiles using the free-surface model are in good agreement with the experimentally determined profiles generated by both the Rushton turbine15 and the paddle impeller.14 There is better agreement between the experimental profiles and the SST model predictions in the inner region of the vortex, whereas in the outer region, the RST model provides better agreement. However, discrepancy exists between the predictions and measurements in the proximity of the impeller shaft at high speeds for both turbulence models. As for the predictions of mean velocity components, both turbulence models are able to capture correctly the general shape of the measured tangential velocity profiles. However, the RST model provides better predictions of the magnitudes of tangential velocity and the locations of the velocity peak, i.e., the boundary between the forced- and free-vortex regions. The SST model overpredicts the width of the forced-vortex region and the velocity in the free-vortex region. Again, the RST-model predicted axial and radial velocities are found to be in better agreement with measurements than the SST model predictions. It should be noted that the axial and radial (except in the impeller stream) velocities are significantly smaller than the tangential velocity, and hence, uncertainty exists regarding the accuracy of measurements of such small velocities using a Pitot-static tube.

Although both the SST and RST models reasonably well predict the shapes of the free liquid surface at different impeller speeds, some features of the flow structure and the mean velocity profiles are better predicted by the RST model. It appears that among the eddy-viscosity based turbulence models, the SST model performs better than the standard k- model in this type of flow. A thorough validation of predictions of the mean velocities and turbulence quantities, which is lacking in the present study because of the nonexistence of data, is necessary in order to be able to assess critically the performance of turbulence models. This is crucially dependent on the acquisition of good quality and detailed experimental data. Acknowledgment This work has been carried out as part of the Chemicals Behaving Badly (Phase II) project via a U.K. Engineering and Physical Sciences Research Council (EPSRC)/Nexia Solutions Limited (NSL) industrial case award. The authors gratefully acknowledge EPSRC and NSL for their financial support, ANSYS UK for providing the CFX-5.7 code, and the User Support Team at ANSYS, in general, and Dr. Ian Jones, in particular, for their help and useful discussions. Nomenclature C ) impeller clearance D ) impeller diameter N ) impeller rotational speed F1 ) blending function in the SST model Fj ) body forces g ) acceleration due to gravity Gij ) coriolis force h ) local free-surface height H ) initial liquid height in vessel k ) turbulent kinetic energy Re ) Reynolds number P ) mean pressure, production of turbulence energy r ) radius rc ) critical radius S ) mean strain rate T ) vessel diameter ui′ ) mean velocity in the xi direction ur ) mean radial velocity uz ) mean axial velocity uθ ) mean tangential velocity Utip ) impeller tip speed xi ) coordinate direction z ) coordinate direction Greek Letters R ) turbulence model constant β ) turbulence model constant  ) rate of dissipation of turbulent kinetic energy φ ) volume fraction µ ) viscosity νt ) kinematic eddy viscosity θ ) circumferential direction F ) density σ ) turbulence model constant ω ) specific dissipation rate (turbulence frequency) Subscripts m ) volume-weighted mixture k ) fluid phase

Ind. Eng. Chem. Res., Vol. 45, No. 8, 2006 2891

l ) liquid phase v ) gas phase Literature Cited (1) Mazzarotta, B. Comminution phenomena in stirred sugar suspensions. Am. Inst. Chem. Eng. Symp. Ser. 1993, 89 (293), 112-117. (2) Alcamo, R.; Micale, G.; Grisafi, F.; Brucato, A.; Ciofalo. M. Largeeddy simulation of turbulent flow in an unbaffled stirred tank driven by a Rushton turbine. Chem. Eng. Sci. 2005; 60 (8-9), 2303-2316. (3) Grisafi, F.; Brucato, A.; Rizzuti, L. Solid-liquid mass transfer coefficients in mixing tanks influence of side wall roughness. Inst. Chem. Eng. Symp. Ser. 1994, 136, 571-578. (4) Amanullah, A.; Tuttiett, B.; Nienow, A. W. Agitator Speed and Dissolved Oxygen Effects in Xanthan Fermentations. Biotechnol. Bioeng. 1998, 57 (2), 198-210. (5) Smit, L.; During, J. Vortex geometry in stirred vessels. In Proceedings 7th European Congress on Mixing, Bruges, Belgium, 1991. (6) Brucato, A.; Ciofalo, M.; Grisafi, F.; Micale, G. Numerical prediction of flow fields in baffled stirred vessels. A comparison of alternative modelling approached. Chem. Eng. Sci. 1998, 53 (21), 3653-3684. (7) Mavros, P.; Mann, R.; Vlaev, S. D.; Bertrand, J. Experimental visualisation and CFD simulation of flow patterns induced by a novel energy-saving dual-configuration impeller in stirred vessels. Chem. Eng. Res. Des. 2001, 79A, 857-866. (8) Montante, G. G.; Lee, K. C.; Brucato, A.; Yianneskis, M. Numerical simulations of the dependency of flow pattern on impeller clearance in stirred vessels. Chem. Eng. Sci. 2001, 56 (12), 3751-3770. (9) Armenante, P.; Chou, C. C.; Hemrajani, R. R. Comparison of experimental and numerical fluid velocity distribution profiles in an unbaffled mixing vessel provided with pitched-blade turbine. Inst. Chem. Eng. Symp. Ser. 1994, 136, 349-356. (10) Dong, L.; Johansen, S. T.; Engh, T. A. Flow induced by an impeller in an unbaffled tanksI. Experimental. Chem. Eng. Sci. 1994, 49 (4), 549560. (11) Armenante, P.; Luo, C.; Chou, C. C.; Fort, I.; Medek, J. Velocity profiles in a closed, unbaffled vessel: Comparison between experimental LDV data and numerical CFD predictions. Chem. Eng. Sci. 1997, 52 (20), 3483-3492. (12) Brunazzi, E.; Galletti, C.; Paglianti, A.; Pintus, S. Laser-Doppler measurements of turbulent-flow parameters: Comparison between baffled and unbaffled stirred tanks. Icheap-6 2003, 3, 1487-1492. (13) Vella, C.; Grisafi, F.; Micale, G.; Rizzutti, L.; Brucato, A. Near impeller flow field in an unbaffled stirred tank. In 11th European Conference on Mixing; Institute of Chemical Engineering: Bamberg, Germany, 2003; pp 629-639. (14) Nagata, S. Mixing: Principle and Applications; John Wiley & Sons: New York, 1975. (15) Ciofalo, M.; Brucato, A.; Grisafi, F.; Torraca, N. Turbulent flow in closed and free-surface unbaffled tanks stirred by radial impeller. Chem. Eng. Sci. 1996, 51 (14), 3557-3573. (16) Wu, H.; Patterson, G. K. Laser-Doppler measurements of turbulentflow parameters in a stirred mixer. Chem. Eng. Sci. 1989, 44, 2207-2221. (17) Lee, C. K.; Yianneskis, M. Turbulence properties of the impeller stream of a Rushton turbine. AIChE J. 1998, 44 (1), 13-24. (18) Hockey, R. M. Ph.D. Thesis, Imperial College, London, 1990. (19) Montante, G.; Lee, K. C.; Brucato, A.; Yanneskis, M. Experimental study of double- to single-loop transition in stirred vessels. Can. J. Chem. Eng. 1999, 77, 649-659, (20) Dong, L.; Johansen, S. T.; Engh, T. A. Flow induced by an impeller in an unbaffled tanksII. Numerical modelling. Chem. Eng. Sci. 1994, 49 (20), 3511-3518. (21) Shekhar, S. M.; Jayanti, S. CFD Study of Power and Mixing Time for Paddle Mixing in Unbaffled Vessels. Chem. Eng. Res. Des 2002, 80A, 482-498. (22) Mununga, L.; Hourigan, K.; Thompson. M. Numerical Study of the Effect of Blade Size on Pumping Effectiveness of a Paddle Impeller in an Unbaffled Mixing Vessel. In International Conference on CFD in the Minerals and Process Industries, Melbourne, Australia, 2003; pp 299304.

(23) CFX-5.7, User Manual; ANSYS: Canonsburg, PA, 2004. (24) Serra, A.; Campolo, M.; Soldati, A. Time-dependent finite-volume simulation of the turbulent flow in a free-surface CSTR. Chem. Eng. Sci. 2001, 56 (8), 2715-2720. (25) Menter, F. R. Two-equation eddy-viscosity turbulence models for engineering applications. AIAA J. 1994, 32 (8), 1598-1605. (26) Zwart, P. J.; Scheuerer, M.; Bogner, M. Numerical Modelling of Free-Surface and Cavitating Flows. Presented at ASTAR International Workshop on Advances Numerical Methods for Multidimensional Simulation of Two-Phase Flow, Germany, 2003. (27) Eggels, J. G. M. Direct and large-eddy simulation of turbulent fluid flow using the lattice-Boltzmann scheme. Int. J. Heat Fluid Flow 1996, 17 (3), 307-323. (28) Derksen, J.; Van den Akker, H. E. A. Large eddy simulation on the flow driven by a Rushton turbine. AIChE J. 2000, 46 (7), 1317-1331. (29) Menter, F. R.; Kuntz, M.; Langtry, R. Ten Years of Industrial Experience with the SST Turbulence Model. Turbul. Heat Mass Transfer 4 2003, 625-632. (30) Wilcox, D. C. Reassessment of the scale-determining equation for advanced turbulence models. AIAA J. 1988, 26, 1299-1310. (31) Sotiropoulos, F.; Patel, V. C. Turbulence anisotropy and near-wall modelling in predicting three-dimensional shear-flows. AIAA J. 1995, 33, 504-514. (32) Hanjalic, K. Advanced turbulence closure models: A view of current status and future prospects. Int. J. Heat Fluid Flow 1994, 15 (3), 178-203. (33) Wilcox, D. C. Turbulence Modelling for CFD; DCW Industries, Inc.: La Canada, CA, 1993. (34) Speziale, C. G.; Sarkar, S.; Gatski, T. B. Modelling the PressureStrain Correlation of Turbulence: An Invariant Dynamical Systems Approach. J. Fluid Mech. 1991, 277, 245-272. (35) Basara, B.; Younis, B. A. Assessment of the SSG pressure-strain model in two-dimensional turbulent separated flows. Appl. Sci. Res. 1995, 55, 39-61. (36) Launder, B. E.; Reece, G. J.; Rodi, W. Progress in the development of a Reynolds-Stress Turbulence Closure. J. Fluid Mech. 1975, 68, 537566. (37) Schneider, G. E.; Raw, M. J. A skewed, positive influence coefficient upwinding procedure for control-volume-based finite-element convection-diffusion computation. Numer. Heat Transfer 1986, 9 (1), 1-26. (38) Barth, T. J.; Jesperson, D. C. The Design and Application of Upwind Schemes on Unstructured Meshes; AIAA Paper 89-0366; American Institute of Aeronautics and Astronautics: Reston, VA, 1989. (39) Leonard, B. P. Survey of Finite Differences of opinion on Numerical Muddling of the Incomprehensible Defective Confusion Equation. Comput. Methods Appl. Mech. Eng. 1979, 19, 59-98. (40) Raw, M. Robustness of Coupled Algebraic Multigrid for the NaVier-Stokes Equations; AIAA Paper 96-0297; American Institute of Aeronautics and Astronautics: Reston, VA, 1996. (41) Rhie, C. M.; Chow, W. L. A Numerical Study of the Turbulence flow past an isolated airfoil with Trailing Edge Separation. AIAA J. 1983, 11, 1525-1532. (42) Esch, T.; Menter, F. R. Heat Transfer Predictions Based on TwoEquation Turbulence Models with Advanced Wall Treatment. Turbul. Heat Mass Transfer 4 2003, 633-639. (43) Haque, J. N.; Mahmud, T.; Roberts, K. J.; Jones, I.; Rhodes, D. CFD Simulation of Hydrodynamics in Unbaffled Agitated Vessels with FreeSurface. In 7th World Congress of Chemical Engineering, Glasgow, U.K., 2005.

ReceiVed for reView September 12, 2005 ReVised manuscript receiVed December 25, 2005 Accepted February 9, 2006 IE051021A