Prediction of Flow Pattern in Stirred Tanks: New Constitutive Equation

flow impellers, namely a pitched-blade downflow turbine (PBTD, Np = 2.1) and a hydrofoil impeller (Np = 0.34) by the impeller boundary condition a...
1 downloads 0 Views 334KB Size
Ind. Eng. Chem. Res. 2001, 40, 1755-1772

1755

Prediction of Flow Pattern in Stirred Tanks: New Constitutive Equation for Eddy Viscosity Nandkishor K. Nere, Ashwin W. Patwardhan, and Jyeshtharaj B. Joshi* Department of Chemical Technology, University of Mumbai, Matunga, Mumbai-400 019, India

The CFD simulations were carried out for the flow produced by two axial flow impellers, namely a pitched-blade downflow turbine (PBTD, Np ) 2.1) and a hydrofoil impeller (Np ) 0.34) by the impeller boundary condition approach. The tank was fully baffled, and the flow regime was turbulent. An attempt has been made to develop a new constitutive equation for Cµ in the description of eddy viscosity given by a standard k- model. The resulting predictions with the new eddy viscosity relation are compared with the experimental data. Also the comparison is sought with the predictions of impeller boundary condition approach using k- model with the standard and modified turbulence parameters (Ranade et al. Chem. Eng. Commun. 1989, 81, 225), zonal model of Sahu et al. (Ind. Eng. Chem. Res. 1998, 37, 2116) and the sliding mesh simulation using a standard k- model. For all these cases, energy balance has been established and compared with the experimental data. This paper also presents a critical review of the computational fluid dynamic (CFD) studies pertaining to the prediction of the flow produced by the axial flow impellers. Introduction Axial flow impellers are the most widely used agitators in chemical process industries. Stirred vessels provided especially with pitched-blade turbines (PBT) are often found in a variety of process industrial applications, such as solid suspension (Zweitering,1 Wiedmann et al.,2 Raghava Rao et al.,3 and Sanger and Deckwer),4 blending, dispersions (Albal et al.5), intensification of heat and mass transfer (Skelland and Seksaria,6 Hixson and Baum,7 Calderblank and MooYoung,8 and Edwards and Wilkinson),9 etc. Better understanding of the hydrodynamics prevailing in the stirred vessel leads to the reliable design of an impeller. Although the use of pitched-blade turbines is common, modeling of the flow pattern is far from complete. A close look at the published literature indicates that the modeling effort has been focused on a very small region (near the impeller), and most of the bulk region has remained neglected. Further, in the impeller region as well, good predictions have been established for some of the flow parameters (axial velocity and to some extent radial velocity) and substantial additional work is needed for the other parameters (tangential velocity, turbulent kinetic energy and its dissipation rate, eddy viscosity, etc.). These aspects have been addressed in the present work. Previous Work on CFD Simulation of Axial Flow Impellers Table 1 summarizes the previous work on the CFD simulation of the flow pattern generated by axial flow impellers. It gives the details as regards to the impeller and turbulence model used along with the grid size, the system investigated and the geometric description of the impeller-vessel configuration along with the region of * Corresponding author. Telephone: 00-91-22-4145616. Fax: 00-91-22-4145614. E-mail: [email protected].

Figure 1. Geometrical details of the stirred vessel.

comparison for predicted radial (u), axial (v), and tangential (w) velocity components, turbulent kinetic energy (k) and its dissipation rate () with the experimental data. Here z/R denotes the dimensionless axial distance from the impeller center-plane while r/R indicates the dimensionless radial coordinate. The impeller center has been considered as the origin as shown in Figure 1. The axial distance below the impeller is positive while that above the impeller is negative. Also the table gives the details regarding the quality of predictions and the conclusions drawn by the authors. It can be concluded from Table 1 that most of the simulations have been carried out for the flow in laminar region (wherein the radial flow dominates). Also, the validation of prediction of the mean and the turbulent flow characteristics throughout the vessel has been carried out only for the flow in unbaffled vessels where tangential flow predominates. Axial flow impellers have been given less consideration as compared to the radial flow impellers. Very few attempts have been made to predict the flow produced by axial flow impellers in baffled stirred vessels at high Reynolds number. Though the various modern impeller modeling approaches and the turbulence models are able to better predict the flow pattern produced by the Rushton turbine (radial flow impeller) and pitched-blade turbine at low Re, they have shown very limited success for the

10.1021/ie0004951 CCC: $20.00 © 2001 American Chemical Society Published on Web 03/10/2001

1756

Ind. Eng. Chem. Res., Vol. 40, No. 7, 2001

Table 1. Previous Work on CFD Simulation of the Flow Produced by Axial Flow Impellers in Stirred Vessel authors, impeller-turbulence model, grid size

system investigated

region of comparison

remarks

Pericleous and Patel (1987), source sink approach

straight-blade turbine, u r/R ) 0 to 1 disk turbine, z/R ) - 0.4 to 0.4 and pitched-blade v r/R ) 0 to 1 turbine. z/R ) - 0.4 to 0.4 Geometry used by Bertrand et al. (1980). H ) T ) 400 mm C/T ) 1/2, D/T ) 1/3

Ranade et al. (1989), black box approach, k- model, (max grid size) 30 × 46 × 9

pitched-blade turbine, turbulent regime, T ) 300 mm, C/T ) 1/3, D/T ) 1/3, H/T ) 1.0

u r/R ) 0.2, 0.4 They have studied the influence of the grid size, impeller z/R ) 0.7 to -0.4 boundary conditions, and the values of the model v r/R ) 0.2, 0.4 parameters on the predicted flow. They have shown z/R ) 0.7 to -0.4 that the model predictions successfully reproduce the w r/R ) 0.85 and 0.95 three-dimensionality and other essential with θ ) -45 characteristics of the flow. They have concluded that to 45 the pair of Cµ and C2 can be varied to minimize the z/R ) 0.2 to 0.33 gap between the prediction and experimental data. They k r/R ) 0 to 0.9 have carried out simulation with higher Cµ and lower z/R ) 0.533 to -0.366 C2 (parameters suggested by Abujelala and Lilley (1982)). Comparison of axial velocity was good, but tangential velocity and turbulent kinetic energy was good only near the impeller region. Contour plots for the normalized resultant velocity, k, , and the length scale were presented.

Fokema et al. (1994), black box approach, standard k- model

pitched-blade turbine, T ) 0.15 m, D/T ) 1/2 H/T ) 1, C/T ) 0.25 and 0.465

u R/R ) 0 to 1 We have carried out sensitivity analysis of the impeller z/R ) 0.177 to 0.424 boundary conditions. They have clearly shown that the k r/R ) 0 to 1, radial and axial velocities strongly affect the predicted z/R ) 0.0071 to 0.177 mean flow field. Hence, a correct specification of these  r/R ) 0 to 1 two velocities is essential. Further, the predicted z/R ) 0.0071 to 0.177 turbulent kinetic energy and the turbulent kinetic energy dissipation rate were found to be sensitive to the respective values of the boundary conditions. A significant result from their work was the insensitivity of the mean flow boundary conditions on the turbulence characteristics and vice versa. However, this result needs to be verified by a detailed comparison of all the flow characteristics throughout the tank.

Murthy et al. (1994), sliding mesh

axial turbine geometry not specified

Only pressure contour plots were presented. No comparison was sought.

Armenante et al. (1994), black box approach, k- model, ASM

PBTD, unbaffled closed vessel, T ) 0.29 m, H ) 0.32 m, D ) 0.102 m, C/T ) 0.22, N ) 5 rps

For straight-blade turbine, the predicted radial jet and the tangential velocity at the impeller center plane was found to agree well with the experimental measurements. For radial flow impellers, the radial velocity agreed within 10-20%, but the prediction of axial velocity was found to be off sometimes by as much as a factor of 2. However, for axial flow impellers, the axial velocity was found to be in good agreement, and radial velocity was found to be in poor agreement. The streamline plots were given to see the flow field qualitatively. The turbulence characteristics were not compared.

u r/T ) 0 to 0.9 z/R ) 0 to -1.32 v r/T ) 0 to 0.9 z/R ) 0 to -1.32 w r/T ) 0 to 0.9 z/R ) 0 to 0.33, -1.32

Only limited comparison was presented. Tangential velocity comparison was excellent. ASM was shown to yield superior predictions. Very poor axial velocity predictions were observed above the impeller near the center (up to 2r/T ) 0.25). No comparison for k and  was sought.

Sturesson et al. (1995), pitched-blade turbine black box (4-bladed), fully approach, turbulent regime k- model with (Re ) 104) low Reynolds flat bottom, no. model and T ) H ) 0.289 m wall functions D/T ) 1/3 near wall C/T ) 1/3

u r/R ) 0.484 z/R ) 0.506 to 0.68 v r/R ) 0.484 z/R ) 0.68 to 0.713 w r/R ) 0.484 z/R ) 0.68 to 0.713

Comparison of the local velocities with the experimental measurements was given near the vessel wall and the vessel bottom. Their studies concentrated on the near wall region and the vessel bottom. All the velocity components were severely underpredicted. Virtually no difference between the prediction with low Re models and wall functions was observed. Despite the presence of the solid boundary the mean and fluctuating velocities scale with impeller tip speed for the three impeller speeds studied. For N ) 3 rps, dimensionless values were a little less than that for 6 and 9 rps.

Sahu and Joshi (1995), several axial flow black box impellers, approach, k- turbulent regime, model, 28 × 33 geometry same as that used by Ranade and Joshi (1989)

u z/R ) -0.36 to 0.56 r/R ) 0 to 0.96 v z/R ) - 0 36 to 0.56 r/R ) 0 to 0.96 w z/R ) - 0 36 to 0.56 r/R ) 0 to 0.96 k z/R ) - 0.36 to 0.56 r/R ) 0 to 0.96

Extensive studies as regards to the effect of grid size, guess values and the parametric sensitivity analysis were carried out. Axial velocity compared well throughout while radial velocity compared well below the impeller, but poor above the impeller. The turbulent kinetic energy comparison was poor.

Ind. Eng. Chem. Res., Vol. 40, No. 7, 2001 1757 Table 1 (Continued) authors, impeller-turbulence model, grid size Harvey et al. (1995), rotating frame of reference, (max) 6 × 105

region of comparison

system investigated pitched-blade turbine, laminar regime, N ) 1.66 rps, T ) 0.145 m, D ) 0.05 m, C/T ) 0.5 (Re ) 21)

u v w

Bakker et al. (1996), black box approach k-, RNG k- and RSM turbulence model laminar: 41000 turbulent: 39000

pitched-blade turbine, T ) 0.145 m, D/T ) 0.35, C/T ) 0.46, N ) 1.66 and 8.33 rps (Re ) 21 and 21500)

u

Armenante and Chou (1996), black box approach, k- and algebraic stress model (ASM)

pitched-blade turbine, 1 or 2 impellers, turbulent regime, D ) 0.102 m, T ) 0.29 m, H/T ) 1.24, C/T ) 2/5, N ) 5 rps

u

v

v w k

remarks

r/R ) 0.1 z/R ) 0 to 0.6, 0.57 r/R ) 0.1 z/R ) 0 to 0.6, 0.57 r/R ) 0.1 z/R ) 0 to 0.6, 0.57

For low Re ()21), axial velocity was in good comparison, radial velocity predictions were good near the impeller, and the tangential velocity comparison was very poor (∼25%). Predicted velocity vector plots were not given for the upper third and bottom third. They have suggested that the predictions can be improved with the finer grids. They have observed the existence of two recirculation loops at low Re. The RFR method can be successfully used to find out the flow variables at the boundary of the impeller.

r/R ) 0 to 0.85, z/R ) 1.35 to -1.35 r/R ) 0 to 0.85 z/R ) 1.35 to -1.35

They have shown that CFM software can accurately predict the laminar flow behavior of the pitched-blade turbine if the proper velocity boundary conditions are provided around the entire periphery of the impeller. In the turbulent regime, CFM however predicted some of the features of the flow field, but dramatically underpredicts the total energy dissipation rate in the vessel. They are not currently capable of describing the large-scale instabilities in the flow. Vector plots for the mean flow by all the models were shown. Also distribution of local energy dissipation rate throughout the tank was shown which showed that the maximum energy dissipation occurs in the impeller region. For turbulent regime the quantitative comparison was presented only for axial and radial velocities. Radial velocity predictions were poor (large underpredictions) near the bottom of the tank by each of the turbulence model used (in the bottom half of the region below impeller). A slight overprediction of the radial velocity in the upper half of the region above the impeller was observed. Axial velocity though showed qualitative agreement, it was overpredicted in the region above the impeller. Axial velocity predictions were poor throughout the tank by all the models used as compared to the radial velocity. Axial velocity predictions were not in qualitative agreement near the bottom of the vessel. They recommended sliding mesh technique for the CFD simulation of flow pattern in stirred vessel. Power no. predicted was less by as much as 50%.

r/R ) 0 to 1 z/R ) 0 to 0.6, -1.0 r/R ) 0 to 1 z/R ) 0 to 0.6, -1.0 r/R ) 0 to 1 z/R ) 0 to 0.6, -1.0 r/R ) 0 to 0.83 z/R ) 0 to 0.6

They have carried out the CFD simulations using three different impeller boundary conditions as regards to the no. of surfaces on which boundary values of the axial and tangential velocities and turbulent characteristics were provided. Overall tangential velocity predictions were satisfactory. Predictions of tangential and axial velocity with k- model were poorer than with ASM Poor predictions for radial velocities. The radial velocity prediction with either model was poor especially near the vessel center. Turbulent kinetic was severely underpredicted by both the models throughout the region of the comparison. Also the difference between the predictions of two models in the upper part of the vessel was very small. Axial velocity comparison though was found satisfactory; it was qualitatively poor near the liquid surface. Contour plots for the dissipation rate were presented, and no quantitative comparison was given. Similar were the results as regards the comparison for the flow produced by the double-PBTD. They have also observed that the presence of the upper impeller alters the flow considerably, producing a strong vertical recirculation pattern between the impellers and significantly reducing the circulation flow below the lower impeller.

1758

Ind. Eng. Chem. Res., Vol. 40, No. 7, 2001

Table 1 (Continued) authors, impeller-turbulence model, grid size

system investigated

Daskopoulous and Rushton turbine, Harris (1996), PBTD, black box T ) 300 mm, approach, C/T ) 1/3, inner-outer D/T ) 1/3, technique, H/T ) 1.0 sliding mesh (same as that technique (only used by Ranade for DT) and Joshi (1989)) turbulence models: k-, RNG - (k-) (only for DT), a standard RSM black box: 40 × 20 × 46 (DT) 60 × 58 × 16 (PBTD) IO: 66 × 31 × 50 (PBTD)

region of comparison u r/R ) 0 to 0.9 z/R ) 0 to 0.533, -0.366 v r/R ) 0 to 0.33 z/R ) 0 w r/R ) 0 to 0.33 z/R ) 0 k r/R ) 0 to 0.9 z/R ) 0 to 0.533, -0.366  r/R ) 0 to 0.33 z/R ) 0

remarks They have assessed the position up to 1996. Comparison of the impeller boundary conditions with that recovered from the iteration (predicted) of the inner-outer approach between the meshes for almost all the variables was poor. The comparison was presented for maximum z/R ratio of 0.933 down to the impeller plane. No comparison for the properties above the impeller. They have observed that the available turbulence and computational models capture the qualitative features of flow and moderately successful in quantitative prediction of the mean flow but seriously underpredict turbulence properties, especially away from the impeller. stream. With IO techniques k was underpredicted by a factor of about two while  was found to be underpredicted by a factor of 30 for PBTD configuration. Similar mean axial velocity profiles by sliding mesh and black box approach, while k was severely underpredicted by sliding mesh technique but the underprediction was not as drastic as that obtained with IO. Black box approach with k- and RSM severely underpredict k in impeller stream of PBTD. Good tangential velocity comparison for PBTD vs that for DT.

Xu and McGrath (1996), momentum source model (similar to that of Pericleous and Patel), k- model, and black box, 61 000

45°, 30°, PBTD (4-blade), u r/R ) 0 to 0.55 The source sink approach improved the radial and turbulent regime, z/R ) 0 to 0.41, -0.02 axial velocity prediction and worsens the tangential T ) 0.305 m, v r/R ) 0 to 0.55, velicity as compared to the predictions by black box D/T ) 1/3, z/R ) 0 to 0.41, -0.02 approach. Power and flow no. agreed well and in C/T ) 1/3, w r/R ) 0 to 0.55 general radial velocity was poor near the bottom and N ) 3.16 rps z/R ) 0 to 0.41, -0.02 tangential velocity was poor throughout the region of comparison. Although these validations have proved that the momentum source model can predict the stirred tank flows realistically with limited computer power and requires no experimental data. Its effectiveness would be undermined when the axial flow is not a dominant component across the impeller such as DT or PBTD with small off bottom clearance. It smears out the time-dependent characteristics. Np was overpredicted (∼12%), excellent flow no. comparison. Power input was calculated by directly integrating the tangential momentum.

Ranade and Dommeti (1996), snapshot approach (source-sink), k- model, 35 × 46 × 38

pitched-blade turbine, turbulent regime, geometry same as that of Ranade and Joshi (1989)

Xuereb and Bertrand (1996), black box approach, k- model

twin propeller (MIXEL-TT), non-Newtonian (Ostwald-de Waele model) T ) 1.85 m, for lower impeller, C/T - 1/3, interimpeller spacing ) T, D/T ) 1/2

Comparison of local velocity field was not given. Countor plots of the predicted flow field and the energy dissipated presented gives only a qualitative picture of flow predicted

Harvey and Rogers (1996), MRF and moving grid: steady and unsteady state simulations, 88 980

6-bladed Rushton turbine, pitched blade turbine, system, laminar regime. Geometry same as that used byHarvey et al. (1995)Re ) 5, 33.59)

Governing equations and the solution procedure same as Harvey et al. (1995). Flow field was solved in a rotating frame of reference. Power no. matched well over a range of impeller Reynolds no. Comparison was presented for DT only. The local velocity field agreed reasonably well. Factor between the time for steady and unsteady computation ) 3 to 4), Little overprediction of tangential velocity and good axial and radial velocity comparison was obtained. No comparison for k and  was sought.

u r/R ) 0 to 1 Axial velocity compared well below the impeller, radial z/R ) 0.2 to 0.533 and tangential velocity comparison was good near v r/R ) 0 to 1 the bottom of the vessel and (very) poor in the region z/R ) -0.366 to 0.533 near the impeller and above the impeller; primary w r/R ) 0 to 1 flow no. (within 10%) and power no. (within 25%) z/R ) 0.2 to 0.533 compare well. No comparison was presented for turbulent characteristics

Ind. Eng. Chem. Res., Vol. 40, No. 7, 2001 1759 Table 1 (Continued) authors, impeller-turbulence model, grid size

system investigated

Harvey et al. (1997), multiple impeller same as that system, four 45° used in earlier pitched-blade work, 332 000 turbine with different diameters, H ) 0.5588 m, T ) 0.33147 (only laminar regime, max Re ) 86)

region of comparison u r/R2 ) 0 to 1.8 z/R2 ) 0 to 0.456, -1.175 v r/R2 ) 0 to 1.8 z/R2 ) 0 to 0.456, -1.175 w r/R2 ) 0 to 1.8 z/R2 ) 0 to 0.456, -1.175

Bakker et al. (1997), 45° pitched-blade sliding mesh turbine, (laminar technique, and transition RNG k- regime) for turbulence various Re (17 to model 12), T ) 0.3 m, D/T ) 1/3, C/T ) 1/3, N ) 3.7 rps, Re ) 10000, flat bottom, baffle off bottom clearance ) T/72 u r/R ) 0 to 0.97, z/R ) 0 to 0.185, 0 to -1.05 v r/R ) 0 to 0.97 z/R ) 0 to 0.185, 0 to -1.05 w r/R ) 0 to 0.97 z/R ) 0 to 0.185, 0 to -1.05 k r/R ) 0 to 0.97 z/R ) 0 to 0.185, 0 to -1.05

remarks Comparison was not presented for near liquid surface and the vessel bottom. It was shown that the overall model predictions were within 20% throughout the tank. Axial and tangential velocities compared well; radial velocity predictions were poor. No k and  comparison. They have studied a variety of impeller configurations, effect of relative impeller sizing, spacing, and baffling on flow distribution within the stirred tank. They have observed that global circulation strongly depends on the relative impeller sizing and spacing. Improper spacing can cause compartmentalization of flow. They have proposed the utility of the numerical method for studying complex multiple impeller flow at low Re. (R2 is the radius of largest impeller and the z ) 0 plane lies at the center of the largest impeller located at 0.1882 m above bottom) was carried out. Simulation of the turbulent flow (Re ) 10000) for only flow no. comparison. No comparison for the predicted velocities was presented. Flow no. agrees reasonably well over a range of Re. At low Re, PBTD created the radial flow and NQ significantly decreased as Re was decreased. After Re ) 400 axial flow starts for PBTD. They recommended sliding mesh technique for CFD simulation of flow pattern in stirred vessel.

Armenante et al. (1997), black box approach, k- and algebraic stress model (ASM)

6-bladed 45° pitched blade turbine, turbulent regime, closed, unbaffled with the flat bottom, H ) 0.293 m, H/T ) 1.0, D/T) 1/3, 450, N ) 11.66 rps

Np and Nq. Predictions of tangential and axial velocity with k- mode were poorer than that with ASM. Good comparison for tangential velocity was obtained except in the region near to the vessel center. The axial velocity compared well throughout the region of comparison while poor radial velocity predictions were obtained below the impeller. Turbulent kinetic energy prediction with ASM was good except in the region just below the impeller where it was severely under predicted. The standard k- model yielded poor k predictions. Np by ASM was about half of the experimental Np and was lower by 30% by standard k- model. Tangential velocities were larger than for the other two velocity components

Ljungqvist and Rasmuson (1997), black box approach, k- model, 51 000

45° pitched-blade turbine, u r/R ) 0 to 943 Have carried out open and close surface simulation. turbulent regime, z/R ) 0.1 and 0.64 Tangential velocity and turbulent kinetic energy do not T ) H ) 0.297 m, v r/R ) 0 to 0.47 compare well, very limited comparison, impeller discharge D/T ) 1/3, z/R ) 0.1 and 0.64 is unaffected by the free slip or no slip condition. C/T ) 1/3 w r/R ) 0 to 943 Overall good predictions of the three velocity z/R ) 0.1 and 0.64 components. In the bottom of the vessel all three k r/R ) 0 to 0.47 components of the velocity are under predicted. (Axial z/R ) 0.64 velocity and tangential velocity by as much as ∼40-50% and the radial velocity with deviation equal to 20%.) The axial component is considerably larger than the other two. The presence of the lid affects the tangential and radial components more than the axial, especially within the projection of the impeller.

Naude et al. (1997), black box approach, k- and RNG k- model, 24 × 23 × 55

two LUMPP LB propellers, twin propeller, turbulent regime, D ) 0.095 T ) H ) 0.190 m H ) 0.228 m D/T ) 1/2 C/T (lower impeller) ) 1/4, interimpeller clearance ) D (Re ) 45 125)

u r/R ) 0 to 0.62 Vector plots of the mean velocity were presented for the z/R ) -1.53 midbaffle plane. Very poor comparison can be seen v r/R ) 0.62 from the presented plots. Detailed velocity z/R ) 0.26 to -0.2, comparison was not given, flow no. and power no. and z/R ) agreed reasonably using the k- model. Power no. -1.06 to -1.5 predictions were significantly lower using RNG k- model. Both the impellers cannot be represented by a single boundary condition. They checked the sensitivity of the impeller boundary conditions. Power consumption was underpredicted by 45%. RNG k- model yields the good results in flows in the curved geometry, time dependent flows, low Re flows or laminar like as well as the flows with swirl. The k- model overpredicted the total energy dissipation rate and the RNG k- model is better suited for the near wall flows.

1760

Ind. Eng. Chem. Res., Vol. 40, No. 7, 2001

Table 1 (Continued) authors, impeller-turbulence model, grid size

system investigated

Weetman (1997), black box approach, multireference frame and sliding mesh approach, max. cells: sliding mesh technique: 131 000 multiple reference frame: 31 000 black box: 54 000

fluid foil impellers (A310 and A200) with only black box, turbulent regime, A310: D ) 0.43 m, T ) 1.22 m, N ) 2.5 rps A200: D ) 0.41 m. T ) 1.22 m, N ) 2.14 rps, C/T ) 1.22

Sahu et al. (1998), black box approach, k- model with optimized set of turbulence parameters, 30 × 37 × 15

various axial flow impellers, turbulent regime geometry same as that used by Ranade and Joshi (1989)

region of comparison

Only velocity vector plots were given. 7% with the average outlet velocities with multiple reference frame and sliding mesh technique. Torque predictions deviate by 20-60%. Sliding mesh technique was found to be stable only after no. of cells cross 54000. Table shows good match of the average velocities under the impeller with all three approaches.

u v w k



Brucato et al. (1998), black box approach, IO, sliding mesh technique standard k-, 582470

remarks

r/R ) 0 to 0.95 z/R ) 0 to 0.586, -0.96 r/R ) 0 to 0.95 z/R ) 0 to 0.586, -0.96 r/R ) 0 to 0.95 z/R ) 0 to 0.586, -0.96 r/R ) 0 to 0.95 z/R ) 0 to 0.586, -0.96 r/R ) 0 to 0.95 z/R ) 0 to 0.586, -0.96

6-bladed, single dual Rushton turbine and the standard constant pitch helical impellers, H/T ) 1, D/T ) 1/3, C/T ) 1/2, N ) 3.0 rps

case of axial flow impellers operating at high Reynolds number in baffled stirred vessels. The main objective of this work is to predict all the flow characteristics throughout the vessel (both in the bulk and in the near impeller region) generated by a pitched-blade downflow turbine (PBTD) and a hydrofoil impeller in a fully baffled vessel in a turbulent regime. System Description The geometrical details of the problem investigated are shown in Figure 1. Impellers investigated were the PBTD and a hydrofoil impeller. The impeller speeds were 4.5 r/s and 1.75 r/s, respectively. In both the cases, the existence of turbulent flow regime was ensured. A quadrant of the vessel was used for the numerical simulation. The impeller was placed at H/3 from the vessel bottom. In the CFD model, three-dimensional time averaged Reynolds transport equations were solved with the k-, model for turbulence while impeller action was simulated by imposing the experimental flow field as impeller boundary conditions. The equations were solved using control volume formulation on a staggered grid arrangement (Patankar)32 with the appropriate boundary conditions in a manner similar to that used in the earlier work (Sahu and Joshi16 and Sahu et al.).30

Zonal modeling significantly improved the prediction of velocities and turbulence quantities. Far away from the impeller, tangential velocity predictions are poor both quantitatively and qualitatively. The turbulent kinetic energy was underpredicted in the region away from the impeller. Energy dissipation rate comparison throughout the tank was poor.

Velocity vector plots in planes immediately upstream and downstream of the baffle were given for the flow generated by the axial flow impeller. Predicted power no. was 2.0 (close to the experimental value of 2.1). They have shown that the predictions given by inner-outer approach are comparable to those given by the SMT. (Re ) 30000 to 70000). A mixed axial-radial flow was observed with the center of the recirculation vortex shifted downward with respect to the impeller’s midplane. In the outer region, a significant axial flow was seen in the plane upstream of the baffle but is absent in the downstream plane.

Present Work In view of the success of the zonal model proposed by Sahu et al.,30 the approach has been extended in the present work. In the zonal modeling, the flow domain is divided into various subdomains depending on the physics of local flow, called zones. In each zone a separate set of model k- parameters is chosen. Though significant improvement in the prediction was achieved by employing the zonal distribution of all the three parameters in the r-z direction, the mean radial and tangential velocity lacked in good comparison. Further, the selection of optimum set of parameters (C1, C2, and Cµ ) may be a tedious job without any systematic procedure. Also, to achieve good results, one may end up with multiple combinations of the parametric distribution. These problems become more prominent if zonal optimization is to be made in all the three dimensions. In view of these difficulties and in order to establish a stepwise systematic procedure for parameter selection, it was thought to be desirable to investigate the parametric sensitivity of C1, C2, and Cµ . Here, in each case, while varying any of the parameters, other parameters were kept as that of the standard k- model. Comparison of the predicted energy dissipation rate was sought with the experimental data obtained by the eddy identification strategy of Sahu et al.33

Ind. Eng. Chem. Res., Vol. 40, No. 7, 2001 1761

Figure 2. Sensitivity of the predictions to C1. Key: ([) experimental; (1) C1 ) 1.24; (2) C1 ) 1.34; (3) C1 ) 1.44; (4) C1 ) 1.64; (A) at z/R ) 0.586; (B) at z/R ) 0.24; (C) at z/R ) -0.24; (D) at z/R ) -0.933.

Sensitivity of C1. CFD simulations with four different values of C1, namely 1.24, 1.34, 1.44, and 1.64, were carried out. Figure 2, parts A-D, shows the comparison between the predicted flow characteristics and the experimental data. Near the vessel bottom (z/R ) 0.586, Figure 1), both the radial and axial velocities were found to be insensitive to C1. The tangential velocity predictions showed sensitivity only in the region up to r/R ) 0.4. As C1 was increased, overpredictions were found to increase. Beyond r/R ) 0.4, all the simulations showed closer comparison. As regards to the turbulent kinetic energy, no systematic trend was observed. However, a value of C1 ) 1.34 gave relatively better predictions. Energy dissipation rate compared well with little overprediction beyond r/R ) 0.4 and was underpredicted for all the values of C1 employed in rest of the region. As C1 was increased,  underpredictions (r/R < 0.4) were found to marginally increase for r/R < 0.4. Figure 2B shows a comparison at a second location (at z/R ) 0.24). All the velocity components (u, v, and w) were found to be insensitive to change in C1. Turbulent kinetic energy showed a trend with respect to C1. As the value of C1 was increased, k overprediction increased; however, the changes were marginal. Comparison of  was good beyond r/R ) 0.3 and all the values of C1 employed resulted nearly in same quantitative predictions. Up to r/R ) 0.3, overpredictions could be seen for all the values of C1, however, the extent decreased with an increase in C1. At both the locations above the impeller (z/R ) -0.24 and -0.933) higher value of C1 was marginally better for u; however, it overpredicted (worsened) axial and tangential velocities (Figure 2, parts C and D). Near the impeller top (z/R ) -0.24), the turbulent kinetic energy

was found to be little sensitive to C1 for r/R < 0.6. As C1 was decreased from 1.64 to 1.34, the k comparison was found to improve while a lower value of C1 ()1.24) yielded little overprediction for r/R < 0.5. The comparison of k was closer with C1 ) 1.34. Beyond r/R ) 0.6, substantial deviation was observed for all the values of C1. The predicted k was found to be insensitive near the liquid surface (i.e., z/R ) -0.933). Near the impeller top (at z/R ) -0.24), energy dissipation rate was sensitive to C1 only in the region bounded by r/R < 0.2, while its predictions showed no sensitivity at z/R ) -0.933 throughout the remaining radial coordinate where there was substantial underprediction. Sensitivity of C2. Results of the simulations with four different values of C2, namely 1.52, 1.72, 1.92, and 2.12, and the comparison with the experimental data is shown in Figure 3, parts A-D. At z/R ) 0.586 (Figure 3A), lower values of C2 did not affect u and v predictions significantly. Higher C2 yielded larger overpredictions for r/R > 0.3 while severe underprediction was the result for r/R > 0.3. Tangential velocity was found to be sensitive for r/R < 0.3 where overpredictions increased as C2 was decreased, while it showed sensitivity only to higher value (C2 > 1.92) in rest of the radial region (r/R > 0.3). Comparison for k was poor (largely overpredicted) for r/R up to 0.4 with C2 ) 2.12, while little underpredictions were observed for r/R < 0.2 for lower C2. It was insensitive to C2 for r/R > 0.4. The energy dissipation rate was underpredicted for all the values, and relatively closer comparison was observed with C2 ) 2.12 in the region up to r/R ) 0.4. The predictions were closer and found to be insensitive to lower values of C2 for r/R > 0.4.

1762

Ind. Eng. Chem. Res., Vol. 40, No. 7, 2001

Figure 3. Sensitivity of the predictions to C2. Key: ([) experimental; (1) C2 ) 1.52; (2) C2 ) 1.60; (3) C2 ) 1.72; (4) C2 ) 1.92; (5) C2 ) 2.12; (A) at z/R ) 0.586; (B) at z/R ) 0.24; (C) at z/R ) -0.24; (D) at z/R ) -0.933.

Near the impeller bottom (Figure 3B), no sensitivity of C2 was observed for all the three velocity components (u, v, and w). However, k showed marginal sensitivity towards changes in C2 and lower C2 resulted in closer comparison. The predicted  was not affected by the change in C2 in most of the radial region, and the values were slightly higher than the experimental data. Near the impeller top (Figure 3C), lower values of C2 though improved the comparison for u for r/R > 0.4, they underpredicted u in rest of the radial region. An improvement in u with lower C2 was marginal. While, the effect was exactly reversed with higher values of C2. Very good comparison for u was observed for r/R > 0.7 with all the values of C2 employed. Axial velocity predictions were not affected by the change in C2. Tangential velocity was severely overpredicted with lower C2 for r/R up to 0.4 while, higher C2 had insignificant effect on w. The region of sensitivity was found to be up to r/R ) 0.4. The predictions of k were improved with an increase in C2 and the comparison was closer for higher C2 ()2.12). The energy dissipation rate again showed insensitivity in the region with r/R > 0.15 with severe underprediction. At z/R ) -0.933, u was severely underpredicted for all the values of C2 while the qualitative trend was not followed for C2 ) 1.72 and 1.92. Axial velocity was predicted poorly for all the values of C2 except for 1.92 up to r/R ) 0.5. For the region with r/R > 0.8, qualitatively good comparison was obtained with lower values of C2 as against the standard value. Tangential velocity was found to be very sensitive and poorly predicted (both qualitatively and quantitatively) for all the values of C2. The turbulent kinetic energy was severely overpredicted for all the C2 employed, and no

significant change in predicted k was observed as a result of change in C2. The predicted values of  were found to be insensitive to C2 and were much lower than the experimental data. Sensitivity of Cµ. CFD simulations were carried out by employing five values (0.17, 0.13, 0.09, 0.05, and 0.02) of Cµ and the predictions were compared with the experimental data (Figure 4, parts A-D). As Cµ was increased from 0.09 to 0.17, radial velocity did not show a significant change while predictions were little closer with lower Cµ ( 0.2. While no such trend was observed for r/R < 0.2, where all the Cµ values except 0.13 yielded fairly close comparison.

Ind. Eng. Chem. Res., Vol. 40, No. 7, 2001 1763

Figure 4. Sensitivity of the predictions to Cµ. Key: ([) experimental; (1) Cµ ) 0.02; (2) Cµ ) 0.05; (3) Cµ ) 0.09; (4) Cµ ) 0.13; (5) Cµ ) 0.17; (A) at z/R ) 0.586; (B) at z/R ) 0.24; (C) at z/R ) -0.24; (D) at z/R ) -0.933.

Near the impeller bottom (z/R ) 0.24), k was found to be affected significantly for Cµ ) 0.13 in the region r/R < 0.3, while the predictions given by other values of Cµ were fairly close. One can observe  (Figure 4A,B) as the most sensitive flow characteristic. The energy dissipation rate predictions near the vessel bottom (z/R ) 0.586) were sensitive to the changes in Cµ. It was underpredicted up to r/R ) 0.2 and overpredicted slightly for the rest of the radial region. The values of  were found to increase as Cµ was increased. The energy dissipation rate predictions near the impeller bottom (z/R ) 0.24) were good and insensitive for r/R > 0.6. Larger deviations (overpredictions) were found for r/R < 0.6 as Cµ was increased from 0.02 to 0.09 and relatively close comparison was sought for Cµ ) 0.02. Radial velocity was found to be sensitive in the region bounded by 0.3 < r/R < 0.7 near the impeller top (z/R ) -0.24). Axial velocity predictions were nearly insensitive for all the values of Cµ . Tangential velocity predictions were found to increase prominently for r/R < 0.4 as Cµ was reduced. Also as Cµ was decreased from 0.17 to 0.02, predicted k was found to increase for r/R > 0.4, while little change in k was observed in the region with r/R < 0.4. The energy dissipation rate was found to be very sensitive in the region with r/R < 0.15 wherein predicted  increased as Cµ was increased. Near the liquid surface (z/R ) -0.933), radial velocity comparison was significantly improved for Cµ ) 0.02 especially in the region with r/R > 0.7. In rest of the radial region, u predictions by all other values of Cµ employed were quantitatively the same. Similarly, overpredictions of axial velocity increased as Cµ was decreased for r/R < 0.7 while the deviation in v

decreased as Cµ was reduced for r/R > 0.7. Thus, one could see the close comparison for r/R < 0.7 and r/R > 0.7 with high and low values of Cµ, respectively. One can easily observe the poor comparison for tangential velocity with Cµ ) 0.02 both qualitatively and quantitatively. Only a marginal effect on w predictions was observed with other values of Cµ. No specific trend was observed. None of the simulations was able to give satisfactory comparison. Turbulent kinetic energy was significantly improved by low values of Cµ ()0.02) while all the other values yielded similar results with severe underprediction throughout the radial distance. The energy dissipation rate was found to be insensitive. Large underpredictions were observed for all the values of Cµ. Zonal selection of C1, C2, and Cµ. The constant C1 arises as the proportionality constant in the modeled term for the generation of the energy dissipation rate due to the vortex stretching by turbulence. While C2 appears as the proportionality constant in the model for the destruction of the energy dissipation rate due to tendency of viscosity to smear out the velocity fluctuations in unmodeled form of -equation. In the standard k- model, C1 is chosen so as to match the von Karman constant or the slope of the logarithmic region in the turbulent boundary layer while the value of the coefficient C2 is chosen to match the decay of isotropic turbulence generated by a grid placed in a flow stream (Ferzieger et al.).34 In view of the observed relative insensitivity of the flow predictions on C1 and C2 (Figures 2 and 3), it was thought desirable to keep the standard values of C1 ()1.44) and C2 ()1.92) as employed in the standard k- model.

1764

Ind. Eng. Chem. Res., Vol. 40, No. 7, 2001

As regards to Cµ, it defines the level of eddy viscosity. Its value can be estimated for a limiting case of twodimensional equilibrium shear layer. The conservation equation for turbulent kinetic energy is given by

( ) (

)

∂Ui ∂Uj ∂Ui ∂ νt′ ∂k′ ∂k′ ∂k′ ) + νt′ + - ′ (1) + Ui ∂t ∂xi ∂xi σk ∂xi ∂xj ∂xi ∂xj For two-dimensional system it reduces to

( ) (

)

∂k′ ∂U ∂V ∂V ∂k′ ∂ νt′ ∂k′ +U ) + - ′ (2) + νt′ ∂t ∂x ∂x σk ∂x ∂y ∂x ∂x For the case of equilibrium shear layer, the diffusive and convective transport of the turbulent kinetic energy is negligible and hence, eq 2 reduces to 2

(∂U ∂y )

νt′

) ′

(3)

which when combined with

ν t ) Cµ

k′2 ′

and Boussinesq’s hypothesis, gives



∂U k′2 ) (-u′v′) ν′t ∂y

(4)

and hence

Cµ )

( ) - u′v′ k′

2

(5)

Experimental measurements of k′ and - u′v′ in twodimensional equilibrium shear layer give a value of 0.3 for the ratio (-u′v′/k′). Hence a value of 0.09 is used as the standard value for Cµ in the k- model. Flow in stirred vessels is complex and three-dimensional in nature. There are regions of recirculating flows, wall curvature comes into the picture, and the swirl is imparted to the flow by the rotating element. In presence of these complexities, it is difficult to obtain conditions of equilibrium shear layer, and hence a constant value of Cµ is unlikely to hold. To understand the possible range of Cµ, Kumar35 made systematic measurements of the cross-correlation (-u′v′) and the turbulent kinetic energy throughout the stirred vessel. It was observed that the ratio (-u′v′/k′) was not constant, and a wide range of this ratio was seen to prevail in the stirred vessel. Similar conclusions have been drawn by several investigators. Abujelala and Lilley 36 have reported the inappropriateness of constant Cµ for the recirculating flows and the flows with curved streamlines. Shih et al.37 have prescribed the Cµ in terms of the local turbulence characteristics (k′ and ′), rotational speed of the system and the strain rate in order to satisfy the criterion of the positivity of the normal stresses and the Schwarz’s inequality in the case of highly strained flows. Sturgess and Syed38 indicated the need for variable Cµ whenever there is a transport of Reynolds stresses (nonequilibrium condition as regards to the generation and dissipation of the turbulence). Similar conditions hold in stirred vessels. Again the dependence of Cµ on the stress production of the

Figure 5. A. Contour plot of optimized Cµ in a stirred vessel. B. Contour plot of (k2/) in a stirred vessel.

turbulent kinetic energy and its dissipation rate (Cµ ) f(P/′)) can be observed from the algebraic stress model (Rodi).39 From the foregoing discussion it follows that a constant value of Cµ is not proper in a stirred vessel. Hence, in view of its realistic significance and the observed high parametric sensitivity, it was decided to carry out CFD simulations with Cµ as a parameter while keeping all other parameters standard. Initially the vessel was divided into a finite number of zones, and correspondingly, a value of Cµ was assigned to each zone. CFD simulation with the zonal distribution was carried out and the predictions were compared with the experimental data. Subsequently, the number of zones was progressively increased, and simultaneously a comparison was sought. Though the comparison was seen to be improved as compared to the other models, results were not up to the mark. Further, number of zones was increased and simulations were carried out. A stage was reached where existence of an infinite number of zones was realized. Hence, distinct values of Cµ to each grid point were assigned, and the CFD simulation was carried out. The resulting predictions were found to give a closer comparison. Hence, the distribution of Cµ was iteratively optimized to yield the predictions until good comparison was obtained.

Ind. Eng. Chem. Res., Vol. 40, No. 7, 2001 1765

Figure 6. Comparison of radial velocity with experimental data. Key: ([) experimental; (1) STD k-; (2) MOD k-; (3) ZON k-; (4) NEW k-; (5) SM STD k-.

To get rational basis for the selection of Cµ distribution, contour plots (Figure 5, parts A and B) of Cµ and (k2/) throughout the vessel were plotted. It was seen that linearly proportionate relationship existed between Cµ and (k2/). Accordingly Cµ distribution was imposed and another CFD simulation was carried out using a new set of Cµ to find out the contours of (k2/) and hence in turn the refined set of Cµ. This procedure was repeated until an excellent comparison between the experimental data and the CFD predictions was observed. To confirm the linear relationship, Cµ ) k2/ was incorporated in the eddy viscosity relation (NEW k- model) in the in-house code, and the resulting predictions were compared. The agreement with the experimental mean velocities was excellent throughout the tank except the tangential velocities near the liquid surface. Further remarkable improvement was observed in terms of agreement between the experimental turbulent kinetic energy and its rate of dissipation. In view of the proven 3D nature of flow in the stirred vessel, the flow characteristics are supposed to vary significantly in the azimuthal direction in the region

bounded by the baffles. To account for this, the transfer of tangential momentum (that is effected through the variation in eddy viscosity) should be taken into consideration. One may conclude that the improvement in the prediction of mean tangential velocity and hence the turbulent kinetic energy might have been achieved by taking into consideration the imposed 3D distribution of the turbulence parameter. New Constitutive Equation for Eddy Viscosity In light of the above discussion and the remarkable success in predicting the flow characteristics throughout the vessel we arrive at following empirical relationship between Cµ and the local k and . The k- model gives the relation for eddy viscosity as

ν′tCµ

k′2 k2 and νt ) Cµ ′ 

(6)

where νt, k, and  are dimensionless eddy viscosity, turbulent kinetic energy, and its dissipation rate,

1766

Ind. Eng. Chem. Res., Vol. 40, No. 7, 2001

Figure 7. Comparison of axial velocity with experimental data. Legends are the same as those of Figure 6.

respectively, while primed symbols indicate respective dimensional counterparts. The relation obtained in dimensionless form is

Results and Discussion

(k)2 Cµ ) 

(7)

Using the definitions of the dimensionless variables

k ) k′/Utip2,  ) ′R/Utip3, and νt ) ν′t/(RUtip) we get

Cµ )

1 (k′) RUtip ′

2

(8)

Substitution of eq 8 into eq 6 gives

ν′t )

( )

1 k′2 RUtip ′

eddy viscosity distribution on the local turbulent characteristics, impeller tip speed and the vessel diameter.

2

(9)

From the relationship obtained, we can conclude that the eddy viscosity and hence the flow field scales with the impeller tip speed and the vessel diameter in a turbulent regime. Thus, one can see the dependence of

Various turbulence models, namely standard k- (STD k-), k- model with the modified parameters (Cµ ) 0.125, C1 ) 1.44, and C2 ) 1.60) (MOD k-) employed by Ranade et al.,11 and the zonal model (ZON k-) of Sahu et al.30 were used to carry out the CFD simulation of the flow generated by the pitched-blade turbine by using impeller boundary condition approach. Also the sliding mesh (SM) simulation using a standard k- model (SM STD k-) was performed. For this purpose FLUENT commercial software was used. The results of all these simulations were compared with the experimental data of Ranade and Joshi40 at four different axial locations, namely z/R ) 0.24, 0.58, -0.24, and -0.933 at midbaffle plane (Figures 6-10). The comparison covers practically the entire region in the vessel both near the impeller and in the bulk. Radial Velocity (u). Figure 6 , parts A-D, shows the comparison for the radial velocity at four different axial locations. At z/R ) 0.58 (12 mm from the vessel bottom), radial velocity predictions as a result of all the simulations were close to each other with a relatively

Ind. Eng. Chem. Res., Vol. 40, No. 7, 2001 1767

Figure 8. Comparison of tangential velocity with experimental data. Legends are the same as those of Figure 6.

better comparison by the NEW k- model. In the region r/R ) 0-0.2, predicted u values by STD k-, MOD, and ZON k- models were qualitatively opposite to the experimental trend though the magnitude was small. Only the NEW k- model could capture the qualitative behavior of the radial velocity throughout the region. All the simulations underpredicted u by as much as 20% in the region bounded by 0.8 > r/R > 0.6. On the other hand SM STD k- yielded a severe underprediction throughout the radial coordinate. At the second location near the impeller bottom (z/R ) 0.24, Figure 6B), almost all the models underpredicted u nominally, the deviation being maximum (as much as 25%) at about r/R ) 0.35. The predictions were excellent for r/R > 0.4 by all the models except STD k- model for which it was underpredicted in the region with r/R > 0.5. Up to r/R ) 0.3, the predicted u by the NEW k- model was similar in magnitude and comparatively closer to the experimental data. It can easily be seen from Figure 6, parts A-D, that the predictions with SM STD k- model were poor (largely underpredicted up to r/R ) 0.4 and overpredicted for r/R > 0.4). At the third location (z/R ) -0.24, Figure 6C), experimental u showed a random trend. No simulation was able to replicate this behavior. Almost all the models employed resulted into good comparison near the wall and the vessel center. Little underpre-

diction was observed for r/R up to 0.3 and little overprediction for the region between r/R ) 0.3 and 0.7. For the region between r/R ) 0.3 and 0.7 the closest comparison was obtained with the MOD k- while the largest deviation observed was due to the STD k- model. The STD k- model gave good comparison up to r/R ) 0.3. The comparison was found to be closer for r/R ) 0.3 to 0.7 by the NEW k- model followed by the MOD k- model. The SM STD k- model, on the other hand, severely overpredicted u for r/R < 0.5 while it underpredicted in the rest of the radial region. Near the liquid surface (z/R ) -0.933, Figure 6D), all the models severely underpredicted u for r/R up to 0.75. Throughout the radial region, only MOD and NEW k- models better predicted qualitative behavior. Thus, NEW k- model was able to capture the radial reversal of the flow with good quantitative comparison near the vessel wall (r/R > 0.75). Comparatively closer predictions were given by MOD k- followed by NEW k- model for r/R up to 0.7. The simulation with SM STD k- model yielded nearly flat radial profile, which could be observed as both qualitatively and quantitatively poor. Axial Velocity (v). Figure 7, parts A-D, depicts the comparison of experimental axial velocity with the predictions of various simulations. At z/R ) 0.587, v was well predicted by both the STD k- and ZON k-

1768

Ind. Eng. Chem. Res., Vol. 40, No. 7, 2001

Figure 9. Comparison of turbulent kinetic energy with experimental data. Legends are the same as those of Figure 6.

models for r/R up to 0.2. Small underpredictions were given by all the models with relatively closer comparison by NEW k- model as compared to others in the region with r/R > 0.45. The best comparison was obtained with the ZON k- model for r/R < 0.45. Both the NEW and MOD k- models overpredicted in the region between r/R ) 0-0.2. One can easily observe insignificant v predictions as a result of simulation with SM STD k- model. At z/R ) 0.24, a very good comparison could be observed by all the models. There was virtually no quantitative difference in the predictions except SM STD k- model. It severely underpredicted (by as much as 400%) throughout. Near the impeller top (z/R ) -0.24) in the region with r/R > 0.4, all the models yielded good predictions with marginal underprediction by the STD k- model. The SM STD k- model, though it was able to predict qualitatively, severely underpredicted v by as much as 100% for r/R up to 0.65 and 400% over the rest of the radial distance. At z/R ) -0.933, a flat trend was closely followed by all the turbulence models. The STD k- model failed to predict upflow near the vessel wall while NEW, ZON, and MOD k- models were able to capture upward axial velocity. Axial velocity trend was found to be followed by the NEW k- model. The overprediction in the magnitude was largest by the MOD k- model. Very little overprediction could be observed as a result of all the simulations up to r/R

) 0.8. Relatively closer comparison was obtained with ZON k- while all the other models employed gave quantitatively poor comparison. SM STD k- model resulted in quantitatively very poor axial velocity comparison. Tangential Velocity (w). Comparison of the predicted tangential velocity with the experimental data is illustrated in Figure 8, parts A-D. At z/R ) 0.586, both the STD and MOD k- models severely overpredicted the tangential velocity for r/R up to 0.3. NEW and ZON k- models gave the satisfactory comparison throughout the radial region though they slightly overpredicted near the vessel center. All the models yielded good comparisons for r/R > 0.3, and there was practically no difference quantitatively. At z/R ) 0.24, a fairly good comparison throughout the radial coordinate by all the turbulence models was sought. The maximum underprediction was around 10% at r/R ) 0.3. All the models resulted in slight overpredictions in the region bounded by r/R ) 0.35-0.55. At z/R ) -0.24, predictions given by the STD k- model were good for r/R up to 0.5 while severe underprediction was the result for r/R > 0.5. Both the MOD and NEW k- models resulted in higher overpredictions near the vessel center (r/R < 0.5) while in the latter half, predictions were satisfactory as compared to the one given by STD k- model. ZON k- model though predicted trend satisfactorily, it

Ind. Eng. Chem. Res., Vol. 40, No. 7, 2001 1769

Figure 10. Comparison of turbulent energy dissipation rate with experimental data. Legends are the same as those of Figure 6.

underpredicted throughout the radial coordinate. At z/R ) -0.933 qualitatively opposite tangential velocities were obtained as a result of the simulations using STD and NEW k- models while MOD k- model gave larger values (as much as 7-8-fold). The ZON k- model gave a prediction opposite to one given by the STD k- model. Thus, in all, no model was able to present correct tangential velocity both qualitatively and quantitatively near the liquid surface. SM STD k- failed to give qualitative predictions below the impeller while above the impeller it resulted in large overpredictions. An almost flat profile of tangential velocity was the result of SM simulation indicating no circumferential flow. Turbulent Kinetic Energy (k). Figure 9, parts A-D, shows the comparison for the turbulent kinetic energy obtained as a result of various CFD simulations with the experimental data. At z/R ) 0.586, turbulent kinetic energy was severely overpredicted by STD k- model throughout the radial coordinate except near the vessel center. The deviation was as much as 100%. While the MOD k- was able to well predict k for r/R > 0.25, it severely underpredicted for r/R up to 0.25. Improvement in k predictions by ZON k- model could be seen (with maximum of 50% underpredictions). On the other hand a better comparison could be observed as a result of the NEW k- model wherein deviation is

restricted for the region between r/R ) 0.35 and r/R ) 0.55 where the maximum deviation was found to be less than 40%. At z/R ) 0.24, k was overpredicted more or less by all the models but the deviation was not high ( 0.6 by the STD k- model while a good comparison was sought throughout by all other models. At z/R ) -0.933, all the models severely underpredicted k throughout the radial distance. No model was able to capture the experimental trend. Largest deviation was observed due to MOD k- model (∼500%) followed by STD k- model (∼400%). The simulations using ZON and NEW k- models yielded good k comparison near the vessel center (r/R < 0.2). Though the k predictions given by ZON and NEW k- models were relatively closer to experimental values, underprediction as much as 65% was observed.

1770

Ind. Eng. Chem. Res., Vol. 40, No. 7, 2001

Table 2. Comparison between the Power Numbers for the PBTD and Hydrofoil Impeller impeller model

turbulence model

Np (PBTD) (exptl Np ) 2.10)

Np(hydrofoil impeller) (exptl Np ) 0.34)

IBC IBC IBC IBC SM

STD k- ZON k- MOD k- NEW k- STD k-

2.06 1.64 1.75 2.13 1.28

0.16 0.44 0.25 0.30 -

No model was able to predict k near the vessel wall satisfactorily. At all the axial locations SM STD k- model yielded very high k near the vessel center while k was almost negligibly predicted in the remaining region (r/R > 0.15). Turbulent Energy Dissipation Rate (E). Figure 10, parts A-D, illustrates the comparison for turbulent energy dissipation rate throughout the tank with the one determined by the eddy identification strategy given by Sahu et al.33 Near the vessel bottom (z/R ) 0.586), it was closely predicted by NEW k- while other models were found to yield overprediction up to r/R ) 0.4. In rest of the region, all the models resulted in closer predictions. At the second location (z/R ) 0.24), the NEW k- model was able to give a closer comparison throughout the radial coordinate while all other models resulted in overprediction up to r/R ) 0.3. At z/R ) -0.24, all the turbulence models severely underpredicted . At z/R ) -0.933, the energy dissipation rate was underpredicted by all the models. The underprediction was found to be the smallest with the NEW k- followed by the ZON k- and MOD k- models. At all the four locations considered, the SM STD k- model yielded very poor (largely overpredicted for r/R < 0.15 and severely underpredicted for r/R > 0.15) comparisons. Though one cannot rely on the experimental determination of the energy dissipation rates, it would be a good idea to see the comparison for power number, which in turn could be seen as the gross measure of the overall energy balance. Energy Balance. To establish the energy balance, the total power dissipation rate was calculated by volume integration of the predicted turbulent energy dissipation rate () as

Utip3 P)F R

∫02π ∫0H ∫0R r dr dz dθ

(10)

Experimentally, the total power consumption was measured by a torque table (the accuracy of experimental measurements was validated by measuring power consumption for a standard disk turbine for which the power number of five has been reported). From the predicted energy dissipation rate (eq 10) and the experimental power consumption (torque table), the values of power number (NP) were estimated using the following equation:

P ) NpFN3D5

(11)

Table 2 shows the comparison of power numbers obtained as a result of various simulations carried out with the experimentally measured power number for the system under consideration. The predicted NP by impeller boundary condition approach with STD and NEW k- models agreed quite well with the experimental value which implies good overall energy balance

Figure 11. Comparison of flow characteristics produced by a hydrofoil impeller. Key: ([) experimental (s) STD k-; (- -) NEW k-.

and, in turn, a good comparison for the turbulent energy dissipation rate. The ZON k- model though gave a closer comparison for all the flow characteristics; the power number was underpredicted by about 16%. However, the SM STD k- yielded a power number that was lower by 40%, indicating a poor overall prediction of the turbulent energy dissipation rate. To check the applicability of the NEW k- model, the same was employed for an axial flow impeller having a markedly different power number as compared to the PBTD. For this purpose, a hydrofoil impeller (Np ) 0.34, D ) 0.167 m, T ) 0.5 m C/T ) 1/3, N ) 1.75 rps) was used. Figure 11 shows the comparison of the predictions with the experimental data of Mishra41 at z/R ) 0.52 (near the vessel bottom). Radial and axial velocity predictions were found to be quantitatively almost the same by both the STD and NEW k- models and compared well with the experimental data. Tangential velocity was overpredicted by STD k- model while closer comparison was observed with the NEW k- model for r/R < 0.2. Their predictions were same in rest of the radial region (r/R > 0.2). The STD k- model

Ind. Eng. Chem. Res., Vol. 40, No. 7, 2001 1771

underpredicted turbulent kinetic energy for r/R < 0.3 with the closer predictions by the NEW k- model. Thus, one can observe the remarkable improvement in the tangential velocity and the turbulent kinetic energy, especially near the vessel center region. Also the power number calculated from the predicted energy dissipation rates was found to be 0.16 by the STD k- model while the NEW k- model resulted in a higher power number (0.3) which could be seen as being much closer to the experimental value of 0.34. Thus, it can be concluded that the new constitutive relation is applicable for the two axial flow impellers having widely different power numbers ranging between 0.34 and 2.1. Hence one can expect its applicability to predict the flow pattern and the turbulent characteristics in a vessel stirred by an axial flow impeller. Conclusions A new constitutive relation for Cµ and hence for the eddy viscosity has been proposed. The new model gives better predictions for the radial, axial, and tangential velocities. A more significant result is the remarkably better comparison throughout the vessel as regards to the turbulent kinetic energy and its dissipation rate. Also, power numbers calculated from the predicted turbulent energy dissipation rates agreed well with the experimental values. The new constitutive relation was also found to hold for another axial flow impeller having widely different power numbers of 0.34. Therefore, the new relation is expected to be useful for all the axial flow impellers having power numbers in the range of 0.34 and 2.1. From the relation obtained for eddy viscosity, we can conclude that eddy viscosity prevailing in a vessel stirred by an axial flow impeller scales with the local turbulent characteristics (k′ and ′), impeller speed and the vessel size. Acknowledgment The research work was supported by a grant from Department of Biotechnology, Government of India (Project No. BT/18/11/PRO485/11/97-PID). N.K.N. is grateful to the University Grants Commission for the award of fellowship. Thanks are also due to M/S Fluent Inc., USA, for providing a complementary copy of the software. Nomenclature C ) impeller clearance from the tank bottom, m Cµ, C1, C2 ) turbulence parameters in the k- model D ) impeller diameter, m H ) liquid height in the vessel, m k′ ) turbulent kinetic energy, m2/s2 k ) dimensionless turbulent kinetic energy (k′/Utip2) N ) impeller speed, rps NP ) power number of the impeller, R ) radius of the vessel, m r ) dimensional radius (m) Re ) reynolds number, ND2F/µ, T ) vessel diameter, m Utip ) impeller tip velocity, m/s U, V ) mean velocity components in x- and y-directions respectively, m/s u′, v′ ) fluctuating velocity components in x- and ydirections respectively, m/s

u, v, w ) dimensionless (with Utip) radial, axial, and tangential velocities, respectively z ) dimensional axial coordinate, m Greek Letters ′ ) turbulent energy dissipation rate, m2/s3  ) dimensionless turbulent energy dissipation rate (′R/ Utip3) F ) density of liquid, kg/m3 θ ) tangential coordinate, rad νt′ ) eddy viscosity, m2/s νt ) dimensionless eddy viscosity (νt′/RUtip)

Literature Cited (1) Zwietering, T. N. Suspending of Solid Particles in Liquid by Agitators. Chem. Eng. Sci. 1958, 8, 244. (2) Wiedemann, J. A.; Steiff, A.; Weinspach, P. M. Fluid Dynamics of Stirred Three Phase Reactors. Ger. Chem. Eng. (Engl. Transl.) 1981, 4, 125. (3) Raghava Rao, K. S.; Rewatkar, V. B.; Joshi, J. B. Critical Impeller Speed for Solid Suspension in Mechanically Agitated Contactors. AIChE J. 1988, 34, 1332. (4) Sanger, P.; Deckwer, W. D. Liquid-Solid Mass Transfer in Aerated Suspensions. Chem. Eng. J. 1981, 22, 179. (5) Albal, R. S.; Shah, Y. T.; Schumpe, A.; Carr, N. L. Mass Transfer in Multiphase Agitated Contactors. Chem. Eng. J. 1983, 27, 61. (6) Skelland, A. H. P.; Seksaria, R. Minimum Impeller Speeds for Liquid-Liquid Dispersion in Baffled Vessels. Ind. Eng. Chem. Process Des. Dev. 1978, 17, 56. (7) Hixson, A. W.; Baum, S. J. Mass Transfer Coefficients in Liquid-Solid Agitation Systems. Ind. Eng. Chem. 1941, 33, 478. (8) Caderblank, P. H.; Moo-Young, M. B. The continuous phase heat and mass-transfer properties of dispersions, Chem. Eng. Sci. 1961, 16, 39. (9) Edwards, M. F.; Wilkinson, W. L.; Heat Transfer in Agitated Vessels, Part II, Non-Newtonian Fluids. Chem. Eng. (London). 1972, 265, 328. (10) Pericleous, K. A.; Patel, M. K. The Modelling of Tangential and Axial Agitators in Chemical Reactors. Phys. Chem. Hydrodyn. 1987, 8, 105. (11) Ranade, V. V.; Joshi, J. B.; Marathe, A. G. Flow Generated by Pitched Blade Turbines II: Simulation Using k- Model. Chem. Eng. Commun. 1989, 81, 225. (12) Fokema, M. D.; Kresta, S. M.; Wood, P. E. Importance of Using the Correct Impeller Boundary and Conditions for CFD Simulations of Stirred Tanks. Can. J. Chem. Eng. 1994, 72, 177. (13) Murthy, J. Y.; Mathur, S. R.; Choudhary, D. CFD Simulation of Flows in Stirred Tank Rectors using a Sliding Mesh Technique. IChem. Symp. Ser. 1994, 136, 155. (14) Armenante, P. M.; Chou, C.; Hemrajani, R. R. Comparison of Experimental and Numerical Fluid Velocity Distribution Profiles in an Unbaffled Mixing Vessel provided with a Pitched-Blade Turbine. IChem. Symp. Ser. 1994, 136, 349. (15) Sturesson, C.; Theliander, H.; Rasmusson, A. An Experimental (LDA) and Numerical Study of Turbulent Flow Behaviour in the Near Wall and Bottoms Regions in An Axially Stirred Vessel. AIChE Symp. Ser. 1995, 91, 102. (16) Sahu, A. K.; Joshi, J. B. Simulation of Flow in Stirred Vessels with Axial Flow Impellers: Effects of Various Numerical Schemes and Turbulence Model Parameters. Ind. Eng. Chem. Res. 1995, 34, 626. (17) Harvey III A. D.; Lee C. K.; Rogers, E. S. Steady-State Modelling and Experimental Measurement of A Baffled Impeller Stirred Tank. AIChE J. 1995, 41, 2177. (18) Bakker, A.; Myers, K. J.; Ward, R. W.; Lee C. K. The Laminar and Turbulent Flow Pattern of A Pitched Blade Turbine. Chem. Eng. Res. Des. 1996, 74, 485. (19) Armenante, P. M.; Chou, C. Velocity Profiles in A Baffled Vessel with Single or Double Pitched Blade Turbines. AIChE J. 1996, 42, 42. (20) Daskopoulos, P.; Harris, C. K. Three-Dimensional CFD Simulations of Turbulent Flow in Baffled Stirred Tanks: An Assessment of the Current Position. IChem. Symp. Ser. 1996, 140, 1.

1772

Ind. Eng. Chem. Res., Vol. 40, No. 7, 2001

(21) Xu, Y.; McGrath, G. CFD Predictions in Stirred Tank Flows. Chem Eng. Res. Des. 1996, 74, 471. (22) Ranade, V. V.; Dommeti, S. M. S. Computational Snapshot of Flow Generated by Axial Impellers in Baffled Stirred Vessels. Chem Eng. Res. Des. 1996, 74, 476. (23) Xuereb, C.; Bertrand, J. 3-D Hydrodynamics in A Stirred Tank by A Double-propeller System and Filled with A Liquid Having Evolving Rheological Properties. Chem. Eng. Sci. 1996, 51, 1725. (24) Harvey, A. D., III; and Rogers, E. S.; Steady and Unsteady Computation of Impeller-Stirred Reactors. AIChE J. 1996, 42, 2701. (25) Bakker, A.; Laroche, R. D.; Wang, M. H.; Calabris, R. V. Sliding Mesh Simulation of Laminar Flow in Stirred Reactors. Chem. Eng. Res. Des. 1997, 75, 42. (26) Armenante, P. M.; Luo, C.; Chou, 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, 3483. (27) Ljungqvist, M.; Rasmuson, A. Hydrodynamics of Open and Closed Stirred Vessels. Proc. 9th Eur. Conf. Mix. 1997, 11, 97. (28) Naude, I.; Xuereb, C.; Trifu, M.; Bertrand, J. Numerical Simulation of Industrial Propellers in A Multi-Stage Agitated Vessel. Proc 9th Eur. Conf. Mix. 1997, 11, 211. (29) Weetman, R. J. Automated Sliding Mesh CFD Computations for Fluidfoil Impellers. Proc 9th Eur. Conf. Mix. 1997, 11, 195. (30) Sahu, A. K.; Kumar, P.; Joshi, J. B. Simulation of Flow in Stirred Vessels with Axial Flow Impeller: Zonal Modelling and Optimisation of Parameters. Ind. Eng. Chem. Res. 1998, 37, 2116. (31) Brucato, A.; Ciofalo, M.; Grisafi, F.; Michale, G. Numerical Prediction of Flow Fields in Baffled Stirred Vessels: A Comparison Of alternative Modeling Approaches. Chem. Eng. Sci. 1998, 53, 3653. (32) Patankar, S. V. Numerical Heat Transfer and Fluid Flow; Hemisphere: New York, 1980.

(33) Sahu, A. K.; Kumar, P.; Patwardhan, A. W.; Joshi, J. B. CFD Modelling and Mixing in Stirred Tanks. Chem. Eng. Sci. 1999, 54, 2285. (34) Ferziger, J. H.; Kline, S. J.; Avva, S. N.; Bordalo, S. N.; Tzuoo, K. L. Zonal Modelling of Turbulent Flows- Philosophy and Accomplishments: Near Wall Turbulence. In Zoran Zaric Memorial Conference; Kline, S. J., Afgan, N. H., Eds.; 1988, p 800. (35) Kumar, P. Design of Bioreactors. Ph.D. Thesis, University of Mumbai, Mumbai, India, 1997. (36) Abujelala, M. T.; Lilley, D. G. Limitations and Empirical Extensions of the k- Model as Applied to the Turbulent Confined Swirling Flow. Chem. Eng. Commun. 1984, 31, 223. (37) Shih, T. H.; Liou, W. W.; Shabbir, A.; Zhu, J. A New k- Eddy-Viscosity Model for High Reynolds Number Turbulent FlowssModel Development and Validation. Comput. Fluids 1995, 24, 227. (38) Sturgess, G. J.; Syed, S. A. Widely-Spaced Coaxial Jet, Diffusion-Flame Combustor: Isothermal Flow Calculations Using the Two Equation Turbulence Model. AIAA Paper No. 82-0113, Orlando, Florida, Jan 11-14, 1982. (39) Rodi, W. Turbulence Models and Their Application in Hydraulics; IAHR/AIRH Monograph; A. A. Balkema: Rotterdam, The Netherlands 1993. (40) Ranade, V. V.; Joshi, J. B. Flow Generated by Pitched Blade Turbines I: Measurements using Laser Doppler Anemometer. Chem. Eng. Commun. 1989, 81, 197. (41) Mishra, V. P. Fluid Mechanics of Mechanically Agitated Reactors. Ph.D. Thesis, University of Mumbai, Mumbai, India, 1993.

Received for review May 17, 2000 Accepted January 24, 2001 IE0004951