- Email: [email protected]
On the efficiency of turbulent mixing in rotating stirrers
Chemical Engineering and Processing 50 (2011) 799–809
Contents lists available at ScienceDirect
Chemical Engineering and Processing: Process Intensification journal homepage: www.elsevier.com/locate/cep
On the efficiency of turbulent mixing in rotating stirrers Filippo Trivellato ∗ Department of Civil and Environmental Engineering, University of Trento, I-38123 Trento, Italy
a r t i c l e
i n f o
Article history: Received 7 October 2010 Received in revised form 1 February 2011 Accepted 26 May 2011 Available online 23 June 2011 Keywords: Gini Kolmogorov Mixing Rushton Sterilization
a b s t r a c t Turbulent mixing is crucial in sterilizing engineering, where it is mandatory to produce the finest and the most homogeneous turbulence all over the fluid domain. This study shows how the turbulence imparted by Rushton impellers, which are typically used in several industrial productions, may be hardly improved by varying the impeller’s geometry. An impeller of new design, a perforated paddle that completely removes both limitations of Rushton impellers, is illustrated. The mathematical modeling is based on the RANS equations and the k–ε model for turbulence closure; great care was lavished to understand the accuracy of the numerical solution. The results are presented in terms of a global indicator of mixing performance that combines both the Kolmogorov scale and the Gini coefficient. The numerical solution disclosed the formation of a peripheral turbulent spot and a seemingly slow motion zone next to the shaft of the rotating perforated paddle. The imparted turbulence may be usefully fine, but many revolutions are needed to grant the safe target level of spatial homogeneity of the turbulence. This was confirmed by both a dedicated statistics treatment and the conventional statistical analysis. © 2011 Elsevier B.V. All rights reserved.
1. Introduction Mixing of single phase flows in mechanically agitated vessels is important in countless industrial processes: colloidal mills, emulsifiers, homogenizers and sterilizers for liquid food are among the industrial activities that can possibly benefit from this study. In mechanically agitated vessels at constant temperature, turbulence (or, in other words, the mixing intensity) is controlled by a variety of factors (the physical properties of the fluid under treatment; the power input, e.g. the impeller rotational speed; and the vessel-impeller design). Scores of different impellers are persistently being developed and an overwhelming literature is available that describes every detail of these impellers; this study is confined to rotating impellers which are broadly classified as radial (such as the Rushton turbine) or axial (propeller and pitched paddle turbine): both kinds of impellers revolve around the central shaft, that identifies the symmetry axis of a cylindrical vessel (Fig. 1). The impeller’s scope is to efficiently produce the most turbulent flow by injecting the least energy to the fluid to be treated: but how much turbulent is the imparted mixing? And how much spatially homogeneous is the generated turbulence? The present study aims at understanding whether revolving impellers (both radial and axial) can actually grant the safe target
∗ Tel.: +39 0461 282 631; fax: +39 0461 282 672. E-mail address: [email protected] 0255-2701/$ – see front matter © 2011 Elsevier B.V. All rights reserved. doi:10.1016/j.cep.2011.05.012
level of turbulence that is needed in sterilizing engineering together with a homogeneous intensity all over the flow field. Cavitation, buoyancy effects, chemical reactions and two-phase flows will not be considered herein. Turbulence that is approximately homogeneous and isotropic is typically produced by flows at high Reynolds numbers through screens, grids or perforated plates. Downstream of these devices, strong shear flows are created by the interactions among the jets issuing from the openings and the wakes formed behind the solid part of the devices. Beyond a distance, the flow is fully developed and the turbulence is basically isotropic. In describing homogeneous turbulence at high Reynolds numbers, it is assumed that energy is injected into eddy motions on a scale of order the size of the vessel; energy is then transferred to smaller and smaller length scales (the well known energy cascade) until it is dissipated at the smallest stable viscous scale, the Kolmogorov scale (hereinafter, the -scale). Turbulent flows with small length scales (fine grain turbulence) have the greatest rate of dissipation and their longitudinal decay is the quickest; turbulence intensities downstream of grids in 2D flows were shown [1] to be controlled by the number and size of the holes; however this fact is controversial, as claimed by the recent and comprehensive review on the dissipation features of perforated plates inserted transversally in cylindrical pipes [2]. It is pointed out that the above review is only relevant to 2D flows and its extension to rotating impellers in cylindrical domains will be investigated in this study.
800
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
of the vortical structures at which viscous dissipation occurs; the Kolmogorov scale was chosen because it conveys an immediate physical meaning and it is also helpful in dictating the smallest scale the mesh elements must have to obtain a meaningful numerical solution. 2. Comparison of impellers
Fig. 1. Slow recirculating flow patterns (invariant tori) in Rushton and disk turbines.
To understand the mixing capability of a stirrer, the temporal evolution of a tracer is checked at convenient locations within the liquid domain [3]; the mixing time is the time needed for the tracer to reach the value of the final equilibrium concentration within ±5% [6]) (or ±2% according to [4]); the average of mixing times calculated at different monitoring locations and the spread of mixing times about the averaged value give a measure of the intensity and homogeneity of the whole mixing. This procedure is quite involved and time consuming, it does not provide a detailed knowledge on the actual overall mixing time and the actual overall quality of mixing. In particular, the application of the mixing time concept to Rushton impellers cannot show a variety of features of the flow field (for example, anticipating the present results, the most turbulent spots of Figs. 8 and 9) and can only suggest improvements of impeller’s geometry by the trial and error method; this will be confirmed in Section 2. And yet the chief limitation is that the mixing time concept is intrinsically dependent on the vessel size [5] and scaling-up may turn out questionable. Khare et al. [7] considered the gas hold-up as a measure of the efficiency of gas–liquid contacting, but this will be discussed at length in Section 2. Apart from conceiving a criterion to characterize the mixing intensity, there is also the need to quantify the quality of mixing; this was done by studying the variance decay exponent or by the mix-norm (see, for instance, Ref. [8]). D’Alessandro et al. [9] quantified mixing quality in terms of entropy and Noack et al. [10] by means of Poincare maps; unfortunately both approaches are not particularly suitable for applied engineering modeling, are not of general value and their results can hardly convey handy implications to improve rotating impellers. The mixing parameter adopted in this work will be the Kolmogorov scale, Eq. (4), which is the minimum stable length scale
As it has been plainly stated by Alvarez-Hernandez et al. [11], mixers are among the most expensive and inefficient equipment in production plants; as far as radial and axial turbines are concerned, it was convincingly demonstrated that the formation of recirculating cells (also known as invariant tori, Fig. 1) is an obstacle to fluid transport. Slow recirculating flow patterns are invariably developed in the stagnate (dead) zones of the vessel that must be eliminated as they are a barrier to fine grain mixing; this translates into questioning the real utility of introducing more and more turbine paddles, which is the customary industrial inclination (this topic will be also substantiated in Fig. 9). Since these plants require many revolutions to achieve the safe target level of homogeneity, energy consumption may turn out as a matter of concern. The removal of stagnate zones is therefore mandatory not only to reduce mixing times, but also to obtain the complete homogenization all over the liquid domain, which is of paramount importance in many industrial processes: the most meaningful example in this context is the sterilization treatment in the liquid food industry. To address the above homogenization issue, the effort of the literature was ceaselessly devoted at improving the impeller’s geometry, mostly by the trial and error method; unfortunately, this approach is far from being successful, as the following elaboration proves. The performance of six impellers (four disks and two pitched paddles) were compared by Khare and Niranjan [7], who measured the gas hold-up εf , that was deemed by the authors as the overall measure of the efficiency of gas–liquid contacting. The hold-up data were analyzed by separating the contributions of large and tiny bubbles from the total gas hold-up. Fig. 2 shows two distinct pictures of the final steady state values of the total hold-up εf as a function of the power dissipation per unit volume of the initially undisturbed liquid. The experiments were intended to assess the most efficient geometry of the impeller, namely the impeller that is capable of producing the highest hold-up at the lowest input power. Fig. 2 clearly shows that the peak values of hold-ups do not differ significantly among the six tested impellers, being all the peak values comprised within the
Fig. 2. Steady state gas hold-up εf produced by six different impellers: (a) large bubbles; (b) tiny bubbles. Reprinted from Khare et al. [7] with kind permission from Elsevier.
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
801
Fig. 3. Present elaboration of Fig. 2: (a) large bubbles; (b) tiny bubbles.
narrow range 0.067–0.076 for large bubbles, and 0.13–0.17 for tiny bubbles; at the same time, the overall behaviour of the six tested impellers is unquestionably different. Apart from the fact that the gas hold-up is a basic measure of mixing efficiency, the graphs of Fig. 2 give no indication about the time required to attain the claimed steady state peak values; in other words, it is not known the amount of energy employed to reach the above steady values of hold-ups and energy consumption is unmistakably relevant to any industrial process. The following present elaboration is useful in addressing this energy issue: Fig. 3 was derived in this study from the data of Fig. 2 by taking the specific hold-up ε¯ f in the ordinate axis (the hold-up per unit input power per unit volume of treated fluid). For both large and tiny bubbles, all the impellers’ performances of Fig. 2 nicely collapse in one single curve, showing that there are vanishing differences among the tested impellers from the energy consumption point of view. Even assuming that εf might be considered as an indicator of mixing performance, apparently the energy required to obtain the target mixing intensity does not differ among the six tested impellers, even though they have different peak power values. Therefore variations of impellers’ geometry may only produce modest, if any, advancements in terms of energy savings. One further problem with impeller stirred vessels is that shear rates are strongly localized at the edge of the rotating paddles and their decay is extremely quick, being already vanishing a short distance away from the rotating paddle; this is readily seen in Fig. 4, where realistic wakes downstream of a six paddle Rushton impeller were calculated by Oshinowo et al. [12]. Many studies agreed that most of the energy (80% according to Koh et al. [13]; 60% by Wu and Patterson [14]; 70% by Cutter [15];
89% by Kresta et al. [16]; 62% by De Silva [17]) is dissipated near the impeller. To provide proper mixing to the most distant zones of the tank, a great deal of energy needs to be imparted; nonetheless, this is still far from being sufficient, especially on scale-up [5], so that stagnant zone of fluid, as a matter of fact, always exists far away from the paddle in spite of the enormous shear gradient delivered by the impeller. This is not the desirable flow field to be realized in a mixing jar for bacteria sterilization. It is also worth to mention that some industrial processes do not allow exceedingly high rotation rates due to the possible rupture of cell membranes or to possible damage of micro-organisms and sludge flocs (incidentally, high shear granulation is used in pharmaceutical, detergent and food industries where granule breakage is purposely performed [18]). The above review of the literature dealing with revolving stirrers has disclosed the following points that deserve further investigation: (i) modifications of the Rushton impellers’ geometry may yield negligible consequences on the mixing performance owing to both the creation of stagnate zones and the extreme localization of shear rates; (ii) mixing produced by grids (or perforated paddles) needs a more handy characterization; and (iii) the spatial homogeneity of the turbulence is as important as its maximum intensity. The present work aims at addressing the three above mentioned issues by the numerical simulation of a new paddle, presented in Section 3 along with a discussion on its accuracy, and at proposing the -scale as a handy indicator to qualify mixing performance (Section 4); in view of the expected strong heterogeneity, a new statistical criterion will be tested and compared to the conventional statistical analysis to address the spatial homogeneity issue (Section 4.1). 3. Numerical simulation 3.1. Batch apparatus
Fig. 4. Localization of shear rates at the outer edge of Rushton paddles rotating clockwise in a baffled tank [7]. Reproduced by kind permission of A. Bakker).
The numerical domain simulated in this study (Fig. 5) comprised a batch-operated cylindrical apparatus with perforated paddles rotating around the shaft (the symmetry axis of the cylinder) and sweeping all the internal volume. This apparatus removes both of the drawbacks discussed above (slow recirculating cells and localization of shear rates). The internal diameter of the cylindrical vessel, D = 106 mm, is representative of an actual experimental apparatus; the paddles’ diameter was di = 100 mm, leaving a 3 mm gap between the rotat-
802
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
Fig. 5. (a) The simulated batch apparatus (2 paddles, in this case); (b) the discretized domain (one quarter of a cylinder, in this case).
ing paddles and the fixed lateral wall; the same 3 mm clearance was also imposed between the paddles and both the ceiling and bottom floor; the paddle thickness was 1 mm and the diameter of the rotating shaft was 20 mm. The height of the paddle was some one quarter of the cylinder radius, namely 13 mm. The role of the holes (more specifically, their number, diameter and distribution) is controversial, as it has been thoroughly discussed by Malavasi et al. [2]. Although the number and diameters of the holes were systematically varied in the present simulations, the holes were always arranged in a staggered pattern that remained unchanged in all the simulations and is depicted in Fig. 6. The holes were modeled as cylindrical in shape with a sharp edge, since this is the customary and easiest way of construction. It is well known (see, for example, Refs. [1,2]) that the edges of orifices have a role in turbulence generation; however, in view of the unusefully difficult task of creating a complicated mesh on a rounded edge, it is firmly believed that even a less efficient edge would serve the purpose of this study: a systematic slight underestimation of turbulence parameters is expected, but this does not affect the general validity of the speculations that will be drawn from the present results. From the practical point of view, to minimize the drag of revolving plates (thus saving input power) it is mandatory to have a great number of holes, but this can pose structural problems to the plate which has to be thin. 3.2. Mesh generation and boundary conditions Taking advantage of the spatial periodicity of the cylindrical domain, only the region of fluid between any two consecutive paddles was actually discretized: this means one quarter of the cylinder or one eight and so on, depending on the adopted number of the rotating paddles. The numerical model of the stirred tank is based on an unstructured finite volume mesh partitioned in two regions: (i) a rotating mesh in the impeller region, in contact with the paddle
and moving with it; (ii) a fixed mesh elsewhere (more specifically, in the narrow gap comprised among the rotating paddles and the vessel walls) as depicted in Fig. 5. The two meshes slide past each other in an unsteady manner; the mobile mesh and the fixed mesh are separated by a mobile cylindrical interface, that provides the appropriate treatment of the continuity between the two meshes. At the mobile interface, a conservative interpolation is used for both mass and momentum equations using a set of fictitious control volumes; all the calculated fluid motion strictly arises from the paddles rotation. The penalty with the sliding mesh is the CPU time, which is an order of magnitude higher than the comparable steady state calculation. Within the rotating domain, one further interface (called volume interface) needs also to be introduced, but this is only for numerical purposes. Fig. 6 shows the mesh refinement on the paddle, on the sharp-edged holes and on the lateral gaps; the minimum and the maximum size of the cells varied from 0.01 mm through 0.5 mm. The cells’ volumes varied from 10−4 mm3 through 0.22 mm3 . The total number of cells varied depending on the simulation and it was typically of order 106 . In all the present numerical experiments, the stirred liquid was water at 20 ◦ C. The fixed lateral walls of the tank were assigned standard wall function boundary conditions; both the top and bottom walls were treated like symmetry boundaries; the surface of the moving zone was assigned an interface boundary condition; and the surface of the impeller was assigned a rotating wall boundary condition. The quality of the mesh is of paramount importance to obtain a credible numerical solution; this is particularly true for the nodes close to wall boundaries due to the steep velocity gradients; following [19], all the nodes at the wall lie within the viscous sublayer (u* y/ < 4, where u* is the friction velocity, y is the distance of the node from the wall and is the kinematic viscosity); again this is a good indication of the satisfactory resolution of the wall
Fig. 6. Mesh refinement on the paddle and on the lateral gaps.
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
803
where S is the modulus of the mean rate of strain tensor Sij ; Sk and Sε are source terms. The constants of the k–ε model that were used in the present study are the universally adopted values: c1ε = 1.44,
Fig. 7. The mobile interface poses a Courant limitation on time marching.
boundary layer and hence of the overall accuracy of the numerical simulation. The most critical boundary condition is the interface connecting the fixed mesh to the rotating mesh; the interface prompts a Courant limitation on the marching time step to negotiate the appropriate exchange of the physical information (Fig. 7): more specifically, it was found that the rotating black inner cell of Fig. 7 must be in contact with the fixed grey outer cell for at least 3–4 time steps, depending on the imposed angular velocity. 3.3. Turbulence closure The flow field of the rotating apparatus is unsteady, 3D and is governed by the Reynolds-Averaged Navier–Stokes (RANS) equations, which are far too known to be reported here; the standard k–ε turbulence model has been widely applied in the literature to close the RANS equations. Different versions of the k–ε model (low Reynolds number model; Reynolds-stress transport model) were also tested in the present work; however, vanishing differences were found among these models and this is consistent with the present literature consensus. Deglon and Meyer [20] proved that the partially satisfactory predictions due to the use of k–ε turbulence models may be due to numerical errors rather than to inherent inadequacies of the turbulent model: thus the standard k–ε closure model was herein adopted because it proved to be the best compromise between accuracy of the results and computational cost. Adopting the convention on repeated indices, the k – equation and the ε – equation for incompressible Newtonian fluid are reported below: ∂ ∂ ∂ (k) + (kui ) = ∂t ∂xi ∂xj ∂ ∂ ∂ (εui ) = (ε) + ∂t ∂xi ∂xj
+
t k
t + ε
∂k ∂xj
∂ε ∂xj
+ Gk − ε + Sk
ε2 ε + Sε + c1ε Gk − c2ε k k (1)
where is the molecular viscosity of the fluid (a physical property) and t is the turbulent viscosity (which is instead a function of the flow field, i.e. of the imparted mixing): t = c
k2
c = 0.09,
k = 1.0,
and
ε = 1.3
The RANS equations were solved using the commercial code Ansys Fluent 12.0. A first order scheme was used to discretize the spatial derivative terms in the governing equations, while a second order scheme was used for the discretization of time derivatives. The under-relaxation parameters for velocity, pressure, turbulent kinetic energy and dissipation energy were 0.3, 0.6, 0.6 and 0.6 respectively. All the computations were performed until the residuals of mean velocities, pressure and turbulence quantities were less than 10−5 . A typical numerical run, involving one million cells, took roughly one day of CPU time to achieve the converged solution in a 2.3 GHz quad-core processor.
3.4. Numerical uncertainty Although there is no shortage of CFD simulations of mixing problems, a thorough discussion on the actual numerical accuracy is hardly found in the literature; therefore, comprehensive sensitivity tests were performed in this study to check whether the present numerical solution was truly unaffected by the adopted mesh; these sensitivity tests are briefly reported in this section. The procedure followed to estimate the numerical errors generated by the present simulation is based on the valuable ASME guideline [21]. The underlying perspective here is to understand whether the present study is capable of predicting with trustworthy accuracy the outcome of mixing numerical simulations in cases in which the measure of turbulent data is not an easy task, as in the present case. Uncertainty in inlet boundary conditions is well known to be a potentially significant contributor to the overall uncertainty; the degree of sensitivity of the present numerical solution to perturbations in inlet conditions was checked by varying the inlet turbulence level from 1% to 30%. The inlet boundary condition was imposed on the holes of the paddle and on the lateral gaps; the flow field is strongly accelerating downstream of the inlet boundary, and then it diverges, slowing down and interacting with the nearby jets. The resulting shear flow was observed to be not affected (less than 1%) by the above variation of the inlet turbulence level; this is reasonable, since accelerating zones are known to dissipate turbulence and this fact considerably mitigates the uncertainty due to the inlet turbulence. The method herein followed for the evaluation of discretization errors is the Grid Convergence Index (GCI). The first step of the procedure is the definition of the representative grid size h:
h=
1 Vi N N
1/3 (2)
i=1
ε
The term incorporating the Boussinesq hypothesis, Gk , describes the generation of turbulent kinetic energy due to the base flow gradients and is given by: Gk = t S 2 ,
c2ε = 1.92,
S 2 = 2Sij Sij =
∂2 Ui Uj ∂Ui ∂Ui + ∂xj ∂xj ∂xi xj
where Vi is the ith cell volume and N is the total number of cells used for the computations. Three significantly different grids were selected and values of the Kolmogorov scale (which is the key variable here) were calculated. In this study, the grid refinement factor r = hcoarse /hfine was always greater than 1.3 (a value based on experience [21]). Roughly geometrically similar cells were always used.
804
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
The apparent order p of the method was calculated by using the following expression: p=
|ln|ε32 /ε21 | + q(p)| ln(r21 )
q(p) = ln
p
r21 − s
(3)
p
r32 − s
s = 1 · sgn
ε 32
ε21
where h1 < h2 < h3 , r21 = h2 /h1 , r32 = h3 /h2 , ε32 = 3 − 2 , ε21 = 2 − 1 , k stands for the solution derived by means of the kth grid (k = 1, 2, 3); sgn is the function signum. The absolute value in Eq. (3) is necessary to ensure extrapolation toward h = 0 [22]. Negative values of the ratio ε32 /ε21 were found in the present simulations; this is an indication of oscillatory convergence, a typical problem occurring when dealing with geometrically complicated fluid domains; however, it is noted that these oscillations were not significant, being their maximum values less than 3%. A fair agreement was observed between the apparent order of the scheme with the formal order used and this is an indication of the grids being in the asymptotic range. 21 were then calculated along with The extrapolated values ext the approximate relative error ea21 : p
21 =
ext
r21 1 − 2 p
r21 − 1
1 − 2 ea21 =
1
The fine-grid convergence index GCI21 fine is given by: GCI21 fine =
1.25ea21 p
r21 − 1
The local order of accuracy p calculated from Eq. (3) was between 0.2 and 3.2, with an average value of pave = 1.5 ; oscillatory convergence was found in 10% of the nodes. The maximum discretization error was some 30%, which is admittedly a large error; however, it must be emphasized that this value was produced close to the boundary layer and translates in a maximum uncertainty of the local velocity of about 0.004 m/s; from a practical point of view, this error cannot possibly be of any relevance as far as the overall accuracy is concerned. Following [23], the iteration error was estimated by
niter,i ∼
in+1 − in i − 1
where n is the iteration number and i is the principal eigenvalue of the solution matrix of the linear system, which can be approximated as: i ∼
|| in+1 − in || || in − in−1 ||
where || · || is an appropriate norm (e.g., the infinity norm L∞ ). The uncertainty ıiter in iteration convergence can then be estimated as: ıiter ∼
described method was found to be typically almost one order of magnitude smaller than the discretization error; so there are reasons to trust the accuracy of the present numerical simulation.
|| niter,i || ave − 1
where ave is the average value of i over a reasonable number of iterations. The iteration error calculated according to the above
4. Results The following parameters were herein tested: (i) number of paddles: 2, 4 and 8 (i.e., 4, 8 and 16 half-paddles); (ii) number of holes in each half-paddle: 0 (i.e., solid plate), 6, 12, 20 and 30; and (iii) diameters of the sharp-edged holes: 0.5, 1, 1.5 and 2 mm. The angular velocity of the shaft varied: Appendix A reports a selection of numerical runs performed at 50 and 100 rad/s. The maximum number of possible combinations of the above set of parameters is huge; however, during the study, it was clear that many combinations were not worth to be tested because of vanishing differences among them. Only 220 selected runs were actually performed; and yet, since most of them turned out to produce nothing but a slight variation in the results, for brevity’s reasons only meaningful results will be actually illustrated hereinafter. The interested reader may find on Appendix A all of the details useful to reproduce the present numerical experiments. The purpose of this section is to quantify turbulent mixing rather than giving sharp indications to pursue the optimal construction of the vessel. The results are presented in terms of the Kolmogorov scale, the -scale:
=
3 ε
1/4 (4)
where is the kinematic viscosity and ε is the turbulent dissipation energy, whose evolution is governed by Eq. (1). The -scale is herein proposed as an easy, and yet effective, parameter to deal with, that dictates indications on the size the mesh elements must have to obtain a meaningful numerical solution. Fig. 8 gives the picture of the -scale distribution in the midplane intersecting the central row of the holes (i.e., at half height of the paddle; the flow is counterclockwise, from below to left). The column of Fig. 8a shows the effect that the variation of the hole diameter produces on the -scale at an angular velocity of 100 rad/s and with 2 paddles (4 half-paddles); the column of Fig. 8b shows the effect of the variation of the holes’ number at 100 rad/s and 2 paddles. There are no stagnate zones. Recalling that di is the impeller diameter and N is the number of revolutions per second, the global Reynolds number: Re =
Ndi2
(5)
is representative of the turbulent regime, being greater than 4 × 104 in all of the present numerical runs. As a first observation, it seems that no dramatic diversification of the results is apparent. The second impression is that the most beneficial and serviceable zone, i.e. the brightest spot with the lowest -scale, is always localized by centrifugal forces at the top peripheral region; the area of this zone is limited in extent, and never exceeding 15% of the entire domain. The turbulence produced in the narrow gap existing among the rotating paddles and the lateral walls of the tank is fine grain turbulence similar to the turbulence produced behind the holes of the rotating paddles. The central region, close to the shaft, behaves basically like a slow motion zone, and this is consistent with what already observed by Bartok et al. [24] in a similar rotating apparatus; it was stated there, and confirmed here, that the mixing pattern in the
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
805
Fig. 8. Distribution of the -scale as a function of the holes’ diameter (column a) and of the holes’ number (column b). Counter clockwise flow from below to left; darkest spots: 1.4 × 10−4 m; brigthest spots: 10−6 m.
central core is modest and almost independent from the angular velocity. Perforated plates may operate either as turbulence suppressors or as turbulence generators. The turbulence generated downstream of a perforated plate is the combined result of three distinct phenomena: the passage of the upstream turbulence through the paddles’ openings; the shear turbulence produced downstream by the interactions of the jets exiting the screen; and the turbulence produced by the perforated plate itself. In the generator role, the plate has a coarser pattern of holes and is operated at high Reynolds numbers. It is known that the turbulence scale imposed initially at the screen controls the turbulence intensity (and below some minimum initial scale, the turbulence always decays); for example, for sufficiently fine plates, the turbulence is of small scale and decays
after 500 hole diameters. At roughly 40 hole diameters, the turbulence is nearly homogeneous and isotropic. The above decay law downstream of grids is only known for 2D flows; it is emphasized that no results are available for axialsymmetric rotating domains. Fig. 9 shows the effects of the variation of the paddles number and of the angular velocity; increasing the number of paddles does not produce finer turbulence, a result possibly unexpected a priori since increasing the number of blades is a popular method used in industries. Interestingly enough, Fig. 9 shows that a 50 rad/s angular velocity produces a noticeable zone of finer turbulence and more extended than in the 100 rad/s run; even the gap between the rotating paddle and the fixed wall of the vessel also plays a better role in producing fine turbulence. That set-up combination turned out to be the most efficient in terms of fine turbulence
806
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
Fig. 9. -Scale distribution as a function of the paddles’ number and the angular velocity. Counter clockwise flow from below to left; darkest spots: 1.4 × 10−4 m; brightest spots: 7 × 10−6 m.
within the tested range of experiments (but nothing can be said at this stage about its homogeneity); it is recalled that providing the optimal setting configuration is beyond the reach of the study. The analysis of all the results (including those not presented explicitly here) showed that the revolving velocity of the paddles (in other words, the power input) is an important parameter in generating the smallest -scale. This result is consistent with what was obtained by Sanchez Perez et al. [25], who enforced the idea that the shear rate in turbulent regime is a function of the angular velocity. Schubauer et al. [26] found that a single screen causing a strong pressure drop is a more satisfactory turbulent promoter than using several screens in series. This evidence is clearly seen in Fig. 9 when comparing the 2- and 8-paddle results: close to the shaft, there is a central zone of slow motion flow; this zone does not produce fine turbulence and is not valuable for sterilizing mixing. Thus the present apparatus, even though it completely removes the limitations described in the introduction (slow recirculating regions and localization of shear rates), still has a questionable performance due to the formation of both a peripheral turbulent spot and a central core of plug zone. And it seems clear that even the optimal combination of the parameters (which was not pursued here) would produce nothing but a slight improvement of the mixing intensity. This translates into the necessity to make many
revolutions to grant the required safe level: this issue is particularly mandatory in liquid food sterilization. In addition, the ideal mixing should produce in the flow field as much a uniform -scale as possible and a criterion is needed to qualify the spatial variability of the -scale. This is of particular relevance in the tested device which, as proved in Figs. 8 and 9, produces a marked spatial heterogeneity of the -scale. 4.1. Statistical treatment and Gini coefficient The numerical simulations described in Section 4 provided the -scale values in each cell of the domain; the range comprised between the minimum through the maximum value of was divided in 200 classes; Fig. 10 shows the frequency distribution of the local -scale values, where after the peak value, the function is typically monotonically decreasing. But the monotonic decrease was not observed in all of the present numerical simulations owing to the strong spatial heterogeneity of the -scale. Secondary peak values did occur that were located in between the first peak and the tail end of the distribution; these secondary peaks were at times as important as the first one. This fact has been tackled here by introducing the Lorenz’s function [27], which was originally intended to describe inequalities of the wealth distribution and is often used in biodiversity studies; it seems that the Lorenz function has not been previously applied in mixing engineering.
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
807
Fig. 10. Example of the frequency distribution of the local -scale.
Firstly, the data sample of the -scale was sorted in ascending order; then the cumulative frequency distribution was calculated by assigning the Weibull frequency to the ith value: Fi =
can be hardly smoothed away by any combination of the setting parameters.
i M+1
where M is the sample size. The cumulated amount of the ith value is then calculated as:
i
i =
4.2. The overall performance For conveniency’s sake, it may be useful to introduce the parameter as a new indicator of the overall performance:
j=1 j M j=1 j
= G¯
Fig. 11 shows an example of the i − Fi plotting: if all the values of the data sample were equal (in other words, if the -scale were perfectly distributed) the graphical representation would be the line of perfect equality, the line OA; by contrast, if only one value were different from zero, the line OBA would be the line of perfect inequality, meaning that the total amount of the -scale is wholly represented by just a single data of the sample. If data values are positive (as in the present -scale sample), the actual Lorenz curve must lie somewhere in between the two extreme lines of perfect equality and perfect inequality: the closer the Lorenz curve to the line OA, the more homogeneous the -scale. The Gini’s coefficient G [28] was conceived to quantify the inequality of data samples and it is defined as the ratio between the shaded area of Fig. 11 and the area of triangle OAB: G=
Fig. 12. The overall performance (number of holes: 6, 20 and 30; dh = 0.5, 1 and 2 mm).
shaded area area triangle OAB
The Gini coefficient is comprised between 0 and 1: if the -scale values were perfectly distributed, the Lorenz function would collapse on the line of perfect equality, the bisectrix OA, and G takes its minimum value (Gmin = 0); by contrast, the most unequal distribution is given by the line OBA and the Gini coefficient takes its maximum value (GMax = 1). The G coefficient calculated by the present simulations was comprised between 0.74 and 0.97, as reported in Appendix A, showing a great deal of spatial heterogeneity that
Fig. 11. Lorenz function and Gini coefficient.
(6)
where ¯ is the weighted mean of (where the weights are the volume of the cells). The parameter is intended to briefly summarize, in only one number, the global performance of the stirring device: the minimum value of is the desirable value for mixing engineering, since it combines at the same time both the finest mixing and the most homogeneous distribution of it. The minimum values of , in other words, detect the best combinations of the parameters. Amazingly enough, the range of is somewhat narrow (see Appendix A), indicating there is no dramatic difference in the results whichever the combination of the parameters. Fig. 12 displays the overall performance of the imparted mixing as a function of the Reynolds number, defined by Eq. (5); for clarity’s sake, only the maximum (2-paddle) and the minimum (4-paddle) values are presented. Fig. 12 shows the asymptotic behaviour starting at Re ∼ 105 (depending on the paddle number), meaning that the holes’ number and diameter are not parameters as relevant as it could be initially thought; by contrast, the number of paddles and the angular velocity are important parameters in this apparatus. It is emphasized that increasing the number of paddles does not produce finer turbulence, but it does promote a better homogenization; the combination of both effects are included in the parameter , that may therefore turn out lower for a 4-paddle impeller. Rather obviously, the simulations with solid paddles (i.e., with no openings) turned out to perform inefficiently in terms of mixing (details can be found in Appendix A). The analysis of the -scale values was also performed by means of the conventional statistics parameters: the weighted mean ¯ (where the weights are the volume of the cells); the standard deviation based on the weighted mean; the skewness coefficient, Sk ; and the kurtosis coefficient Ku . Again their definition is far too known to be reported here. Summing up briefly the main results, the statistical analysis based on the conventional parameters confirmed the strong heterogeneity of the flow field and it is worth to report here only the major differences from the above described (Section 4.1) dedicated statistics. More specifically:
808
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
(i) both the standard deviation and the parameter revealed the same behaviour on the involved parameters (paddles numbers, angular velocity and diameter of holes); (ii) most of the above simulations showed a negative skewness coefficient, indicating a significant concentration of -scale values lower than the average value (Appendix A); (iii) the 2-paddle turbine has the lowest kurtosis coefficient, meaning a rounded peak and short, thin tail of the distribution. Even though the present results are relevant to closed cylindrical containers with rotating perforated paddles, however the dependence of the pressure drop (not shown here for brevity’s sake) on relevant parameters – such as the global Reynolds number, the paddle thickness and the holes number and diameter – were found to be in qualitatively agreement with what described by Malavasi et al. [2], who tackled the different case of the perforated plate transversally inserted in a cylindrical conduit. Although the present simulations have been performed in a batch mixing apparatus, it is clear that the obtained results are just as valid, at least qualitatively, as in a semi-continuous or even a continuous mixing apparatus [29]. Also, all the literature studies which have been carried out on liquid and gas mass transfer in mechanically agitated reactors indicated that the operating ambient pressure does not influence the mass transfer coefficient; thus the present results can be extended even to micro-organisms inactivation based on fluids at supercritical stage, that is at operating pressures of the order of hundreds of bars. 5. Conclusions The modest efficiency of impeller agitated jars was already perceived – more or less explicitly – by a number of authors, that showed there is room to improve the current state of the art on rotating devices for mixing engineering, whose actual efficiency in terms of fine turbulence and its spatial homogeneity may be questioned. This study has emphasized, from the fluidodynamic point of view, the reasons why Rushton rotating paddles have a disappointing mixing performance: to be more specific, the formation of stagnate/slow recirculating zones and the highly localized, even though exorbitant, shear rates. The global mixing efficiency is hardly improved by the impeller’s geometry, as Fig. 3 proves. In spite of the removal of these limitations by the proposed new agitator (Fig. 5a shows the best configuration, where rotating paddles sweep away any stagnate zone), the numerical simulation showed that no combination of parameters is capable of avoiding the formation of a central slow motion zone close to the shaft and the creation of a peripheral turbulent spot (clearly visible in Figs. 8 and 9); although these spots are of serviceable turbulence, potentially useful for applications, still the spatial homogeneity of the -scale is not adequate; thus many revolutions are needed to achieve the safe target level of homogeneity. It is concluded that rotating impellers seem intrinsically not efficient since the rotation itself amplifies the limitations of the apparatus: so there are reasons to reconsider the actual convenience of rotating stirrers in terms of energy consumption. Although cavitation, buoyancy effects, chemical reactions and two-phase flows were not considered in the numerical simulation, it is felt that the present conclusions should be barely affected by these phenomena. According to [1], the number and size of the holes control turbulence intensities downstream of 2D grids; this study showed that the above 2D result cannot be extended to revolving perforated plates in cylindrical domains and this is consistent with the dissipation features of perforated plates illustrated by Malavasi et al. [2].
As the concluding remark, it seems appropriate to mention oscillating translating grids as an alternative technique in place of rotating grids. The translational vibrating grid is a well established technology in water treatment plants, particularly in coagulation and flocculation plants (see, among many others, Ref. [30]): in this case, the grid is plane and typically vibrating at small amplitudes and large frequencies. Liem et al. [30] demonstrated the great potential of vibrating grids to produce fine and homogeneous turbulence; their exhaustive review [30] also confirmed that rotating impellers produce significant spatial variations of turbulent quantities and that different types of impellers did not reveal any major difference among each other; this point has been furtherly corroborated in this study (Fig. 3). Acknowledgments The author is grateful to Dr. G. Bazzanella for performing Ansys simulations. The comments of the three anonymous Referees are gratefully acknowledged. The kind permission of Elsevier and of Professor Andrè Bakker to reproduce Figs. 2 and 4 is appreciated. Appendix A. Details of numerical simulations Pn = paddles’ number; Hn = holes’ number; dh = holes’ diameter; ω = angular velocity; ¯ = weighted Kolmogorov scale; G = Gini coefficient; = G · ¯ (Eq. (6)); = standard deviation; Sk = skewness; Ku = kurtosis; n/a = not applicable. Pn
Hn
dh [mm]
ω [s−1 ]
¯ [m]
2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
30 30 30 30 30 30 20 20 20 20 20 20 12 12 12 12 12 12
0.5 0.5 1.0 1.0 2.0 2.0 0.5 0.5 1.0 1.0 2.0 2.0 0.5 0.5 1.0 1.0 2.0 2.0
100 50 100 50 100 50 100 50 100 50 100 50 100 50 100 50 100 50
66 107 66 108 59 91 64 110 66 107 62 95 69 116 66 109 64 102
0.88 0.82 0.87 0.83 0.85 0.86 0.87 0.82 0.88 0.84 0.88 0.87 0.85 0.78 0.87 0.81 0.89 0.86
4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4
30 30 30 30 30 30 20 20 20 20 20 20 12 12 12 12 12 12
0.5 0.5 1.0 1.0 2.0 2.0 0.5 0.5 1.0 1.0 2.0 2.0 0.5 0.5 1.0 1.0 2.0 2.0
100 50 100 50 100 50 100 50 100 50 100 50 100 50 100 50 100 50
38 63 40 63 33 53 41 68 40 65 38 56 42 71 41 68 40 65
0.94 0.93 0.95 0.93 0.97 0.93 0.94 0.92 0.94 0.93 0.93 0.94 0.94 0.92 0.94 0.93 0.94 0.93
8 8
30 30
2.0 3.0
100 100
270 180
0.95 0.94
2 2 4 4
0 0 0 0
n/a n/a n/a n/a
100 50 100 50
74 122 48 80
0.82 0.74 0.95 0.92
G
[m]
[m]
Sk
58.1 87.7 57.4 89.6 50.1 78.3 55.7 90.2 58.1 89.9 54.6 82.7 58.7 90.5 57.4 88.3 56.7 87.7
20.5 43.8 23.8 40.8 23.8 33.8 21.9 43.2 22.3 39.2 23.8 35.0 24.7 45.0 22.4 39.0 22.2 36.7
−1.80 −0.59 −0.23 −0.54 0.27 −0.48 −0.28 −0.53 −0.30 −0.46 0.05 −0.48 0.05 −0.39 −0.25 −0.48 −0.29 −0.59
5.78 3.58 3.37 3.26 2.81 3.33 3.91 3.45 3.60 3.46 3.32 3.59 3.68 3.27 3.53 3.50 3.70 3.62
35.7 58.6 38.0 58.6 32.0 49.3 38.5 62.6 37.6 60.5 35.3 52.6 39.5 65.3 38.5 63.2 37.6 60.5
8.4 16.0 7.8 14.5 6.2 13.4 9.0 17.6 8.0 14.4 7.8 13.3 9.4 18.3 8.2 15.6 7.9 15.5
−1.11 −1.48 0.33 −0.68 0.17 −0.46 −0.54 −1.20 −0.10 −0.82 −0.18 −0.64 0.11 −0.76 0.29 −0.58 0.37 −0.74
8.54 6.69 12.03 6.25 12.67 6.94 8.79 5.75 8.83 6.50 8.40 6.85 9.33 5.09 9.27 5.72 13.14 6.29
48 34
−1.88 −1.23
14.23 7.26
23.5 40.5 9.8 17.9
0.32 −0.07 0.84 −0.34
3.32 3.13 10.17 4.85
257 169 60.7 90.3 45.6 73.6
Ku
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
References [1] R. Liu, D.S.K. Ting, Turbulent flow downstream of a perforated plate: sharpedged orifice vs. finite-thickness holes, J. Fluids Eng. 129 (2007) 1164–1171. [2] S. Malavasi, G. Messa, S. Macchi, The pressure drop coefficient through sharpedged perforated plates, in: 32◦ Hydraulics and Hydraulic Structures National Congress, Palermo, Italy, 2010. [3] V. Buwa, A. Dewan, A.F. Nassar, F. Durst, Fluid dynamics and mixing of singlephase flowing a stirred vessel with a grid disc impeller: experimental and numerical investigations, Chem. Eng. Sci. 61 (2006) 2815–2822. [4] T. Kumaresan, B.J. Jyeshtharaj, Effect of impeller design on the flow pattern and mixing in stirred tanks, Chem. Eng. J. 115 (2006) 173–193. [5] V.W. Uhl, J.B. Gray, Mixing: Theory and Practice, Orlando Academic Press, 1986. [6] B.S. Choi, B. Wan, S. Philyaw, K. Dhanasekharan, T.A. Ring, Residence time distributions in a stirred tank: comparison of CFD predictions with experiment, Ind. Eng. Chem. Res. 43 (20) (2004) 6548–6556. [7] A.S. Khare, K. Niranjan, An experimental investigation into the effect of impeller design on gas hold-up in a highly viscous Newtonian liquid, Chem. Eng. Sci. 54 (1999) 1093–1100. [8] L. Cortellezzi, A. Adrover, M. Giona, Feasibility, efficiency and transportability of short duration optimal mixing protocols, J. Fluid Mech. 597 (2008) 199–231. ´ Control of mixing in fluid flow: a maximum [9] D. D’Alessandro, M. Dahleh, I. Mezc, entropy approach, IEEE Trans. Automat. Contr. 44 (1999) 1852–1863. ´ G. Tadmor, A. Banaszuk, Optimal mixing in recirculation [10] B.R. Noack, I. Mezic, zones, Phys. Fluids 16 (2004) 867–888. [11] M.M. Alvarez-Hernandez, T. Shinbrot, J. Zalc, F.J. Muzzio, Practical chaotic mixing, Chem. Eng. Mixing 57 (2002) 3749–3753. [12] L. Oshinowo, A. Bakker, E. Marshall, Simulating Mixing Time with Computational Fluid Dynamics, vol. JA111, Fluent Inc. Press, 2000, pp. 1–5. [13] P.T.L. Koh, J.G.R. Andrews, P.T.H. Uhlherr, Flocculation in stirred tanks, Chem. Eng. Sci. 39 (1984) 975–985. [14] H. Wu, G.K. Patterson, Laser-Doppler measurements of turbulent flow parameters in a stirred mixer, Chem. Eng. Sci. 44 (1989) 2207–2221.
809
[15] L.A. Cutter, Flow and turbulence in a stirred tank, AIChE J. 12 (1966) 35–45. [16] S.M. Kresta, K.J. Bittorf, D.J. Wilson, Internal annular wall jets: radial flow in a stirred tank, AIChE J. 47 (2001) 2390–2401. [17] I.P.D. De Silva, Turbulence dissipation in stirred jars, J. Eng. Mech. ASCE 132 (11) (2006) 1260–1268. [18] R.M. Smith, L.X. Liu, J.D. Litster, Breakage of drop nucleated granules in a breakage only high shear mixer, Chem. Eng. Sci. (2010) (published online). [19] A. Arbel, A. Shklyar, Numerical simulations of turbulent flow through fine screen, in: MASCOT09, IMACS Workshop, IAC-CNR, Rome, Italy, 2009. [20] D.A. Deglon, C.J. Meyer, CFD modeling of stirred tanks: numerical considerations, Miner. Eng. 19 (2006) 1059–1068. [21] I.B. Celik, U. Ghia, P.J. Roache, C.J. Freitas, H. Coleman, P.E. Raad, Procedure for estimation and reporting of uncertainty due to discretization in CFD applications, J. Fluids Eng. ASME 130 (7) (2008) 1–4. [22] I. Celik, O. Karatekin, Numerical experiments on application of richardson extrapolation with nonuniform grids, J. Fluids Eng. ASME 119 (1997) 584–590. [23] J.H. Ferziger, M. Peric, Further discussion of numerical errors in CFD, Int. J. Numer. Methods Fluids 23 (1996) 1263–1274. [24] W. Bartok, C.E. Heath, M.A. Weiss, Mixing in a jet stirred reactor, AIChE J. 6 (4) (1960) 685–687. [25] J.A. Sanchez Perez, E.M. Rodriguez Porcel, J.L. Casas Lopez, J.M. Fernandez Sevilla, Y. Chisti, Shear rate in stirred tank and bubble column bioreactors, Chem. Eng. J. 124 (2006) 1–5. [26] G.B. Schubauer, W.G. Spangerberg, P.S. Klebanoff, Aerodynamic characteristics of damping screen, in: Tech. Note 2001, NACA, Washington, DC, 1950. [27] M.O. Lorenz, Methods of measuring the concentration of wealth, Q. Publ. Am. Stat. Assoc. 9 (70) (1905) 209–219. [28] C. Gini, Measurements of inequality and income, Econ. J. 31 (1921) 124–126. [29] J. Ding, X. Wang, X. Zhou, N. Ren, W. Guoa, CFD optimization of continuous stirred-tank (CSTR) reactor for biohydrogen production, Bioresour. Technol. 101 (18) (2010) 7005–7013. [30] L.E. Liem, D.W. Smith, S.J. Stanley, Turbulent velocity in flocculation by means of grids, J. Environ. Eng. ASCE 125 (3) (1999) 224–233.
Contents lists available at ScienceDirect
Chemical Engineering and Processing: Process Intensification journal homepage: www.elsevier.com/locate/cep
On the efficiency of turbulent mixing in rotating stirrers Filippo Trivellato ∗ Department of Civil and Environmental Engineering, University of Trento, I-38123 Trento, Italy
a r t i c l e
i n f o
Article history: Received 7 October 2010 Received in revised form 1 February 2011 Accepted 26 May 2011 Available online 23 June 2011 Keywords: Gini Kolmogorov Mixing Rushton Sterilization
a b s t r a c t Turbulent mixing is crucial in sterilizing engineering, where it is mandatory to produce the finest and the most homogeneous turbulence all over the fluid domain. This study shows how the turbulence imparted by Rushton impellers, which are typically used in several industrial productions, may be hardly improved by varying the impeller’s geometry. An impeller of new design, a perforated paddle that completely removes both limitations of Rushton impellers, is illustrated. The mathematical modeling is based on the RANS equations and the k–ε model for turbulence closure; great care was lavished to understand the accuracy of the numerical solution. The results are presented in terms of a global indicator of mixing performance that combines both the Kolmogorov scale and the Gini coefficient. The numerical solution disclosed the formation of a peripheral turbulent spot and a seemingly slow motion zone next to the shaft of the rotating perforated paddle. The imparted turbulence may be usefully fine, but many revolutions are needed to grant the safe target level of spatial homogeneity of the turbulence. This was confirmed by both a dedicated statistics treatment and the conventional statistical analysis. © 2011 Elsevier B.V. All rights reserved.
1. Introduction Mixing of single phase flows in mechanically agitated vessels is important in countless industrial processes: colloidal mills, emulsifiers, homogenizers and sterilizers for liquid food are among the industrial activities that can possibly benefit from this study. In mechanically agitated vessels at constant temperature, turbulence (or, in other words, the mixing intensity) is controlled by a variety of factors (the physical properties of the fluid under treatment; the power input, e.g. the impeller rotational speed; and the vessel-impeller design). Scores of different impellers are persistently being developed and an overwhelming literature is available that describes every detail of these impellers; this study is confined to rotating impellers which are broadly classified as radial (such as the Rushton turbine) or axial (propeller and pitched paddle turbine): both kinds of impellers revolve around the central shaft, that identifies the symmetry axis of a cylindrical vessel (Fig. 1). The impeller’s scope is to efficiently produce the most turbulent flow by injecting the least energy to the fluid to be treated: but how much turbulent is the imparted mixing? And how much spatially homogeneous is the generated turbulence? The present study aims at understanding whether revolving impellers (both radial and axial) can actually grant the safe target
∗ Tel.: +39 0461 282 631; fax: +39 0461 282 672. E-mail address: [email protected] 0255-2701/$ – see front matter © 2011 Elsevier B.V. All rights reserved. doi:10.1016/j.cep.2011.05.012
level of turbulence that is needed in sterilizing engineering together with a homogeneous intensity all over the flow field. Cavitation, buoyancy effects, chemical reactions and two-phase flows will not be considered herein. Turbulence that is approximately homogeneous and isotropic is typically produced by flows at high Reynolds numbers through screens, grids or perforated plates. Downstream of these devices, strong shear flows are created by the interactions among the jets issuing from the openings and the wakes formed behind the solid part of the devices. Beyond a distance, the flow is fully developed and the turbulence is basically isotropic. In describing homogeneous turbulence at high Reynolds numbers, it is assumed that energy is injected into eddy motions on a scale of order the size of the vessel; energy is then transferred to smaller and smaller length scales (the well known energy cascade) until it is dissipated at the smallest stable viscous scale, the Kolmogorov scale (hereinafter, the -scale). Turbulent flows with small length scales (fine grain turbulence) have the greatest rate of dissipation and their longitudinal decay is the quickest; turbulence intensities downstream of grids in 2D flows were shown [1] to be controlled by the number and size of the holes; however this fact is controversial, as claimed by the recent and comprehensive review on the dissipation features of perforated plates inserted transversally in cylindrical pipes [2]. It is pointed out that the above review is only relevant to 2D flows and its extension to rotating impellers in cylindrical domains will be investigated in this study.
800
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
of the vortical structures at which viscous dissipation occurs; the Kolmogorov scale was chosen because it conveys an immediate physical meaning and it is also helpful in dictating the smallest scale the mesh elements must have to obtain a meaningful numerical solution. 2. Comparison of impellers
Fig. 1. Slow recirculating flow patterns (invariant tori) in Rushton and disk turbines.
To understand the mixing capability of a stirrer, the temporal evolution of a tracer is checked at convenient locations within the liquid domain [3]; the mixing time is the time needed for the tracer to reach the value of the final equilibrium concentration within ±5% [6]) (or ±2% according to [4]); the average of mixing times calculated at different monitoring locations and the spread of mixing times about the averaged value give a measure of the intensity and homogeneity of the whole mixing. This procedure is quite involved and time consuming, it does not provide a detailed knowledge on the actual overall mixing time and the actual overall quality of mixing. In particular, the application of the mixing time concept to Rushton impellers cannot show a variety of features of the flow field (for example, anticipating the present results, the most turbulent spots of Figs. 8 and 9) and can only suggest improvements of impeller’s geometry by the trial and error method; this will be confirmed in Section 2. And yet the chief limitation is that the mixing time concept is intrinsically dependent on the vessel size [5] and scaling-up may turn out questionable. Khare et al. [7] considered the gas hold-up as a measure of the efficiency of gas–liquid contacting, but this will be discussed at length in Section 2. Apart from conceiving a criterion to characterize the mixing intensity, there is also the need to quantify the quality of mixing; this was done by studying the variance decay exponent or by the mix-norm (see, for instance, Ref. [8]). D’Alessandro et al. [9] quantified mixing quality in terms of entropy and Noack et al. [10] by means of Poincare maps; unfortunately both approaches are not particularly suitable for applied engineering modeling, are not of general value and their results can hardly convey handy implications to improve rotating impellers. The mixing parameter adopted in this work will be the Kolmogorov scale, Eq. (4), which is the minimum stable length scale
As it has been plainly stated by Alvarez-Hernandez et al. [11], mixers are among the most expensive and inefficient equipment in production plants; as far as radial and axial turbines are concerned, it was convincingly demonstrated that the formation of recirculating cells (also known as invariant tori, Fig. 1) is an obstacle to fluid transport. Slow recirculating flow patterns are invariably developed in the stagnate (dead) zones of the vessel that must be eliminated as they are a barrier to fine grain mixing; this translates into questioning the real utility of introducing more and more turbine paddles, which is the customary industrial inclination (this topic will be also substantiated in Fig. 9). Since these plants require many revolutions to achieve the safe target level of homogeneity, energy consumption may turn out as a matter of concern. The removal of stagnate zones is therefore mandatory not only to reduce mixing times, but also to obtain the complete homogenization all over the liquid domain, which is of paramount importance in many industrial processes: the most meaningful example in this context is the sterilization treatment in the liquid food industry. To address the above homogenization issue, the effort of the literature was ceaselessly devoted at improving the impeller’s geometry, mostly by the trial and error method; unfortunately, this approach is far from being successful, as the following elaboration proves. The performance of six impellers (four disks and two pitched paddles) were compared by Khare and Niranjan [7], who measured the gas hold-up εf , that was deemed by the authors as the overall measure of the efficiency of gas–liquid contacting. The hold-up data were analyzed by separating the contributions of large and tiny bubbles from the total gas hold-up. Fig. 2 shows two distinct pictures of the final steady state values of the total hold-up εf as a function of the power dissipation per unit volume of the initially undisturbed liquid. The experiments were intended to assess the most efficient geometry of the impeller, namely the impeller that is capable of producing the highest hold-up at the lowest input power. Fig. 2 clearly shows that the peak values of hold-ups do not differ significantly among the six tested impellers, being all the peak values comprised within the
Fig. 2. Steady state gas hold-up εf produced by six different impellers: (a) large bubbles; (b) tiny bubbles. Reprinted from Khare et al. [7] with kind permission from Elsevier.
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
801
Fig. 3. Present elaboration of Fig. 2: (a) large bubbles; (b) tiny bubbles.
narrow range 0.067–0.076 for large bubbles, and 0.13–0.17 for tiny bubbles; at the same time, the overall behaviour of the six tested impellers is unquestionably different. Apart from the fact that the gas hold-up is a basic measure of mixing efficiency, the graphs of Fig. 2 give no indication about the time required to attain the claimed steady state peak values; in other words, it is not known the amount of energy employed to reach the above steady values of hold-ups and energy consumption is unmistakably relevant to any industrial process. The following present elaboration is useful in addressing this energy issue: Fig. 3 was derived in this study from the data of Fig. 2 by taking the specific hold-up ε¯ f in the ordinate axis (the hold-up per unit input power per unit volume of treated fluid). For both large and tiny bubbles, all the impellers’ performances of Fig. 2 nicely collapse in one single curve, showing that there are vanishing differences among the tested impellers from the energy consumption point of view. Even assuming that εf might be considered as an indicator of mixing performance, apparently the energy required to obtain the target mixing intensity does not differ among the six tested impellers, even though they have different peak power values. Therefore variations of impellers’ geometry may only produce modest, if any, advancements in terms of energy savings. One further problem with impeller stirred vessels is that shear rates are strongly localized at the edge of the rotating paddles and their decay is extremely quick, being already vanishing a short distance away from the rotating paddle; this is readily seen in Fig. 4, where realistic wakes downstream of a six paddle Rushton impeller were calculated by Oshinowo et al. [12]. Many studies agreed that most of the energy (80% according to Koh et al. [13]; 60% by Wu and Patterson [14]; 70% by Cutter [15];
89% by Kresta et al. [16]; 62% by De Silva [17]) is dissipated near the impeller. To provide proper mixing to the most distant zones of the tank, a great deal of energy needs to be imparted; nonetheless, this is still far from being sufficient, especially on scale-up [5], so that stagnant zone of fluid, as a matter of fact, always exists far away from the paddle in spite of the enormous shear gradient delivered by the impeller. This is not the desirable flow field to be realized in a mixing jar for bacteria sterilization. It is also worth to mention that some industrial processes do not allow exceedingly high rotation rates due to the possible rupture of cell membranes or to possible damage of micro-organisms and sludge flocs (incidentally, high shear granulation is used in pharmaceutical, detergent and food industries where granule breakage is purposely performed [18]). The above review of the literature dealing with revolving stirrers has disclosed the following points that deserve further investigation: (i) modifications of the Rushton impellers’ geometry may yield negligible consequences on the mixing performance owing to both the creation of stagnate zones and the extreme localization of shear rates; (ii) mixing produced by grids (or perforated paddles) needs a more handy characterization; and (iii) the spatial homogeneity of the turbulence is as important as its maximum intensity. The present work aims at addressing the three above mentioned issues by the numerical simulation of a new paddle, presented in Section 3 along with a discussion on its accuracy, and at proposing the -scale as a handy indicator to qualify mixing performance (Section 4); in view of the expected strong heterogeneity, a new statistical criterion will be tested and compared to the conventional statistical analysis to address the spatial homogeneity issue (Section 4.1). 3. Numerical simulation 3.1. Batch apparatus
Fig. 4. Localization of shear rates at the outer edge of Rushton paddles rotating clockwise in a baffled tank [7]. Reproduced by kind permission of A. Bakker).
The numerical domain simulated in this study (Fig. 5) comprised a batch-operated cylindrical apparatus with perforated paddles rotating around the shaft (the symmetry axis of the cylinder) and sweeping all the internal volume. This apparatus removes both of the drawbacks discussed above (slow recirculating cells and localization of shear rates). The internal diameter of the cylindrical vessel, D = 106 mm, is representative of an actual experimental apparatus; the paddles’ diameter was di = 100 mm, leaving a 3 mm gap between the rotat-
802
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
Fig. 5. (a) The simulated batch apparatus (2 paddles, in this case); (b) the discretized domain (one quarter of a cylinder, in this case).
ing paddles and the fixed lateral wall; the same 3 mm clearance was also imposed between the paddles and both the ceiling and bottom floor; the paddle thickness was 1 mm and the diameter of the rotating shaft was 20 mm. The height of the paddle was some one quarter of the cylinder radius, namely 13 mm. The role of the holes (more specifically, their number, diameter and distribution) is controversial, as it has been thoroughly discussed by Malavasi et al. [2]. Although the number and diameters of the holes were systematically varied in the present simulations, the holes were always arranged in a staggered pattern that remained unchanged in all the simulations and is depicted in Fig. 6. The holes were modeled as cylindrical in shape with a sharp edge, since this is the customary and easiest way of construction. It is well known (see, for example, Refs. [1,2]) that the edges of orifices have a role in turbulence generation; however, in view of the unusefully difficult task of creating a complicated mesh on a rounded edge, it is firmly believed that even a less efficient edge would serve the purpose of this study: a systematic slight underestimation of turbulence parameters is expected, but this does not affect the general validity of the speculations that will be drawn from the present results. From the practical point of view, to minimize the drag of revolving plates (thus saving input power) it is mandatory to have a great number of holes, but this can pose structural problems to the plate which has to be thin. 3.2. Mesh generation and boundary conditions Taking advantage of the spatial periodicity of the cylindrical domain, only the region of fluid between any two consecutive paddles was actually discretized: this means one quarter of the cylinder or one eight and so on, depending on the adopted number of the rotating paddles. The numerical model of the stirred tank is based on an unstructured finite volume mesh partitioned in two regions: (i) a rotating mesh in the impeller region, in contact with the paddle
and moving with it; (ii) a fixed mesh elsewhere (more specifically, in the narrow gap comprised among the rotating paddles and the vessel walls) as depicted in Fig. 5. The two meshes slide past each other in an unsteady manner; the mobile mesh and the fixed mesh are separated by a mobile cylindrical interface, that provides the appropriate treatment of the continuity between the two meshes. At the mobile interface, a conservative interpolation is used for both mass and momentum equations using a set of fictitious control volumes; all the calculated fluid motion strictly arises from the paddles rotation. The penalty with the sliding mesh is the CPU time, which is an order of magnitude higher than the comparable steady state calculation. Within the rotating domain, one further interface (called volume interface) needs also to be introduced, but this is only for numerical purposes. Fig. 6 shows the mesh refinement on the paddle, on the sharp-edged holes and on the lateral gaps; the minimum and the maximum size of the cells varied from 0.01 mm through 0.5 mm. The cells’ volumes varied from 10−4 mm3 through 0.22 mm3 . The total number of cells varied depending on the simulation and it was typically of order 106 . In all the present numerical experiments, the stirred liquid was water at 20 ◦ C. The fixed lateral walls of the tank were assigned standard wall function boundary conditions; both the top and bottom walls were treated like symmetry boundaries; the surface of the moving zone was assigned an interface boundary condition; and the surface of the impeller was assigned a rotating wall boundary condition. The quality of the mesh is of paramount importance to obtain a credible numerical solution; this is particularly true for the nodes close to wall boundaries due to the steep velocity gradients; following [19], all the nodes at the wall lie within the viscous sublayer (u* y/ < 4, where u* is the friction velocity, y is the distance of the node from the wall and is the kinematic viscosity); again this is a good indication of the satisfactory resolution of the wall
Fig. 6. Mesh refinement on the paddle and on the lateral gaps.
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
803
where S is the modulus of the mean rate of strain tensor Sij ; Sk and Sε are source terms. The constants of the k–ε model that were used in the present study are the universally adopted values: c1ε = 1.44,
Fig. 7. The mobile interface poses a Courant limitation on time marching.
boundary layer and hence of the overall accuracy of the numerical simulation. The most critical boundary condition is the interface connecting the fixed mesh to the rotating mesh; the interface prompts a Courant limitation on the marching time step to negotiate the appropriate exchange of the physical information (Fig. 7): more specifically, it was found that the rotating black inner cell of Fig. 7 must be in contact with the fixed grey outer cell for at least 3–4 time steps, depending on the imposed angular velocity. 3.3. Turbulence closure The flow field of the rotating apparatus is unsteady, 3D and is governed by the Reynolds-Averaged Navier–Stokes (RANS) equations, which are far too known to be reported here; the standard k–ε turbulence model has been widely applied in the literature to close the RANS equations. Different versions of the k–ε model (low Reynolds number model; Reynolds-stress transport model) were also tested in the present work; however, vanishing differences were found among these models and this is consistent with the present literature consensus. Deglon and Meyer [20] proved that the partially satisfactory predictions due to the use of k–ε turbulence models may be due to numerical errors rather than to inherent inadequacies of the turbulent model: thus the standard k–ε closure model was herein adopted because it proved to be the best compromise between accuracy of the results and computational cost. Adopting the convention on repeated indices, the k – equation and the ε – equation for incompressible Newtonian fluid are reported below: ∂ ∂ ∂ (k) + (kui ) = ∂t ∂xi ∂xj ∂ ∂ ∂ (εui ) = (ε) + ∂t ∂xi ∂xj
+
t k
t + ε
∂k ∂xj
∂ε ∂xj
+ Gk − ε + Sk
ε2 ε + Sε + c1ε Gk − c2ε k k (1)
where is the molecular viscosity of the fluid (a physical property) and t is the turbulent viscosity (which is instead a function of the flow field, i.e. of the imparted mixing): t = c
k2
c = 0.09,
k = 1.0,
and
ε = 1.3
The RANS equations were solved using the commercial code Ansys Fluent 12.0. A first order scheme was used to discretize the spatial derivative terms in the governing equations, while a second order scheme was used for the discretization of time derivatives. The under-relaxation parameters for velocity, pressure, turbulent kinetic energy and dissipation energy were 0.3, 0.6, 0.6 and 0.6 respectively. All the computations were performed until the residuals of mean velocities, pressure and turbulence quantities were less than 10−5 . A typical numerical run, involving one million cells, took roughly one day of CPU time to achieve the converged solution in a 2.3 GHz quad-core processor.
3.4. Numerical uncertainty Although there is no shortage of CFD simulations of mixing problems, a thorough discussion on the actual numerical accuracy is hardly found in the literature; therefore, comprehensive sensitivity tests were performed in this study to check whether the present numerical solution was truly unaffected by the adopted mesh; these sensitivity tests are briefly reported in this section. The procedure followed to estimate the numerical errors generated by the present simulation is based on the valuable ASME guideline [21]. The underlying perspective here is to understand whether the present study is capable of predicting with trustworthy accuracy the outcome of mixing numerical simulations in cases in which the measure of turbulent data is not an easy task, as in the present case. Uncertainty in inlet boundary conditions is well known to be a potentially significant contributor to the overall uncertainty; the degree of sensitivity of the present numerical solution to perturbations in inlet conditions was checked by varying the inlet turbulence level from 1% to 30%. The inlet boundary condition was imposed on the holes of the paddle and on the lateral gaps; the flow field is strongly accelerating downstream of the inlet boundary, and then it diverges, slowing down and interacting with the nearby jets. The resulting shear flow was observed to be not affected (less than 1%) by the above variation of the inlet turbulence level; this is reasonable, since accelerating zones are known to dissipate turbulence and this fact considerably mitigates the uncertainty due to the inlet turbulence. The method herein followed for the evaluation of discretization errors is the Grid Convergence Index (GCI). The first step of the procedure is the definition of the representative grid size h:
h=
1 Vi N N
1/3 (2)
i=1
ε
The term incorporating the Boussinesq hypothesis, Gk , describes the generation of turbulent kinetic energy due to the base flow gradients and is given by: Gk = t S 2 ,
c2ε = 1.92,
S 2 = 2Sij Sij =
∂2 Ui Uj ∂Ui ∂Ui + ∂xj ∂xj ∂xi xj
where Vi is the ith cell volume and N is the total number of cells used for the computations. Three significantly different grids were selected and values of the Kolmogorov scale (which is the key variable here) were calculated. In this study, the grid refinement factor r = hcoarse /hfine was always greater than 1.3 (a value based on experience [21]). Roughly geometrically similar cells were always used.
804
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
The apparent order p of the method was calculated by using the following expression: p=
|ln|ε32 /ε21 | + q(p)| ln(r21 )
q(p) = ln
p
r21 − s
(3)
p
r32 − s
s = 1 · sgn
ε 32
ε21
where h1 < h2 < h3 , r21 = h2 /h1 , r32 = h3 /h2 , ε32 = 3 − 2 , ε21 = 2 − 1 , k stands for the solution derived by means of the kth grid (k = 1, 2, 3); sgn is the function signum. The absolute value in Eq. (3) is necessary to ensure extrapolation toward h = 0 [22]. Negative values of the ratio ε32 /ε21 were found in the present simulations; this is an indication of oscillatory convergence, a typical problem occurring when dealing with geometrically complicated fluid domains; however, it is noted that these oscillations were not significant, being their maximum values less than 3%. A fair agreement was observed between the apparent order of the scheme with the formal order used and this is an indication of the grids being in the asymptotic range. 21 were then calculated along with The extrapolated values ext the approximate relative error ea21 : p
21 =
ext
r21 1 − 2 p
r21 − 1
1 − 2 ea21 =
1
The fine-grid convergence index GCI21 fine is given by: GCI21 fine =
1.25ea21 p
r21 − 1
The local order of accuracy p calculated from Eq. (3) was between 0.2 and 3.2, with an average value of pave = 1.5 ; oscillatory convergence was found in 10% of the nodes. The maximum discretization error was some 30%, which is admittedly a large error; however, it must be emphasized that this value was produced close to the boundary layer and translates in a maximum uncertainty of the local velocity of about 0.004 m/s; from a practical point of view, this error cannot possibly be of any relevance as far as the overall accuracy is concerned. Following [23], the iteration error was estimated by
niter,i ∼
in+1 − in i − 1
where n is the iteration number and i is the principal eigenvalue of the solution matrix of the linear system, which can be approximated as: i ∼
|| in+1 − in || || in − in−1 ||
where || · || is an appropriate norm (e.g., the infinity norm L∞ ). The uncertainty ıiter in iteration convergence can then be estimated as: ıiter ∼
described method was found to be typically almost one order of magnitude smaller than the discretization error; so there are reasons to trust the accuracy of the present numerical simulation.
|| niter,i || ave − 1
where ave is the average value of i over a reasonable number of iterations. The iteration error calculated according to the above
4. Results The following parameters were herein tested: (i) number of paddles: 2, 4 and 8 (i.e., 4, 8 and 16 half-paddles); (ii) number of holes in each half-paddle: 0 (i.e., solid plate), 6, 12, 20 and 30; and (iii) diameters of the sharp-edged holes: 0.5, 1, 1.5 and 2 mm. The angular velocity of the shaft varied: Appendix A reports a selection of numerical runs performed at 50 and 100 rad/s. The maximum number of possible combinations of the above set of parameters is huge; however, during the study, it was clear that many combinations were not worth to be tested because of vanishing differences among them. Only 220 selected runs were actually performed; and yet, since most of them turned out to produce nothing but a slight variation in the results, for brevity’s reasons only meaningful results will be actually illustrated hereinafter. The interested reader may find on Appendix A all of the details useful to reproduce the present numerical experiments. The purpose of this section is to quantify turbulent mixing rather than giving sharp indications to pursue the optimal construction of the vessel. The results are presented in terms of the Kolmogorov scale, the -scale:
=
3 ε
1/4 (4)
where is the kinematic viscosity and ε is the turbulent dissipation energy, whose evolution is governed by Eq. (1). The -scale is herein proposed as an easy, and yet effective, parameter to deal with, that dictates indications on the size the mesh elements must have to obtain a meaningful numerical solution. Fig. 8 gives the picture of the -scale distribution in the midplane intersecting the central row of the holes (i.e., at half height of the paddle; the flow is counterclockwise, from below to left). The column of Fig. 8a shows the effect that the variation of the hole diameter produces on the -scale at an angular velocity of 100 rad/s and with 2 paddles (4 half-paddles); the column of Fig. 8b shows the effect of the variation of the holes’ number at 100 rad/s and 2 paddles. There are no stagnate zones. Recalling that di is the impeller diameter and N is the number of revolutions per second, the global Reynolds number: Re =
Ndi2
(5)
is representative of the turbulent regime, being greater than 4 × 104 in all of the present numerical runs. As a first observation, it seems that no dramatic diversification of the results is apparent. The second impression is that the most beneficial and serviceable zone, i.e. the brightest spot with the lowest -scale, is always localized by centrifugal forces at the top peripheral region; the area of this zone is limited in extent, and never exceeding 15% of the entire domain. The turbulence produced in the narrow gap existing among the rotating paddles and the lateral walls of the tank is fine grain turbulence similar to the turbulence produced behind the holes of the rotating paddles. The central region, close to the shaft, behaves basically like a slow motion zone, and this is consistent with what already observed by Bartok et al. [24] in a similar rotating apparatus; it was stated there, and confirmed here, that the mixing pattern in the
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
805
Fig. 8. Distribution of the -scale as a function of the holes’ diameter (column a) and of the holes’ number (column b). Counter clockwise flow from below to left; darkest spots: 1.4 × 10−4 m; brigthest spots: 10−6 m.
central core is modest and almost independent from the angular velocity. Perforated plates may operate either as turbulence suppressors or as turbulence generators. The turbulence generated downstream of a perforated plate is the combined result of three distinct phenomena: the passage of the upstream turbulence through the paddles’ openings; the shear turbulence produced downstream by the interactions of the jets exiting the screen; and the turbulence produced by the perforated plate itself. In the generator role, the plate has a coarser pattern of holes and is operated at high Reynolds numbers. It is known that the turbulence scale imposed initially at the screen controls the turbulence intensity (and below some minimum initial scale, the turbulence always decays); for example, for sufficiently fine plates, the turbulence is of small scale and decays
after 500 hole diameters. At roughly 40 hole diameters, the turbulence is nearly homogeneous and isotropic. The above decay law downstream of grids is only known for 2D flows; it is emphasized that no results are available for axialsymmetric rotating domains. Fig. 9 shows the effects of the variation of the paddles number and of the angular velocity; increasing the number of paddles does not produce finer turbulence, a result possibly unexpected a priori since increasing the number of blades is a popular method used in industries. Interestingly enough, Fig. 9 shows that a 50 rad/s angular velocity produces a noticeable zone of finer turbulence and more extended than in the 100 rad/s run; even the gap between the rotating paddle and the fixed wall of the vessel also plays a better role in producing fine turbulence. That set-up combination turned out to be the most efficient in terms of fine turbulence
806
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
Fig. 9. -Scale distribution as a function of the paddles’ number and the angular velocity. Counter clockwise flow from below to left; darkest spots: 1.4 × 10−4 m; brightest spots: 7 × 10−6 m.
within the tested range of experiments (but nothing can be said at this stage about its homogeneity); it is recalled that providing the optimal setting configuration is beyond the reach of the study. The analysis of all the results (including those not presented explicitly here) showed that the revolving velocity of the paddles (in other words, the power input) is an important parameter in generating the smallest -scale. This result is consistent with what was obtained by Sanchez Perez et al. [25], who enforced the idea that the shear rate in turbulent regime is a function of the angular velocity. Schubauer et al. [26] found that a single screen causing a strong pressure drop is a more satisfactory turbulent promoter than using several screens in series. This evidence is clearly seen in Fig. 9 when comparing the 2- and 8-paddle results: close to the shaft, there is a central zone of slow motion flow; this zone does not produce fine turbulence and is not valuable for sterilizing mixing. Thus the present apparatus, even though it completely removes the limitations described in the introduction (slow recirculating regions and localization of shear rates), still has a questionable performance due to the formation of both a peripheral turbulent spot and a central core of plug zone. And it seems clear that even the optimal combination of the parameters (which was not pursued here) would produce nothing but a slight improvement of the mixing intensity. This translates into the necessity to make many
revolutions to grant the required safe level: this issue is particularly mandatory in liquid food sterilization. In addition, the ideal mixing should produce in the flow field as much a uniform -scale as possible and a criterion is needed to qualify the spatial variability of the -scale. This is of particular relevance in the tested device which, as proved in Figs. 8 and 9, produces a marked spatial heterogeneity of the -scale. 4.1. Statistical treatment and Gini coefficient The numerical simulations described in Section 4 provided the -scale values in each cell of the domain; the range comprised between the minimum through the maximum value of was divided in 200 classes; Fig. 10 shows the frequency distribution of the local -scale values, where after the peak value, the function is typically monotonically decreasing. But the monotonic decrease was not observed in all of the present numerical simulations owing to the strong spatial heterogeneity of the -scale. Secondary peak values did occur that were located in between the first peak and the tail end of the distribution; these secondary peaks were at times as important as the first one. This fact has been tackled here by introducing the Lorenz’s function [27], which was originally intended to describe inequalities of the wealth distribution and is often used in biodiversity studies; it seems that the Lorenz function has not been previously applied in mixing engineering.
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
807
Fig. 10. Example of the frequency distribution of the local -scale.
Firstly, the data sample of the -scale was sorted in ascending order; then the cumulative frequency distribution was calculated by assigning the Weibull frequency to the ith value: Fi =
can be hardly smoothed away by any combination of the setting parameters.
i M+1
where M is the sample size. The cumulated amount of the ith value is then calculated as:
i
i =
4.2. The overall performance For conveniency’s sake, it may be useful to introduce the parameter as a new indicator of the overall performance:
j=1 j M j=1 j
= G¯
Fig. 11 shows an example of the i − Fi plotting: if all the values of the data sample were equal (in other words, if the -scale were perfectly distributed) the graphical representation would be the line of perfect equality, the line OA; by contrast, if only one value were different from zero, the line OBA would be the line of perfect inequality, meaning that the total amount of the -scale is wholly represented by just a single data of the sample. If data values are positive (as in the present -scale sample), the actual Lorenz curve must lie somewhere in between the two extreme lines of perfect equality and perfect inequality: the closer the Lorenz curve to the line OA, the more homogeneous the -scale. The Gini’s coefficient G [28] was conceived to quantify the inequality of data samples and it is defined as the ratio between the shaded area of Fig. 11 and the area of triangle OAB: G=
Fig. 12. The overall performance (number of holes: 6, 20 and 30; dh = 0.5, 1 and 2 mm).
shaded area area triangle OAB
The Gini coefficient is comprised between 0 and 1: if the -scale values were perfectly distributed, the Lorenz function would collapse on the line of perfect equality, the bisectrix OA, and G takes its minimum value (Gmin = 0); by contrast, the most unequal distribution is given by the line OBA and the Gini coefficient takes its maximum value (GMax = 1). The G coefficient calculated by the present simulations was comprised between 0.74 and 0.97, as reported in Appendix A, showing a great deal of spatial heterogeneity that
Fig. 11. Lorenz function and Gini coefficient.
(6)
where ¯ is the weighted mean of (where the weights are the volume of the cells). The parameter is intended to briefly summarize, in only one number, the global performance of the stirring device: the minimum value of is the desirable value for mixing engineering, since it combines at the same time both the finest mixing and the most homogeneous distribution of it. The minimum values of , in other words, detect the best combinations of the parameters. Amazingly enough, the range of is somewhat narrow (see Appendix A), indicating there is no dramatic difference in the results whichever the combination of the parameters. Fig. 12 displays the overall performance of the imparted mixing as a function of the Reynolds number, defined by Eq. (5); for clarity’s sake, only the maximum (2-paddle) and the minimum (4-paddle) values are presented. Fig. 12 shows the asymptotic behaviour starting at Re ∼ 105 (depending on the paddle number), meaning that the holes’ number and diameter are not parameters as relevant as it could be initially thought; by contrast, the number of paddles and the angular velocity are important parameters in this apparatus. It is emphasized that increasing the number of paddles does not produce finer turbulence, but it does promote a better homogenization; the combination of both effects are included in the parameter , that may therefore turn out lower for a 4-paddle impeller. Rather obviously, the simulations with solid paddles (i.e., with no openings) turned out to perform inefficiently in terms of mixing (details can be found in Appendix A). The analysis of the -scale values was also performed by means of the conventional statistics parameters: the weighted mean ¯ (where the weights are the volume of the cells); the standard deviation based on the weighted mean; the skewness coefficient, Sk ; and the kurtosis coefficient Ku . Again their definition is far too known to be reported here. Summing up briefly the main results, the statistical analysis based on the conventional parameters confirmed the strong heterogeneity of the flow field and it is worth to report here only the major differences from the above described (Section 4.1) dedicated statistics. More specifically:
808
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
(i) both the standard deviation and the parameter revealed the same behaviour on the involved parameters (paddles numbers, angular velocity and diameter of holes); (ii) most of the above simulations showed a negative skewness coefficient, indicating a significant concentration of -scale values lower than the average value (Appendix A); (iii) the 2-paddle turbine has the lowest kurtosis coefficient, meaning a rounded peak and short, thin tail of the distribution. Even though the present results are relevant to closed cylindrical containers with rotating perforated paddles, however the dependence of the pressure drop (not shown here for brevity’s sake) on relevant parameters – such as the global Reynolds number, the paddle thickness and the holes number and diameter – were found to be in qualitatively agreement with what described by Malavasi et al. [2], who tackled the different case of the perforated plate transversally inserted in a cylindrical conduit. Although the present simulations have been performed in a batch mixing apparatus, it is clear that the obtained results are just as valid, at least qualitatively, as in a semi-continuous or even a continuous mixing apparatus [29]. Also, all the literature studies which have been carried out on liquid and gas mass transfer in mechanically agitated reactors indicated that the operating ambient pressure does not influence the mass transfer coefficient; thus the present results can be extended even to micro-organisms inactivation based on fluids at supercritical stage, that is at operating pressures of the order of hundreds of bars. 5. Conclusions The modest efficiency of impeller agitated jars was already perceived – more or less explicitly – by a number of authors, that showed there is room to improve the current state of the art on rotating devices for mixing engineering, whose actual efficiency in terms of fine turbulence and its spatial homogeneity may be questioned. This study has emphasized, from the fluidodynamic point of view, the reasons why Rushton rotating paddles have a disappointing mixing performance: to be more specific, the formation of stagnate/slow recirculating zones and the highly localized, even though exorbitant, shear rates. The global mixing efficiency is hardly improved by the impeller’s geometry, as Fig. 3 proves. In spite of the removal of these limitations by the proposed new agitator (Fig. 5a shows the best configuration, where rotating paddles sweep away any stagnate zone), the numerical simulation showed that no combination of parameters is capable of avoiding the formation of a central slow motion zone close to the shaft and the creation of a peripheral turbulent spot (clearly visible in Figs. 8 and 9); although these spots are of serviceable turbulence, potentially useful for applications, still the spatial homogeneity of the -scale is not adequate; thus many revolutions are needed to achieve the safe target level of homogeneity. It is concluded that rotating impellers seem intrinsically not efficient since the rotation itself amplifies the limitations of the apparatus: so there are reasons to reconsider the actual convenience of rotating stirrers in terms of energy consumption. Although cavitation, buoyancy effects, chemical reactions and two-phase flows were not considered in the numerical simulation, it is felt that the present conclusions should be barely affected by these phenomena. According to [1], the number and size of the holes control turbulence intensities downstream of 2D grids; this study showed that the above 2D result cannot be extended to revolving perforated plates in cylindrical domains and this is consistent with the dissipation features of perforated plates illustrated by Malavasi et al. [2].
As the concluding remark, it seems appropriate to mention oscillating translating grids as an alternative technique in place of rotating grids. The translational vibrating grid is a well established technology in water treatment plants, particularly in coagulation and flocculation plants (see, among many others, Ref. [30]): in this case, the grid is plane and typically vibrating at small amplitudes and large frequencies. Liem et al. [30] demonstrated the great potential of vibrating grids to produce fine and homogeneous turbulence; their exhaustive review [30] also confirmed that rotating impellers produce significant spatial variations of turbulent quantities and that different types of impellers did not reveal any major difference among each other; this point has been furtherly corroborated in this study (Fig. 3). Acknowledgments The author is grateful to Dr. G. Bazzanella for performing Ansys simulations. The comments of the three anonymous Referees are gratefully acknowledged. The kind permission of Elsevier and of Professor Andrè Bakker to reproduce Figs. 2 and 4 is appreciated. Appendix A. Details of numerical simulations Pn = paddles’ number; Hn = holes’ number; dh = holes’ diameter; ω = angular velocity; ¯ = weighted Kolmogorov scale; G = Gini coefficient; = G · ¯ (Eq. (6)); = standard deviation; Sk = skewness; Ku = kurtosis; n/a = not applicable. Pn
Hn
dh [mm]
ω [s−1 ]
¯ [m]
2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
30 30 30 30 30 30 20 20 20 20 20 20 12 12 12 12 12 12
0.5 0.5 1.0 1.0 2.0 2.0 0.5 0.5 1.0 1.0 2.0 2.0 0.5 0.5 1.0 1.0 2.0 2.0
100 50 100 50 100 50 100 50 100 50 100 50 100 50 100 50 100 50
66 107 66 108 59 91 64 110 66 107 62 95 69 116 66 109 64 102
0.88 0.82 0.87 0.83 0.85 0.86 0.87 0.82 0.88 0.84 0.88 0.87 0.85 0.78 0.87 0.81 0.89 0.86
4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4
30 30 30 30 30 30 20 20 20 20 20 20 12 12 12 12 12 12
0.5 0.5 1.0 1.0 2.0 2.0 0.5 0.5 1.0 1.0 2.0 2.0 0.5 0.5 1.0 1.0 2.0 2.0
100 50 100 50 100 50 100 50 100 50 100 50 100 50 100 50 100 50
38 63 40 63 33 53 41 68 40 65 38 56 42 71 41 68 40 65
0.94 0.93 0.95 0.93 0.97 0.93 0.94 0.92 0.94 0.93 0.93 0.94 0.94 0.92 0.94 0.93 0.94 0.93
8 8
30 30
2.0 3.0
100 100
270 180
0.95 0.94
2 2 4 4
0 0 0 0
n/a n/a n/a n/a
100 50 100 50
74 122 48 80
0.82 0.74 0.95 0.92
G
[m]
[m]
Sk
58.1 87.7 57.4 89.6 50.1 78.3 55.7 90.2 58.1 89.9 54.6 82.7 58.7 90.5 57.4 88.3 56.7 87.7
20.5 43.8 23.8 40.8 23.8 33.8 21.9 43.2 22.3 39.2 23.8 35.0 24.7 45.0 22.4 39.0 22.2 36.7
−1.80 −0.59 −0.23 −0.54 0.27 −0.48 −0.28 −0.53 −0.30 −0.46 0.05 −0.48 0.05 −0.39 −0.25 −0.48 −0.29 −0.59
5.78 3.58 3.37 3.26 2.81 3.33 3.91 3.45 3.60 3.46 3.32 3.59 3.68 3.27 3.53 3.50 3.70 3.62
35.7 58.6 38.0 58.6 32.0 49.3 38.5 62.6 37.6 60.5 35.3 52.6 39.5 65.3 38.5 63.2 37.6 60.5
8.4 16.0 7.8 14.5 6.2 13.4 9.0 17.6 8.0 14.4 7.8 13.3 9.4 18.3 8.2 15.6 7.9 15.5
−1.11 −1.48 0.33 −0.68 0.17 −0.46 −0.54 −1.20 −0.10 −0.82 −0.18 −0.64 0.11 −0.76 0.29 −0.58 0.37 −0.74
8.54 6.69 12.03 6.25 12.67 6.94 8.79 5.75 8.83 6.50 8.40 6.85 9.33 5.09 9.27 5.72 13.14 6.29
48 34
−1.88 −1.23
14.23 7.26
23.5 40.5 9.8 17.9
0.32 −0.07 0.84 −0.34
3.32 3.13 10.17 4.85
257 169 60.7 90.3 45.6 73.6
Ku
F. Trivellato / Chemical Engineering and Processing 50 (2011) 799–809
References [1] R. Liu, D.S.K. Ting, Turbulent flow downstream of a perforated plate: sharpedged orifice vs. finite-thickness holes, J. Fluids Eng. 129 (2007) 1164–1171. [2] S. Malavasi, G. Messa, S. Macchi, The pressure drop coefficient through sharpedged perforated plates, in: 32◦ Hydraulics and Hydraulic Structures National Congress, Palermo, Italy, 2010. [3] V. Buwa, A. Dewan, A.F. Nassar, F. Durst, Fluid dynamics and mixing of singlephase flowing a stirred vessel with a grid disc impeller: experimental and numerical investigations, Chem. Eng. Sci. 61 (2006) 2815–2822. [4] T. Kumaresan, B.J. Jyeshtharaj, Effect of impeller design on the flow pattern and mixing in stirred tanks, Chem. Eng. J. 115 (2006) 173–193. [5] V.W. Uhl, J.B. Gray, Mixing: Theory and Practice, Orlando Academic Press, 1986. [6] B.S. Choi, B. Wan, S. Philyaw, K. Dhanasekharan, T.A. Ring, Residence time distributions in a stirred tank: comparison of CFD predictions with experiment, Ind. Eng. Chem. Res. 43 (20) (2004) 6548–6556. [7] A.S. Khare, K. Niranjan, An experimental investigation into the effect of impeller design on gas hold-up in a highly viscous Newtonian liquid, Chem. Eng. Sci. 54 (1999) 1093–1100. [8] L. Cortellezzi, A. Adrover, M. Giona, Feasibility, efficiency and transportability of short duration optimal mixing protocols, J. Fluid Mech. 597 (2008) 199–231. ´ Control of mixing in fluid flow: a maximum [9] D. D’Alessandro, M. Dahleh, I. Mezc, entropy approach, IEEE Trans. Automat. Contr. 44 (1999) 1852–1863. ´ G. Tadmor, A. Banaszuk, Optimal mixing in recirculation [10] B.R. Noack, I. Mezic, zones, Phys. Fluids 16 (2004) 867–888. [11] M.M. Alvarez-Hernandez, T. Shinbrot, J. Zalc, F.J. Muzzio, Practical chaotic mixing, Chem. Eng. Mixing 57 (2002) 3749–3753. [12] L. Oshinowo, A. Bakker, E. Marshall, Simulating Mixing Time with Computational Fluid Dynamics, vol. JA111, Fluent Inc. Press, 2000, pp. 1–5. [13] P.T.L. Koh, J.G.R. Andrews, P.T.H. Uhlherr, Flocculation in stirred tanks, Chem. Eng. Sci. 39 (1984) 975–985. [14] H. Wu, G.K. Patterson, Laser-Doppler measurements of turbulent flow parameters in a stirred mixer, Chem. Eng. Sci. 44 (1989) 2207–2221.
809
[15] L.A. Cutter, Flow and turbulence in a stirred tank, AIChE J. 12 (1966) 35–45. [16] S.M. Kresta, K.J. Bittorf, D.J. Wilson, Internal annular wall jets: radial flow in a stirred tank, AIChE J. 47 (2001) 2390–2401. [17] I.P.D. De Silva, Turbulence dissipation in stirred jars, J. Eng. Mech. ASCE 132 (11) (2006) 1260–1268. [18] R.M. Smith, L.X. Liu, J.D. Litster, Breakage of drop nucleated granules in a breakage only high shear mixer, Chem. Eng. Sci. (2010) (published online). [19] A. Arbel, A. Shklyar, Numerical simulations of turbulent flow through fine screen, in: MASCOT09, IMACS Workshop, IAC-CNR, Rome, Italy, 2009. [20] D.A. Deglon, C.J. Meyer, CFD modeling of stirred tanks: numerical considerations, Miner. Eng. 19 (2006) 1059–1068. [21] I.B. Celik, U. Ghia, P.J. Roache, C.J. Freitas, H. Coleman, P.E. Raad, Procedure for estimation and reporting of uncertainty due to discretization in CFD applications, J. Fluids Eng. ASME 130 (7) (2008) 1–4. [22] I. Celik, O. Karatekin, Numerical experiments on application of richardson extrapolation with nonuniform grids, J. Fluids Eng. ASME 119 (1997) 584–590. [23] J.H. Ferziger, M. Peric, Further discussion of numerical errors in CFD, Int. J. Numer. Methods Fluids 23 (1996) 1263–1274. [24] W. Bartok, C.E. Heath, M.A. Weiss, Mixing in a jet stirred reactor, AIChE J. 6 (4) (1960) 685–687. [25] J.A. Sanchez Perez, E.M. Rodriguez Porcel, J.L. Casas Lopez, J.M. Fernandez Sevilla, Y. Chisti, Shear rate in stirred tank and bubble column bioreactors, Chem. Eng. J. 124 (2006) 1–5. [26] G.B. Schubauer, W.G. Spangerberg, P.S. Klebanoff, Aerodynamic characteristics of damping screen, in: Tech. Note 2001, NACA, Washington, DC, 1950. [27] M.O. Lorenz, Methods of measuring the concentration of wealth, Q. Publ. Am. Stat. Assoc. 9 (70) (1905) 209–219. [28] C. Gini, Measurements of inequality and income, Econ. J. 31 (1921) 124–126. [29] J. Ding, X. Wang, X. Zhou, N. Ren, W. Guoa, CFD optimization of continuous stirred-tank (CSTR) reactor for biohydrogen production, Bioresour. Technol. 101 (18) (2010) 7005–7013. [30] L.E. Liem, D.W. Smith, S.J. Stanley, Turbulent velocity in flocculation by means of grids, J. Environ. Eng. ASCE 125 (3) (1999) 224–233.