Three-Dimensional Deformation Dynamics of Trailing Vortex

Dec 31, 2008 - dynamics resulting from vortical structures in stirred vessels are presented. ... stirred vessel, and their implications for enhancemen...
0 downloads 0 Views 5MB Size
Download PDF
8148

Ind. Eng. Chem. Res. 2009, 48, 8148–8158

Three-Dimensional Deformation Dynamics of Trailing Vortex Structures in a Stirred Vessel Yann Bouremel,* Michael Yianneskis, and Andrea Ducci Experimental and Computational Laboratory for the Analysis of Turbulence (ECLAT), DiVision of Engineering, King’s College London, Strand WC2R 2LS, U.K.

Novel three-dimensional methodologies suitable for the characterization and quantification of the deformation dynamics resulting from vortical structures in stirred vessels are presented. The analytical approaches are developed in the context of the three-dimensional deformation encountered by a spherical fluid element in the vicinity of the trailing vortices around a Rushton impeller and allow the identification of the rate and type of dynamics (stretching, compression, biaxial strain, or biaxial compression) encountered by such an element in the vortex flow field, as well as the resulting surface deformation at different parts of the vortex. The results facilitate a more in-depth understanding of the deformation processes involved in fluid mixing in a stirred vessel, and their implications for enhancement of mixing in such processes are discussed. 1. Introduction Mixing processes in stirred vessels are encountered in the majority of chemical industries, and their optimization is of great importance to improve energy efficiency, product yield, and selectivity. Substantial efforts have been made in recent years to improve understanding of the fluid mechanics of mixing at different spatial scales, both laminar and turbulent, but such understanding is still far from complete; in addition, process optimization often depends critically on process-specific data that are very useful, but generalization and/or extrapolation of such knowledge to different vessel designs and scales is hampered by the complex dynamics involved. Turbulent mixing has been addressed through either largescale studies, which approximate the flow patterns through the vessel, or statistical turbulence and structural turbulence studies that are based on averaged fluctuating quantities and focus more specifically on the turbulent kinetic energy dissipated by viscosity into internal energy; on the other hand, laminar mixing involves an increase of interfacial areas and diffusion.1 Shear and normal deformations create new surface areas that diffuse and mix. Experimental studies such as those of Ali et al.2 and Chan et al.3 on drop dispersion inside reactors stirred by pitchedblade and radial discharge impellers have shown that ligament stretching as well as turbulent fragmentation are responsible for the deformation/break of drops inside the trailing vortices. Experimental, numerical, and theoretical work on drop deformation and burst in shear flow have been reviewed by Rallison.4 Erwin5 considered how the flow surface deformation rate can be expressed in terms of the strain tensor principal values and principal directions, comparing the energy for different types of strain to deform the fluid. More recently, Ottino6 formulated equations for deformations in one, two, and three dimensions, the length, area, and volumetric deformations, respectively. An approach often employed to characterize mixing is that of the Baldyga and Bourne micromixing model,7 the engulfmentdeformation-diffusion model, where the deformation of structures embedded in a turbulent fluid is considered. Much effort has been expended on turbulent mixing, but the small scales involved in the dissipation of turbulent energy8,9 * To whom correspondence should be addressed. E-mail: [email protected]. Tel.: +44 (0)20 7848 2522. Fax: +44 (0)20 7848 2932.

and the nonuniform distribution of such energy in stirred vessels10 present considerable obstacles to its accurate determination. On the other hand, there is a plethora of data already available in the published literature on the velocity distributions in stirred tanks11-14 that can be readily harnessed to extract information on nonturbulent aspects of mixing and in particular the strain dynamics of vortical structures that are dissipative in nature. The present work aims to make use of such information and to develop appropriate methodologies in order to better comprehend how fluid elements are deformed in the vortical flowfields in stirred vessels. The deformation of a fluid element/ blob of fluid15 in the trailing vortices16 emanating from a Rushton impeller is employed as a paradigm, in the context of a spherical element subjected to stretching and compression near the impeller blades. Better knowledge of the strain and deformation dynamics should enable to locate regions of high stretching/ compression and/or deformation inside a stirred reactor where mixing might be performed more efficiently, as regions of high deformation create large interfacial areas that in turn facilitate diffusion of concentration fields. The methodologies were first developed on a simpler flowsa vortex ringsand in two dimensions.17,18 Clearly the flows concerned are fully three-dimensional, and in this paper the approaches are extended to encompass deformation of an element in all three directions so as to enable a more complete picture of the deformation processes, for the generic case where the element/blob concerned is of a size of the same order as the vortical structures it encounters. The work aims thus to provide a foundation on which more complex approaches, for example, including effects of interfacial tension, can be developed in future. The methodologies are developed by making use of relevant laser Doppler anemometry (LDA) data12,19 and building on appropriate analytical and theoretical considerations for simpler flows.6,20 In the following section the stirred vessel configuration in which the data were obtained and the principal measurement system characteristics are described. Subsequently, in sections 3 and 4, the deformation rate tensor and surface deformation rate approaches are developed and explained. These are then applied to the flow field of interest, and the results obtained in terms of the trailing vortices, the rates of stretching and compression, and the surface deformation and dissipation rate

10.1021/ie801481v CCC: $40.75  2009 American Chemical Society Published on Web 12/31/2008

Ind. Eng. Chem. Res., Vol. 48, No. 17, 2009

8149

are presented and discussed. The final section provides an outline of the main conclusions and observations from the work and discusses their implications for further development and optimization of mixing processes. 2. Flow Configuration and Experimental Apparatus The mean and rms velocity data presented here are reported in detail in ref 19. The laser Doppler anemometry (LDA) measurements presented therein carried out in the impeller stream of a stirred tank (T ) 150 mm) are utilized in this work to characterize the strain dynamics of the trailing vortices generated by a Rushton turbine. The tank design in which the measurements were made was standard with four equi-spaced baffles of width T/10 and thickness T/75. The impeller of diameter D ) T/3 was set at a height of C ) T/3 from the bottom of the tank. The thickness, width, and height of the impeller blades were T/75, T/12, and T/15, respectively. The working fluid was white oil FC 2012 W of viscosity ν ) 4 × 10-5 m2 s-1. A schematic diagram of a standard configuration tank is shown in Figure 1; the rig was equipped with a trough to minimize the refractive distortion at the tank cylindrical surface. A more extensive description of the vessel is given in refs 12 and 19. The LDA data analyzed in this paper was obtained at an impeller rotational speed of N ) 2670 rpm, corresponding to a Re ) ND2ν-1 ) 2780. The radial, tangential, and axial velocity components (ur, uθ, uz) were obtained in the mid plane between two consecutive baffles for 0.16 < r/D < 0.76 and 0.8 < z/D < 1.2 with a measurement spatial resolution of ∆r ) ∆z ) 0.02D. Full three-dimensional characteristics of the trailing vortex region were obtained by carrying a phase resolved analysis of the data. In the remainder of the paper a cylindrical coordinate system (r, θ, and z) with origin in the center of the base of the vessel will be used. The impeller rotated in a counterclockwise direction as viewed from above so that the tangential coordinate increases when moving away from the trailing blade. 3. Deformation Rate Tensor The rate of deformation of a fluid particle can be estimated from the gradient of the velocity vector (∇u) which can be decomposed into a symmetric part, S, denoted as the strain rate tensor, and an antisymmetric one, Ω, denoted as the rotation rate tensor.20 The definitions of the strain and rotation rate tensors in cylindrical coordinates are given in eqs 1 and 2:

(

S) ∂ur ∂r 1 ∂ uθ 1 ∂ur r + 2 ∂r r r ∂θ 1 ∂ur ∂uz + 2 ∂z ∂r

[ ()

1 ∂ uθ 1 ∂ur r + 2 ∂r r r ∂θ 1 ∂uθ ur + r ∂θ r 1 1 ∂uz ∂uθ + 2 r ∂θ ∂z

[ () ] [ ] [

(

]

] [ [

0 -ωz ωθ 1 ωz 0 -ωr Ω) 2 -ωθ ωr 0 where the vorticity ω is defined in eq 3: ω)

(

) (

) (

)

]

)

1 ∂ur ∂uz + 2 ∂z ∂r 1 1 ∂uz ∂uθ + 2 r ∂θ ∂z ∂uz ∂z (1)

]

(2)

)

∂ur ∂uz 1 ∂uz ∂uθ 1∂ 1 ∂ur i + i + (ru ) i r ∂θ ∂z r ∂z ∂r θ r ∂r θ r ∂θ z (3)

Figure 1. Schematic diagram of a standard configuration stirred vessel with external trough to minimize beam refraction.

A better understanding of the strain dynamics encountered in different flow regions can be obtained by estimating the principal components of the strain rate tensor as shown in eq 4:

(

)

(

-ω3* ω2* 0 S11* 0 0 1 -ω1* ω3* ∇u ) 0 S22* 0 + 0 2 0 0 S33* -ω2* ω1* 0

)

(4)

where S11*, S22* and S33* are the eigenvalues of S and ω1*, ω2*, and ω3* are the components of the vorticity along the local principal axes of the strain rate tensor. Each Sii* and the corresponding eigenvector defines the local intensity and direction of stretching and compression; if Sii* > 0 the fluid element considered is stretched, whereas if Sii* < 0 it is compressed. This reference system conversion allows the main local rates of stretching and/or compression to be quantified and also whether local vortical regions are affected by biaxial strain or biaxial compression to be identified by computing the parameter S11* × S22* × S33*. For an incompressible fluid (∇ · u) regions with S11* × S22* × S33* > 0 are subject to biaxial compression, while regions with S11* × S22* × S33* < 0 are subject to biaxial strain. When considering the principal reference system of the strain rate tensor the amount of kinetic energy dissipated by the viscous stresses to deform a fluid particle can be obtained from eq 5: ε ) 2ν(S11*2 + S22*2 + S33*2)

(5)

To assess the accuracy of the different derivative terms shown in eqs 1 and 2 that are necessary to estimate S and Ω, some preliminary data analysis was carried out to check whether the continuity equation (eq 6) was satisfied. 1 ∂uθ ur ∂ur ∂uz + + + )0 (6) r ∂θ r ∂r ∂z A comparison between the tangential profiles of (1/r)(∂uθ/ ∂θ)N-1 and -[(ur/r) + (∂ur/∂r) + (∂uz/∂z)]N-1 is shown in Figure 2a-c for three positions of coordinates (r/D ) 0.3, z/D ) 0.88), (r/D ) 0.46, z/D ) 1.14), and (r/D ) 0.42, z/D ) 1.08), respectively. The two types of profile show a generally good agreement for all the locations considered, with the standard deviation of the absolute difference between the two curves being approximately 10% of the variation of (1/r)(∂uθ/∂θ)N-1 over 0 < θ < 60°. Higher discrepancies between the two profiles are encountered close to the trailing and leading blades for low and high values of θ, respectively. The maximum difference is present near the blade at (r/D ) 0.46, z/D ) 1.08) for θ ≈ 56°.

8150

Ind. Eng. Chem. Res., Vol. 48, No. 17, 2009

Figure 2. Variation of (1/r)(∂uθ/∂θ)N-1 and -[(ur/r) + (∂ur/∂r) + (∂uz/ ∂z)]N-1 with the tangential coordinate θ; (a) r/D ) 0.3, z/D ) 0.88; (b) r/D ) 0.46, z/D ) 1.08; (c) r/D ) 0.42, z/D ) 1.14.

As a consequence the range of θ investigated was narrowed to 6° e θ e 51°.

Figure 3. Sketch of the deformation of a material surface with n ) i3 in a time interval dt.

4. Surface Deformation Rate The deformation rate of a material surface was used by Erwin5 to classify different types of flow associated to specific deformation characteristics (simple shear, pure shear, and pure strain) in terms of mixing effectiveness. More recently Ottino6 derived an expression for the specific rate of stretching, η, of a material surface, A, shown in eq 7: D ln η ) ∇ · u - S : nn Dt

(7)

where nn denotes the dyadic tensor obtained from a unit vector n perpendicular to the material surface considered. For incompressible flow (∇ · u ) 0) eq 7 can be rewritten as D ln η ) -S : nn ) -S(n) · n ) -Snn Dt

(8)

where Snn is the component along n of the deformation vector S(n) applied on a surface with normal unit vector n. Equation 8 can be readily derived by considering the sketch shown in Figure 3, where a material surface of normal i3 and initial size ∆x1 × ∆x2 is stretched along both the x1 and the x2 directions during an infinitesimal time interval dt. As a first approximation (neglecting terms of the order of dt2), the specific variation of the surface area dA/A can be obtained by adding the two shaded areas shown in Figure 3. The proportionality between dA/A and the component of the deformation tensor normal to the surface considered is given by eq 9: d ln(η) )

(

)

∂u1 ∂u2 dA dA ) ) + dt ) A ∆x1∆x2 ∂x1 ∂x2 ∂u3 dt ) -S33 dt (9) ∂x3

As previously mentioned the present study aims at providing an estimation of the surface deformation rate that would be experienced by a virtual spherical particle positioned in proximity of the trailing vortices. The spherical shape was selected for its perfect symmetry as in this way the final estimate of the rate of surface deformation is independent of the orientation of the particle in space.

Figure 4. Sketch of the deformation rate Snn locally normal to the external surface of a sphere.

The rate of surface deformation was estimated for virtual spherical particles of radius rs centered in each point of measurement. The local deformation tensor present on the surface of the sphere was three-dimensionally interpolated from the values of the deformation tensors obtained in the points of measurement. It should be noted that before carrying the 3D linear interpolation, the local deformation tensor had to be converted from a cylindrical to a Cartesian reference system. The surface of the sphere was meshed according to the spherical coordinate system shown in Figure 4, where R and β are the polar and azimuthal angles, respectively. The rate of surface deformation, D ln η/Dt, of a generic facet positioned on the external surface of the sphere was obtained from eq 8, with n ) [sin(R) cos(β), sin(R) sin(β), cos(R)]. The total surface deformation rate, η′tot, was estimated from the integral of Snn over the external area of the sphere (eq 10).

∫∫

1 π 2π S sin(R) dβ dR (10) 4π 0 0 nn It should be noted that the estimate of η′tot is dependent on the values of ∆R and ∆β used to mesh the external surface of the spherical particle. For example, the contour plots of Snn(R, β) obtained for angular mesh sizes of 3° and 25° are shown in η′tot ) -

Ind. Eng. Chem. Res., Vol. 48, No. 17, 2009

8151

Figure 5. Contour plots of Snn/N for spheres of different azimuthal and polar mesh size; (a) ∆R ) ∆β ) 3°; (b) ∆R ) ∆β ) 25°.

Figure 6. Variation of η′tot and As/(4πrs2) for different sphere mesh sizes, ∆R, ∆β; circles, η′tot at (r/D ) 0.38, z/D ) 1.03, θ ) 10°); squares, η′tot at (r/D ) 0.4, z/D ) 1.03, θ ) 10°); diamonds, As/(4πrs2).

Figure 5a,b, respectively. It is evident that finer meshes provide a continuous variation of Snn along the external surface of the sphere with a final shape that deviates significantly less to the one of a perfect sphere. Moreover it has to be stressed that smaller values of ∆R and ∆β allow a wider range of facet orientations to be considered, which are more likely to be aligned with the vector joining the sphere origin and the points of the measurement grid closest to the origin. This enables the interpolation procedure to be optimized and in Snn(R, β) the characteritics of the deformation tensors which were directly estimated in the measurement points best captured. The variation of η′tot/N with angular mesh size, ∆R and ∆β, is shown on the left-hand ordinate (LHO) of Figure 6, where circle and square symbols are associated to two adjacent points of coordinates (r/D ) 0.38, z/D ) 1.03, θ ) 10°) and (r/D ) 0.4, z/D ) 1.03, θ ) 10°), respectively. For values of ∆R, ∆β < 6°, η′tot showed a variation within 10% of the estimate obtained for the lowest mesh size (i.e., η′tot(3°)). Coarser meshes were more strongly affected by the orientation of the sphere facets with respect to the points of the measurement mesh. In particular it was found that η′tot varied by more than 40% for ∆R, ∆β > 20° for both the points considered. The variation of the sphere surface As with respect to the mesh size is also shown on the right-hand ordinate (RHO) of Figure 6. The approximate estimate of the sphere surface As obtained from the integral over the finite mesh is normalized with the real external surface of a sphere of radius rs. As expected finer meshes yield a final shape closer to a perfect sphere as As/(4πrs2) approaches 1 for decreasing ∆R and ∆β (and is around 0.99 for ∆R, ∆β < 6°). To minimize the

Figure 7. Schematic diagram of the upper and lower trailing vortices by means of iso-vorticity contours behind a blade.

variability of η′tot due to the sphere mesh size the data presented in the rest of the paper were obtained for ∆R, ∆β < 5°. 5. Results and Discussion In the following subsections the presentation and discussion of the data has been organized in four parts, that is, trailing vortex vorticity, rates of strain, surface deformation, and dissipation rate, so as to gain a complete three-dimensional understanding of the vortex and deformation dynamics occurring in the trailing vortex region. In brief, the rationale of the work was to identify and quantify regions of the trailing vortices where the linear and surface deformations reach maximum and minimum values and to assess the efficiency of the system by comparing the local surface deformation and the dissipation rate. 5.1. Trailing Vortices. A three-dimensional view of the trailing vortices generated behind a RT blade is shown in Figure 7, where the shape of the vortices was visualized by estimating iso-vorticity surfaces; the bottom vortex, in red, is associated to ωθ/N ) 27, while the upper one, in blue, to ωθ/N ) -27. In agreement with the findings of refs 11 and 13 the two counterrotating vortices are separated by the impeller disk and their surface is curved backward and radially outward when moving away from the leading blade. It should be noted that the trajectory of the trailing vortices may be affected by small differences in the impeller/vessel configuration. The contour plots of the phase-resolved component of the tangential vorticity, ωθ/N, are shown in Figure 8a,b for axial positions above and below the impeller centerline, respectively. To improve the legibility of the figures and at the same time provide sufficient information about the quantity considered, both plots have been subdivided in six sectors, each representing

8152

Ind. Eng. Chem. Res., Vol. 48, No. 17, 2009

Figure 8. Contour plots of the tangential vorticity component, ωθ, for six different axial positions above (a) and below (b) the impeller centerline.

the phase resolved term (over an angle θ of 60°) for six axial positions (i.e., z/D ) 0.99, 1.01, 1.03, 1.05, 1.07, and 1.09 in Figure 8a and z/D ) 0.89, 0.91, 0.93, 0.95, 0.97, and 0.99 in Figure 8b), with the axial coordinate increasing counterclockwise. From Figure 8a,b it is evident that the upper and lower trailing vortices have similar strength as both of them show an absolute maximum intensity of 60N at near opposite locations with respect to the impeller centerline (i.e., z/D ) 1.05 and z/D ) 0.97). The slight asymmetry with respect to the plane of the disk is due to the slight curvature of the impeller stream and the vortex axis to the horizontal (5-7°, as reported in ref 11). The axis of each vortex was estimated for three selected axial positions where the vorticity reaches higher intensity values, namely, z/D ) 1.03-1.07 and z/D ) 0.93-0.97 for the upper and lower trailing vortices, respectively. The radial coordinate rtv(θ) of the axis was estimated from the vorticity weighted average shown in eq 11:

Figure 9. Contour plots of the principal components of the strain rate tensor, Sii*, for six axial positions above the impeller centerline: (a) S11*; (b) S22*; and (c) S33*.

Ntot

rtv(θ) )

∑ω i)1 Ntot

∑

θi(θ)ri

(11)

ω (θ) i)1 θi

where Ntot is the total number of points used to carry the average. The loci of the axis radial coordinates rtv(θ) of the upper and lower trailing vortices are shown as white lines in Figure 8a,b for all the z/D positions selected. The presence of inner trailing vortices12,14,21 generated at the inner edge of the blade is also clear in Figure 8a,b for axial positions z/D ) 1.03 and 0.93, right above and below the impeller disk, respectively. It should be noted that the inner trailing vortices are counter-rotating with respect to the corresponding external ones. 5.2. Rates of Stretching and Compression. The contour plots of the principal components of the strain rate tensor, S11*, S22*, and S33* are shown in Figure 9a-c, respectively, for six

Figure 10. Contour plots of the triple product of the principal strain rates, S11* × S22* × S33*, for six axial positions above the impeller centerline.

heights above the impeller centerline. The strain rates S11* and S33* are always associated to compression and stretching, respectively, while S22* changes sign, from positive to negative values (stretching to compression), when positions of increasing

Ind. Eng. Chem. Res., Vol. 48, No. 17, 2009

8153

Figure 11. Comparison of the axial profiles of the principal components of the strain rate, Sii*, and of their triple product, S11* × S22* × S33*, with the axial profile of the tangential vorticity component, ωθ, at r/D ) 0.44 and θ ) 40°; (a) S11* (LHO), ωθ (RHO); (b) S22* (LHO), ωθ (RHO); (c) S33* (LHO), ωθ (RHO); (d) S11* × S22* × S33* (LHO), ωθ (RHO).

Figure 12. Correlation between the angle γ1* and the tangential vorticity, ωθ, for all the points of measurement located in the range z/D ) 1.01-1.07.

Figure 13. Joint probability density function of the angles γ1* and δ1* for points with ωθ/N < -25 and z/D ) 1.01-1.07.

z/D are considered. From Figure 9a,b, it is evident that S11* and S33* are the dominant strain rates as they reach peak absolute intensities of 50N, while S22* varies within a narrower range of -12 e |S22*/N| e 12 across the different heights considered (see Figure 9b). As expected all the three strain rates assume local absolute maxima at the border of the vortex core, closer to the inner and outer edges of the leading blade, for each z/D considered. For heights z/D ) 1.03-1.05, where the strength of the upper trailing vortex is highest (see Figure 8a), the strain rate S22* is closer to zero, and as a consequence the deformation

Figure 14. Probability density function of γi* for points with ωθ/N < -25 and z/D ) 1.01-1.07; (a) γ1*; (b) γ3*.

Figure 15. Contour plots of the rate of surface deformation, η′tot/N, for six axial positions above the impeller centerline.

of the local flow is almost bidimensional and mainly controlled by two strain rates of stretching and compression mutually orthogonal to each other. On the contrary, for axial positions z/D ) 0.99-1.01 closer to the impeller centerline, the local deformation characteristics are more complex and governed by all the strain rates. The behavior described is well reflected in Figure 10, where the contour plot of the strain rate triple product is shown for the six axial heights investigated. S11* × S22* × S33* assumes low negative values (≈ -3 × 104N3) for axial positions closer to the impeller centerline (z/D ) 0.99-1.01), confirming that

8154

Ind. Eng. Chem. Res., Vol. 48, No. 17, 2009

Figure 16. Profiles of η′tot and ωθ against the radial distance from the vortex axis, ∆R/D, for selected axial positions below and above the impeller centerline; (a) z/D ) 0.93; (b) z/D ) 0.95; (c) z/D ) 1.03; (d) z/D ) 1.05.

the local deformation rate is indeed three-dimensional and subject to biaxial strain. Opposite deformation characteristics are present at the upper border of the trailing vortex (z/D ) 1.07-1.09). For this range of z/D, S11* × S22* × S33* assumes positive values (≈0.5 × 104N3) indicating that the local flow is subject to biaxial compression. However, it has to be stressed that for this range of z/D the absolute value of S11* × S22* × S33* is at least six times smaller than that one shown at locations closer to the impeller centerline. Analogous graphs representing the contour plots of S11*, S22*, and S33* for the bottom vortex are not shown for brevity of presentation. However, their near-symmetry with respect to the impeller centerline (the reason for the small flow asymmetry has been mentioned earlier) can be clearly seen in Figure 11a-d, where the axial profiles of Sii* and of S11* × S22* × S33* are compared with the axial profile of the tangential vorticity at r/D ) 0.44 and θ ) 40°. As previously pointed out, the triple product shows a sharp minimum at the impeller centerline and two maxima of lower intensity at z/D ) 0.93 and z/D ) 1.09. From Figure 11c,d it can be seen that S22* and S33* intensities drop for z/D values corresponding to the maximum and minimum of the vorticity, while they assume higher absolute values at the border of the vortices. This behavior is partially present in Figure 11b, where S11* shows a major deviation from zero only at the impeller centerline, in between the two counterrotating vortices. Besides providing an estimate of the intensity of the strain rates, an effort was made to determine the directions of stretching and compression associated to S11*, S22*, and S33* in the core of the trailing vortices. It should be noted that the direction of the eigenvectors ii*, related to the strain rate Sii*, varies for each point of measurement, and therefore a statistical approach was selected to determine which directions are more likely to be associated with each strain rate. In the rest of the

section, γi* and δi* are defined as the angles included between the ith principal eigenvector and the tangential and radial directions, respectively. The correlation between the angle γ1* and the tangential vorticity is shown in Figure 12 for all the points of measurement located at z/D ) 1.01-1.07. It should be noted that the radial distance, ∆R/D, of each point from the vortex axis (the vortex axis shown as white line in Figure 8a) corresponds to symbols of hues varying from blue to red. Despite the large scatter, it is evident that the points (γ1*, ωθ) are organized in a Gaussianlike distribution, with locations of highest vorticity being associated to 70° e γ1* e 100°, while lower values of γ1* correspond to points of small ωθ, at the margin of the vortex core (i.e., circles of blue and red hues). It should be stressed that for values of γ1* close to 90° the local eigenvector ith is perpendicular to the θ direction and located in the plane r-z. To identify the direction of i1* in the core of the upper trailing vortex only the points with vorticity ωθ/N < -25 were analyzed further. The joint probability density function of the angles γ1* and δ1* for locations at z/D ) 1.01-1.07 with ωθ/N < -25 is shown in Figure 13. It is evident that for points in the vortex core, i1* is more likely to be oriented in space at angles γ1* ) 80°-90° and δ1* ≈ 60° from the tangential and radial directions, respectively. The orientations of the other two eigenvectors can be determined by considering that they are placed on a plane perpendicular to i1* and are mutually orthogonal to each other. It has to also to be stressed that the plane i2*-i3* contains the θ direction as γ1* ) 80°-90°. As a consequence the orientation of i2* and i3* can be determined by estimating the angles γ2* and γ3* between the eigenvectors and the θ direction. The probability distribution of the γ2* and γ3* for locations at z/D ) 1.01-1.07 with ωθ/N < -25 are shown in Figure 14a,b, respectively. Also in this case it is evident that i2* and i3* are

Ind. Eng. Chem. Res., Vol. 48, No. 17, 2009

8155

Figure 17. Correlation of η′tot/N with S11* × S22* × S33* for six positions above the impeller centerline; (a) z/D ) 0.99; (b) z/D ) 1.01; (c) z/D ) 1.03; (d) z/D ) 1.05; (e) z/D ) 1.07; (f) z/D ) 1.09.

more likely to be oriented with angles of γ2* ≈ 50° and γ3* ≈ -40° from θ. As expected i2* and i3* are orthogonal to each other. 5.3. Surface Deformation and Dissipation Rate. As outlined in section 4 the current results are intended to provide an estimate of the rate of surface deformation that would be experienced by a spherical particle placed in each point of the interrogation volume. In particular the results presented hereafter are related to spherical particles of diameter Ds ) 0.08D which is comparable to the size of the two trailing vortices (Dtv ) 0.1D).11 The contour plot of the total rate of surface deformation, η′tot/N, is shown in Figure 15 for six axial positions above the impeller centerline. The rate of surface deformation assumes mainly negative values (-4 e η′tot/N e -2), associated to a reduction of the sphere external surface, above and below the upper trailing vortex for axial positions z/D ) 0.99-1.01 and z/D ) 1.09, respectively. On the contrary external surface

expansion of the sphere is present for axial positions z/D ) 1.03-1.05, where η′tot is positive and shows local peaks in a range of 1 e η′tot/N e 3. Similarly to the behavior shown by the rates of strain, Sii* (see Figure 9), η′tot/N local absolute peak values are located at the edge of the upper trailing vortex for each axial position considered. However, it has to be remarked that the intensity of the surface deformation is considerably lower (|η′tot/Sii*| < 0.1) than those of the strain rates. A better view of the variation of η′tot/N can be obtained from Figure 16a-e, where the profiles of the surface deformation are compared with those of the tangential component of the vorticity at selected axial positions both above and below the impeller centerline, where η′tot/N shows local peak values. Similarly to the plot of Figure 12, the variables are plotted against the radial distance ∆R/D from the vortex axis. It should be noted that each plot of Figure 16 includes points with 16° e θ e 50° and 0.38 e r/D e 0.7. Direct comparison of Figure

8156

Ind. Eng. Chem. Res., Vol. 48, No. 17, 2009

Figure 18. Contour plots of the principal components of the strain rate tensor, Sii*, for six axial positions above the impeller centerline; (a) S11*; (b) S22*; (c) S33*.

16a,b with Figure 16c,d shows that η′tot/N has a similar behavior for axial positions below (z/D ) 0.93-0.95) and above the centerline (z/D ) 1.03-1.05), with a double peak shape at the edges of the vortex and a minimum in the vortex core. Despite the maxima of η′tot having also similar intensity of ≈2 between analogous locations in the upper and lower trailing vortices, it is evident that in Figure 16a,b the second peak for ∆R/D > 0 is less pronounced than the corresponding one in Figure 16c,d for positions above the impeller centerline. The correlation between the surface deformation and the type of strain dynamics (biaxial compression/strain) locally experienced by the flow is shown in Figure 17a-f, where η′tot is plotted against S11* × S22* × S33* for z/D ) 0.99-1.09. Similarly to Figure 16, each plot of Figure 17 includes points with 16° e θ e 50° and 0.38 e r/D e 0.7. As already mentioned in section 5.2, the strain rate tensor is strongly threedimensional close to the impeller centerline as S11 × S22 × S33 < 0 at z/D ) 0.99 (Figure 16a), while it becomes twodimensional (S11 × S22 × S33 < 0) for z/D in the core of the upper trailing vortex (Figure 16c,d). However, it is worth pointing out that for all the axial positions considered, points with positive η′tot/N are generally correlated to S11 × S22 × S33 ≈ 0, indicating that external surface expansion might be more likely to occur in flow regions where the local deformation strain rate is two-dimensional. The contour plot of the kinetic energy dissipated by the trailing vortex structures to deform a fluid particle is shown in Figure 18a for six axial positions above the impeller centerline. Similarly to the principal strain rates, ε is higher at the edges of the trailing vortex, close to the inner and outer tips of the blade. Higher values of ε occur at axial positions in the bottom half of the trailing vortex (z/D ) 0.99-1.03), with local maximum intensities within 7000 < ε/(νN2) < 10 000 (i.e., 2.5 < ε/N3D2 < 3.6). The present data confirm that the amount of

energy dissipated by the trailing vortex structures is significantly lower than the one dissipated by the smallest dissipative scales in the impeller stream. For example, Ducci and Yianneskis9 measured a maximum ensemble-averaged value of the dissipation rate of turbulent kinetic energy (12N3D2) which is around 4 times higher than the average found here from phase-resolved data. The difference is certainly under-estimated as from a phase-resolved analysis of a single term comprising ε; they suggested that phase-resolved turbulence dissipation rate should be at least 3 times higher than the corresponding one obtained from ensemble-averages (i.e., 36N3D2). Similar observations were also made by Hucher et al.22 who estimated instantaneous variations of the dissipation rate of the turbulence kinetic energy up to 20 times the corresponding ensemble-averaged value. An indication of the deformation efficiency of the flow can be obtained by estimating the ratio η′tot/(εν-1)0.5 which is indicative of the amount of energy locally dissipated by the viscous forces to extend the sphere external surface.6 It should be noted that η′tot/(εν-1)0.5 can also be considered as the ratio between the characteristic time scales of the viscous and surface deformations. The contour plot of the flow efficiency is shown in Figure 18b for six axial positions above the impeller centerline. Similarly to the behavior outlined for η′tot in Figure 15, positive (≈ 0.05) and negative (≈ -0.05) efficiencies are found in the core (z/D ) 1.03-1.05) and at the margin (z/D ) 0.99) of the upper trailing vortex, respectively. From this estimate it is evident that only a small amount (≈5%) of the local energy dissipated by the trailing vortices is actually employed to expand or shrink the external surface of a virtual particle. The variation of the rate of surface deformation, η′tot/N, with the size of the virtual particle is shown in Figure 19a,b for Ds/D ) 0.04-0.1. The radial profiles of η′tot(Ds) for different Ds/D at z/D ) 1.05 and θ ) 50° are shown in Figure 19a, where the radial profile of the tangential vorticity is also plotted to provide a reference of the vortex core location. From Figure 19 it is evident that spherical particles of increasing diameter are subject to higher rates of surface deformation with the maximum of η′tot for Ds/D ) 0.1 being almost 4 times greater than the corresponding maximum obtained for Ds/D ) 0.04. This behavior might be expected as particles of size closer to the trailing vortex diameter are more affected by the rates of strain and compression occurring at this length scale. It should also be noted that as the ratio Ds/D is increased, the radial distance between the two local maxima of each profile decreases. This can be explained by considering that for a given particle position, larger particles are affected by rates of strain over longer distances, and for a location of the particle close to the vortex core, larger particles are simultaneously affected by both peaks of the rate of strain present at the edges of the vortex core. The general increase of the intensity of the surface deformation rate with increasing particle diameter is also confirmed in Figure 19b, where the probability distributions of the absolute value assumed by η′tot/N at z/D ) 1.05 (15° < θ < 50°; 0.275 < r/D < 0.725) are plotted for Ds/D ) 0.04-0.1. It is evident that pdf associated to larger Ds/D is broader with a larger number of points positioned at z/D ) 1.05 experiencing higher values of η′tot. 6. Conclusions The deformation dynamics in the region of the trailing vortices in a vessel stirred by a Rushton impeller have been studied in an effort to develop three-dimensional approaches

Ind. Eng. Chem. Res., Vol. 48, No. 17, 2009

8157

Figure 19. (a) Radial profiles of η′tot/N and ωθ/N for four different diameters of the virtual spherical particle (z/D ) 1.05; θ ) 50°); (b) probability distribution of η′tot/N for four different Ds/D at z/D ) 1.05 (15° < θ < 50°; 0.275