Development of Fluid Particle Breakup and ... - ACS Publications

Subscriber access provided by GAZI UNIV

Article

Development of fluid particle breakup and coalescence closure models for the complete energy spectrum of isotropic turbulence Jannike Solsvik, and Hugo A. Jakobsen Ind. Eng. Chem. Res., Just Accepted Manuscript • DOI: 10.1021/acs.iecr.5b04077 • Publication Date (Web): 17 Dec 2015 Downloaded from http://pubs.acs.org on January 20, 2016

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Industrial & Engineering Chemistry Research is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Development of fluid particle breakup and coalescence closure models for the complete energy spectrum of isotropic turbulence Jannike Solsvik∗ and Hugo A. Jakobsen Department of Chemical Engineering, Norwegian University of Science and Technology (NTNU), Trondheim, Norway ∗ E-mail:

[email protected]

Abstract The constitutive equations proposed in the literature for describing the fluid particle (i.e., bubble and drops) breakage and coalescence phenomena created by turbulence mechanisms are generally limited to the inertial range of scales and infinite Reynolds numbers. A consistent approach for extending these breakage and coalescence kernels to the complete energy spectrum of isotropic turbulence (i.e., the energy containing, inertial and dissipation subranges) and for a larger range of integral scale Reynolds numbers is proposed in this study. The model energy spectrum for the complete range of scales of isotropic turbulence proposed by Pope 1 is employed in this work. A corresponding integral relation for the secondorder longitudinal structure function model can be deduced from the model energy spectrum through the use of a Fourier transform pair via the scalar velocity correlation function. However, this integral relation is cumbersome to solve numerically thus an approximate algebraic model is proposed in this work for computing the second-order longitudinal structure function to reduce the computational costs for solving the complete population balance equation.

1 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The novel modeling procedure is illustrated extending a few of the more popular drop and bubble breakage frequency and collision rate kernel (coalescence) closures to the complete range of scales in the spectrum of isotropic turbulence. The predictions of the extended closure laws are compared to those of the original model versions of the number density of vortices, breakage frequency, daughter size distribution and collision rate which are limited to the inertial range of isotropic turbulence at infinite integral scale Reynolds numbers.

Introduction The fluid particle breakup and coalescence phenomena taking place in turbulent dispersions are generally modeled employing the Kolmogorov turbulence energy spectrum model and the Kolmogorov second-order longitudinal structure function, which are valid for the inertial subrange of turbulence exclusively. The Kolmogorov inertial range models are known to be limited to fully developed turbulent flows in the infinite Reynolds number limit. However, dispersed multiphase flows in industrial chemical reactors, combustors and separators are generally characterized by finite integral scale Reynolds numbers thus the impact of Reynolds numbers might become very important and ought to be considered in the closure models. Moreover, the requirement that the breaking mother particle size corresponds to a vortex size within the inertial subrange of turbulence is not always fulfilled. Furthermore, for finite Reynolds number flows, the width of the inertial subrange of turbulence in which the Kolmogorov models are valid becomes more narrow and the inertial subrange might even disappear. Therefore, to develop more reliable and accurate models for drop and bubble breakup and collision rate the three-dimensional turbulent energy spectrum model and a consistent second-order longitudinal structure function for the complete range of scales are needed. After sufficiently accurate and consistent relations for these functions are established, the existing breakage and collision rate models limited to the inertial subrange of turbulence can be extended to the complete spectrum of turbulence. In a couple of recent articles, the energy spectrum model of Pope 1 has been adopted but with parameter values limited to infinite Reynolds numbers, and the structure function models used were not consistent with the model spectrum and 2 ACS Paragon Plus Environment

Page 2 of 38

Page 3 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


thus not sufficiently accurate. That is, either the second-order longitudinal structure function for the inertial range of turbulence was erroneously used for the complete spectrum 2 or a rough and inconsistent approximation for the entire spectrum was applied instead 3,4 . Moreover, it is noted that the correct effect of Reynolds number variations were disregarded because fixed parameter values for cL and cη were employed in the model energy spectrum thus limiting the spectrum to infinite Reynolds numbers. In the present work an exact relation derived from statistical turbulence theory is used to transform the three-dimensional energy spectrum into the corresponding second-order longitudinal structure function. The accuracy limitations of this modeling procedure is thus restricted to the uncertainties induced by the choice of model energy spectrum employed for describing the dissipative, inertial and energy containing subranges of turbulence. Moreover, it is important to consider the impact of integral scale Reynolds number variations thus the model spectrum parameter values cannot be kept constant because such an approach limits the spectrum to infinite Reynolds numbers. An analytical solution was obtained for the second-order longitudinal structure function considering both the energy-containing and inertial subranges of turbulence. However, for the dissipation subrange no analytical solution exists. Therefore, for the inertial to dissipation subranges of scales a semi-empirical interpolating relation has been adopted from the literature 5 and verified by numerical solutions. Based on these relations, a new approximate formula representing the second-order longitudinal structure function within the complete energy spectrum is proposed. By use of the model energy spectrum of Pope 1 and the new structure function model proposed in this work, a novel relation for the vortex number density for the complete spectrum have been developed. The new number density model predictions are compared to the results obtained by the existing model formulations limited to the inertial subrange. Moreover, the impact of extending the existing breakage and collision rate models limited to the inertial subrange of turbulence, by use of a three-dimensional energy spectrum and a consistent second-order longitudinal structure function valid for the complete range of turbulence scales have been elucidated. To validate the extended models the local experimental data of the breakage and coalescence



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

characteristics must be obtained under known flow conditions which can be associated to particular integral scale Reynolds numbers. Proper model validation thus requires both characterizations of the flow determining the local turbulent energy dissipation rate, turbulent energy, integral length scale in order to compute the local integral scale Reynolds number, and in parallel measuring the single fluid particle breakage events and investigating binary coalescence events. Independent data has been used to provide close to universal model parameter values for the turbulence model provided by Pope 1 and thus for the derivative structure functions, but the other empirical model parameters in the population balance kernels in which all disregarded and eventually unknown dependences are merged must be refitted to experimental data when the turbulence model is exchanged from the Kolmogorov models used in previous work to the extended Pope 1 spectrum model used in this work. It is noted that the other additional empirical parameters in the population balance closures which depend on physiochemical characteristics of the system (e.g., surfactants, mass transfer, dispersed phase physical data and transport coefficients) and possibly random events. Such detailed data is not available yet, even though several research groups are providing local experimental data on the fluid particle size characteristics, but normally without a suitable characterization of the local flow properties. For example, many of the previous works are performed in stirred tanks, nozzles, etc in which the turbulence properties are far from uniform and only average flow characteristics are provided. Future work must thus focus on constructing and building a suitable apparatus for performing the necessary experiments in order to gather sufficient experimental data. The large eddy simulation (LES) method can also be used to validate the statistical turbulence theory closures. In particular, the LES method has been used to validate the vortex number density model and the second order structure function 6 . Identification of individual vortices is thus required in order to verify the number density. Provided the fluid particles are large enough to be resolved by this model, the vortex fluid particle interactions can be investigated as well as the model assumptions of the absolute turbulent particle velocity, the relative particle velocity and the relative particle - vortex velocity being equal to the second-order longitudinal structure function.


Page 4 of 38

Page 5 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


However, if the LES resolution is not fine enough the fluid particles are not resolved. It might thus be better to apply direct numerical simulations (DNS) for model validation instead, resolving all scales of turbulence. The local turbulent energy dissipation rate per unit mass might be more accurately predicted by the DNS model than by LES. This issue is considered to be outside the scope of this article. The article contents is outlined as follows. In the first section, the necessary relationships from statistical turbulence theory are presented. The model energy spectrum of Pope 1 is outlined and the integral scale Reynolds number variation dependency computations are explained. Moreover, the exact statistical turbulence theory relationship allowing the second-order longitudinal structure function model to be determined from the three-dimensional energy spectrum model by quadrature computations is presented. Thereafter, the number density model for turbulent vortices are extended from the inertial subrange to the complete range of scales. The modeling of the breakage and coalescence (collision rate) phenomena created by turbulence is the last topic considered. A discussion of the statistical turbulence theory and the results of the present work are given, and finally conclusions are drawn.

Turbulence theory The energy spectrum model of Kolmogorov For the inertial subrange of scales, Kolmogorov 7 ,8 ,9 ,10 proposed a simple energy spectrum model: E(κ ) = Cε 2/3 κ −5/3

(1)

where ε denotes the turbulent kinetic energy dissipation rate per unit mass and κ represents the wave number. The value of the Kolmogorov energy spectrum parameter is normally set to be 11,12 C ≈ 1.5. It is important to note that this parameter value is determined in the infinite Taylor-scale Reynolds



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 38

number (Reλ ) limit. It is further expected that for finite Reynolds numbers, the parameter value varies somewhat with Reλ 11,12 .

The energy spectrum model of Pope Pope 1 proposed an energy spectrum model describing the energy containing, inertial and dissipating subranges of isotropic turbulence defined by: E(κ ) = Cε 2/3 κ −5/3 fL (κ L) fη (ηκ )

(2)

where fL (κ L) and fη (ηκ ) are semi-empirical non-dimensional functions. The semi-empirical function fL (κ L) determines the shape of the energy spectrum in the energy containing subrange and approaches unity for large κ L. L denotes the integral length scale and is defined by

L=

k3/2 ε

(3)

where k represents the turbulent kinetic energy. Moreover, fη (ηκ ) governs the shape of energy spectrum within the dissipation subrange and approaches unity for small κη , where η is the Kolmogorov length scale. In the inertial subrange, both fL (κ L) and fη (ηκ ) approach unity and thus the Kolmogorov −5/3 spectrum with constant C is recovered. The function fL (κ L) is given by { fL (κ L) =

κL [(κ L)2 + cL ]1/2

} 5 +p0 3

(4)

in which p0 is specified with a value of 2 and cL is a positive adjustable parameter. The function fη (ηκ ) is given by

{ } fη (κη ) = exp −β ([(κη )4 + c4η ]1/4 − cη )

where β = 5.2 and cη is a positive adjustable parameter.


(5)

Page 7 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


When the values of k, ε and ν are fixed, the model spectrum is specified by Equations (2)–(5) (in which three model parameters are determined by C = 1.5, β = 5.2 and p0 = 2). The values of the two adjustable parameters cL and cη are then fitted so that E(κ ) and 2νκ 2 E(κ ) integrate to k and ε , respectively:

∫ ∞

E(κ ) dκ

(6)

2νκ 2 E(κ ) dκ

(7)

k=

ε=

∫ ∞ 0

0

The second-order longitudinal structure function An exact relationship between the three-dimensional energy spectrum and the second-order longitudinal structure function has been established within the field of statistical turbulence theory 13 : 4 ⟨[δ v] ⟩1D,p (r) = 3 2

∫ ∞ 0

[ { }] sin(κ r) cos(κ r) E(κ ) 1 − 3 + dκ (κ r)3 (κ r)2

(8)

where E(κ ) denotes the three-dimensional energy spectrum, and r is the distance between the two measuring points. The subscript p stand for "parallel". For the case of using the Pope 1 model energy spectrum, the numerical quadrature solution of the second-order longitudinal structure function (8) is plotted in Fig. 1. In the sequent the analytical solution obtained for the second-order longitudinal structure function within the energy containing and inertial subranges is elucidated 14 . We consider the model spectrum (2) at very high Reynolds numbers. The very high wavenumber part of the spectrum is not represented, hence we set fη (ηκ ) = 1. The model spectrum is thus given by: −5/6 2/3 5/3

E(κ ) =CcL

ε

L

z2

(9)

11

[1 + z2 ] 6

−1/2

in which a scaled, non-dimensional wavenumber z has been introduced; z = cL


κ L.


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 38

The relation (8) becomes 4 −1/3 ⟨[δ v] ⟩1D,p (r) = CcL k 3 2

−1/2

in which κ r = zs, thus κ r = (cL

∫ ∞

z2

0

[1 + z2 ] 6

11

[ { }] sin(zs) cos(zs) 1−3 + dz (zs)3 (zs)2 −1/2

κ L)s and s = κ r/(cL

(10)

κ L). To determine the model energy

spectrum parameter cL the constraint (6) is employed: ∫ ∞

k= 0

−1/3

E(κ ) dκ = CcL

∫ ∞

z2

0

[1 + z2 ] 6

k

11

dz

(11)

The last integral can be solved analytically based on the Euler’s beta integral on the interval [0, ∞) 15 :

∫ ∞

I1 =

0

z2 11

[1 + z2 ] 6

dz =

3 Γ( 12 )Γ( 13 ) = 1.2620 10 Γ( 56 )

(12)

Hence, the value of cL is determined from the constraint (6): cL = (CI1 )3

(13)

The relation (8) thus yields: 4 k ⟨[δ v] ⟩1D,p (r) = k − 4 3 I1 2

[

∫ ∞

z2

0

[1 + z2 ] 6

11

] sin(zs) cos(zs) + dz (zs)3 (zs)2

(14)

The integral has been solved by symbolic MATLAB and the analytic solution yields:

I2 =

[

∫ ∞

z2

0

[1 + z2 ] 6

11

] sin(zs) cos(zs) I1 + dz = [T1 (r) + T2 (T3 (r)T4 (r) − T5 (r))] (zs)3 (zs)2 3

where 2 F T1 (r) = [s(r)]2

(( ) 1 ( ) ) 1 2 3 [s(r)]2 − , 3 2 4

T2 = 3

3/2

( ) 2 Γ 3


(15)

(16)

(17)

Page 9 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


T3 (r) = 27 × 2 1 T4 (r) = F 352π

1/3

2/3

[s(r)]

( ) 2 Γ 3

) (( ) 11 ( ) 7 6 17 [s(r)]2 , 3 6 4

22/3 T5 (r) = K 4 ([s(r)]) 2π [s(r)]2/3 3

(18)

(19)

(20)

where F, K and Γ denote the hypergeometric, Bessel, and gamma functions, respectively. Hence, the analytical result for the second-order longitudinal structure function describing the inertial and energy-containing subranges of turbulence is given by 4 ⟨[δ v]2 ⟩1D,p (r) = k {1 − [T1 (r) + T2 (T3 (r)T4 (r) − T5 (r))]} 3

(21)

The solution is shown in Fig. 1. Kolmogorov 7 ,8 proposed a second-order longitudinal structure function model valid for the inertial subrange of scales exclusively: ⟨[δ v]2 ⟩1D,p (r) = C1 (ε r)2/3

(22)

in which the value of the Kolmogorov parameter is approximated by C1 ≈ 2.0 at sufficiently large Reynolds number 11,12 . Employing the three-dimensional model energy spectrum of Pope 1 as basis, no exact algebraic expression is available for the second-order longitudinal structure function ⟨[δ v]2 ⟩1D,p (r) considering the dissipation subrange. However, there exist a few approximate semi-empirical algebraic relations for this structure function designed upon experimental, numerical and theoretical analyzes. There is of course no exact algebraic expression available for the structure function within the complete spectrum of turbulence as no analytical solution has been found for the dissipation subrange. In this work a formula for the second-order longitudinal structure function denoting a combi-



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 38

nation of the analytical solution valid for the energy containing and inertial subranges of scales given by (21) and the part of the semi-empirical formula proposed by Sawford and Hunt 5 valid for the dissipation subrange, is employed: ( 2 )2/3 r 4 {1 − [T1 (r) + T2 (T3 (r)T4 (r) − T5 (r))]} ⟨[δ v] ⟩1D,p (r) = k 2 3 rd + r2 2

(23)

where rd is a function of the integral-scale Reynolds number, ReL . The length scale relation rd = (15C1 )3/4 η governs the transition from the dissipation to the inertial subrange. The Kolmogorov structure function model parameter value is approximated by C1 ≈ 2.0. With minor differences, this model coincides with the numerical quadrature solution in Fig. 1.

Number density of turbulent vortices The number density of vortices may generally be expressed as a function of a consistent set of three-dimensional energy spectrum and second-order longitudinal structure function models 16 :

n˙ r (r) =

24E(2π /r) r5 ⟨[δ v]2 ⟩1D,p (r)

(24)

The conventional number density models found in the literature are limited to the inertial range of the energy spectrum. However, since the inertial subrange of turbulence might become very narrow as the integral scale Reynolds number (ReL ) decreases to finite values as observed in engineering equipment, the number density model must be expanded to the complete spectrum 2–4 . In this work the energy spectrum model of Pope 1 , considering both the energy containing, inertial and dissipation subranges of isotropic turbulence, was adopted. Moreover, a consistent second-order longitudinal structure function formula is required for the complete spectrum of turbulence. In previous attempts to extend the number density model to the entire spectrum of turbulence, the Kolmogorov second-order longitudinal structure function model which is limited to the inertial subrange was employed even in cases were the energy spectrum model was extended to consider 10 ACS Paragon Plus Environment

Page 11 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


the entire spectrum of turbulence 2 . Other researchers did employ a rough inconsistent structure function approximation instead in an attempt to extend the structure function model as well to a wider spectrum 3,4 . The number density model proposed by Ghasempour et al. 2 for the turbulent vortices supposedly valid for the entire energy spectrum is given in Table 1. The conventional number density model based on the Kolmogorov models valid for the inertial subrange of turbulence solely is given in the same table for comparison. The novel number density model proposed in this work for the turbulent vortices in the complete energy spectrum is also given in Table 1. In this case the generalized number density function (24) has thus been extended to the complete energy spectrum in a consistent manner using both an extended energy spectrum 1 and an extended second-order structure function (23). The number density models in Tab. 1 are plotted in Fig. 2. The given second order longitudinal structure function model predictions are in good qualitative agreement with recent DNS data 17 .

Fluid particle breakup and coalescence The drop and bubble breakage and coalescence phenomena are determining the dispersed phase size distribution in a variety of multiphase flow systems. Constitutive equations for the fluid particle breakage and coalescence phenomena are needed adopting the population balance equation 18,19 for describing the behavior of dispersed flow systems. Numerous closure models for fluid particle breakage and coalescence can be found in the literature. Reviews of these closure models are provided by, e.g., Solsvik et al. 20 and Liao and Lucas 21 ,22 . The breakage and coalescence closure models reported in the literature are derived considering the inertial subrange of turbulence exclusively. In the following, the breakup models proposed by Coulaloglou and Tavlarides 23 and Luo and Svendsen 24 and the coalescence model by Prince and Blanch 25 are extended to the entire spectrum of turbulence scales. The results obtained by the extended models are compared with the corresponding results from the original kernels which are limited to the inertial subrange of turbulence.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 38

Breakup model by Coulaloglou and Tavlarides The following simplifying assumptions were considered in the phenomenological drop breakup model proposed by Coulaloglou and Tavlarides 23 : (i) the drops are in a turbulent flow field which is locally isotropic, (ii) the drop size is within the inertial subrange (i.e. η ≪ d ≪ L and E(κ ) ∝ −5/3), (iii) viscous effects are negligible, and (iv) the drop deforms due to local pressure fluctuations (statistically treated analogous to the velocity correlation function Batchelor 26 ). The breakup criteria used in the model derivation is that an oscillating deformed drop will break if the turbulent kinetic energy transmitted to the drop by turbulent eddies exceeds the drop surface energy. The breakup frequency b of a drop of diameter d was thus given by [

] 1 1 ∆N(d) b(d) = [fraction of drops breaking] = breakup time τb N(d)

(25)

Considering Equation (25) the drop breakup time based on the inertial subrange of turbulence was given by 27

τb =

d d ≈√ ≈ c1 d 2/3 ε −1/3 2 v(d) ¯ ⟨[δ v] ⟩1D,p

(26)

where v(d) ¯ represents the mean translational velocity of drops with diameter d, roughly approximated by the second-order longitudinal structure function. The fraction of drops breaking was deduced from the two-dimensional Maxwell–Boltzmann speed distribution 20 :

∫∞ Ec

[ ] Ec ∆N(d) P(E) dE = exp − ¯ = N(d) E

(27)

In Equation (27) the critical turbulent kinetic energy was taken as the drop surface energy: E c = c2 σ d 2 Based on the formula for the drop translational kinetic energy defined by E¯ = 12 ρd


(28) (π

) 3 v¯2 , the d 6

Page 13 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


mean turbulent kinetic energy of a drop with diameter d was approximated by E¯ = c3 ρd d 3 ⟨[δ v]2 ⟩1D,p = c5 ρd ε 2/3 d 11/3

(29)

where the mean translational velocity of drops with diameter d, v(d), ¯ was roughly approximated by the second-order longitudinal structure function. Substitution of Equations (26) and (27) into (25); after use of Equations (28) and (29), yields the following expression for the drop breakup frequency:

b(d) = c6 d

−2/3 1/3

ε

[ exp −

c7 σ ρd ε 2/3 d 5/3

] (30)

Considering the drop breakup model proposed by Coulaloglou and Tavlarides 23 , the expressions for τb (26) and E¯ (29) are based on the inertial subrange of turbulence. A drop breakup time consistent for the complete turbulence energy spectrum is given by d c8 d τb = c7 ( )1/2 )1/2 = ( ( )2/3 ⟨[δ v]2 ⟩1D,p 4 d2 {1 − [T1 (d) + T2 (T3 (d)T4 (d) − T5 (d))]} 3 k r2 +d 2

(31)

d

in which the expression for ⟨[δ v]2 ⟩1D,p is taken from Eq. (23). Furthermore, the mean turbulent kinetic energy consistent with the entire energy spectrum yields )2/3 ( d2 3 2 34 ¯ {1 − [T1 (d) + T2 (T3 (d)T4 (d) − T5 (d))]} (32) E = c9 ρd d ⟨[δ v] ⟩1D,p = c10 ρd d k 2 3 rd + d 2 It is noticed that rd is a function of ReL .

Breakup model by Luo and Svendsen The breakup model by Coulaloglou and Tavlarides 23 includes two unknown parameter that need to be determined by fitting to experimental data. Luo and Svendsen 24 intended to develop a more fundamental fluid particle breakage model for turbulent dispersions employing basic kinetic theory 13 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 38

concepts. The breakage criterion proposed in the model derivation states that a fluid particle breaks if the turbulent kinetic energy of the bombarding vortices becomes greater than a critical value. The ratio of the partial breakup density function, Ω(d, d ′ ), to the number of fluid particles per unit volume, n, might be given as 20 : Ω(d, d ′ ) = n

∫

π Pb (d, d ′ , r)vrel (d + r)2 n˙ r dr 4

(33)

in which the fluid particle–vortex relative velocity vrel is approximated by the second-order lon( )1/2 gitudinal structure function relation; vrel ≈ ⟨[δ v]2 ⟩1D,p (r) . The number density was modeled in accordance with Eq. (24), with a correction factor (1 − αd ) in order to take into account the existence of the dispersed phase 28 . The probability for a drop or bubble of size d to break and form a daughter particle of size d ′ because of a collision with a vortex of size r was expressed by (

c f (d, d ′ )πσ d 2 Pb (d, d ′ , r) = exp − 1 ( π 3 ) ρc 2 6 r ⟨[δ v]2 ⟩1D,p (r)

) (34)

where the mean translational vortex velocity has been approximated by the second order longitudinal structure function, and c f (d, d ′ ) denotes the surface area increase coefficient given by ( ′ 3 )2/3 ( [ ′ 3 ])2/3 c f (d, d ′ ) = [dd 3] + 1 − [dd 3] − 1. Furthermore, the energy spectrum function E(κ ) employed is the Kolmogorov function (1). The total breakage density of drops or bubbles of size d is defined by

Ω(d) =

1/3 d/2 ∫

Ω(d, d ′ )

3[d ′ ]2 ′ dd d3

(35)

0

The breakage frequency b(d) is given in terms of Ω(d) by 20

b(d) =

Ω(d) 2n

(36)

The normalized daughter particle size distribution (volume based) was determined from the fol14 ACS Paragon Plus Environment

Page 15 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


lowing relation 20 : ∗ PDSD,V (Vd ,Vd ′ ) = d

Ω(Vd ,Vd ′ ) Ω(Vd )

(37)

The details for relating the volume based and diameter based distribution functions are outlined elsewhere 20 . The original breakage model by Luo and Svendsen 24 was restricted to the inertial subrange of turbulence using the Kolmogorov structure function and energy spectrum models. The extended version of the breakup model by Luo and Svendsen 24 proposed in this work for the complete spectrum of turbulence is achieved by adopting the relation (23) for the structure function and definition (2) for the energy spectrum in Eq. (33).

Coalescence model by Prince and Blanch Prince and Blanch 25 developed a bubble coalescence kernel based on the assumptions of isotropic turbulence and that the bubble sizes lies in the inertial subrange of turbulence. The bubble collision rate resulting from turbulent motion was expressed as a function of bubble size, concentration and velocity:

θi j = ni n j Si j vrel,i j

(38)

Here, Si j is the collision cross-sectional area of the bubbles defined by ] [ π di d j 2 Si j = + 4 2 2

(39)

The mean relative velocity between the two bubbles of sizes i, j is given by

vrel,i j = (v¯2i + v¯2j )1/2

(40)

The mean translational velocity of bubbles of size i was approximated by the Kolmogorov secondorder longitudinal structure function: 2/3

v¯2i = ⟨[δ v]2 ⟩1D,p (di ) = k1 ε 2/3 di 15

ACS Paragon Plus Environment

(41)


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 38

Substitution of Equations (39), (40) and (41) into (38) yields the following expression for θi j : ( )2 ( 2/3 2/3 )1/2 θi j = k2 ni n j ε 1/3 di + d j di + d j

(42)

The coalescence rate model by Prince and Blanch 25 is thus restricted to the inertial subrange of turbulence. To extend this collision rate model to the complete energy spectrum, expression (41) is replaced by (23). The result is: [ ]2 θi j = k 3 n i n j d i + d j

[( ( k

) )2/3 di2 {1 − [T1 (di ) + T2 (T3 (di )T4 (di ) − T5 (di ))]} rd2 + di2 1/2  ( )2/3 2 dj { } 1 − [T1 (d j ) + T2 (T3 (d j )T4 (d j ) − T5 (d j ))]  + k 2 2 rd + d j (43)

where rd is a function of ReL .

Results and discussion In the simulation study, the following parameters and relations are used if not otherwise specified:

ε = 1, L = 1, k = (Lε )2/3 , ν = k1/2 L/ReL , and η = (ν 3 /ε )1/4 .

Energy spectra The energy spectrum model of Pope 1 describes the complete spectrum of isotropic turbulence. In several publications employing the energy spectrum model proposed by Pope 1 the parameters cL and cη were erroneously fixed at the values representative for infinite integral scale Reynolds number flows. In order to fulfill the constrains (6) and (7), the parameters cL and cη in the energy spectrum model of Pope 1 (Eq. 2) must be refitted for each Reynolds number considered. The energy spectrum model of Pope 1 is thus rather computationally demanding. To reduce the com16 ACS Paragon Plus Environment

Page 17 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


putational load a possible approximation might be to make empirical Reynolds number dependent relations fitted to parameter variation profiles covering the range of Reynolds numbers of interest and solve the model energy spectrum function using parameter values from these parameterizations instead of refitting the parameter values for each flow condition by use of the given constraints. There is consensus in the literature that the energy spectrum model of Pope 1 might be employed for the entire energy spectrum of isotropic turbulence. The spectrum model has been validated comparing the model predictions with experimental data covering several flow types and conditions 1 . The performance of the spectrum model was considered very good. The spectrum model also compare well with direct numerical simulation (DNS) results as discussed by Han et al. 4 .

Structure function The population balance model is considered to be employed as part of the combined multifluidpopulation balance model framework, thus the numerical quadrature solution of the second-order longitudinal structure function relation (8) is considered too computational demanding for CFD simulations in which the Reynolds number might have substantial spatial and temporal variations. For this reason, in this work the numerical quadrature solution method for determining the structure function (8) in terms of the energy spectrum model of Pope 1 was exclusively applied to verify the accuracy of the semi-empirical algebraic relations. All the models considered for the second-order longitudinal structure function desribing the complete range of scales include attempts to take into account the impact of Reynolds number variations. Moreover, bearing in mind that there is an exact relationship between the three-dimensional energy spectrum model and the second-order longitudinal structure function model (8), the structure function can be computed accurately provided a sufficiently accurate energy spectrum model is available. In this article the energy spectrum model of Pope 1 was employed because the Reynolds number dependence of the model has been validated and found to represent experimental data accurately. For these reasons, it is assumed that the structure function computed based on this 17 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

spectrum model by use of a sufficiently accurate quadrature rule from the exact relationship (8) yields a fairly accurate approximation. However, the numerical quadrature solution is rather computationally demanding as part of any PBE computations, thus we did propose to use a simplified algebraic formula representing the structure function and intend to use the numerical solution to verify the accuracy of the simplified model. In Fig. 1 the Reynolds number dependency of the quadrature approximation of the secondorder longitudinal structure function (8) valid for the entire energy spectrum is shown. The corresponding results obtained by the simplified model are not plotted in this figure because merely minor differences exist between the numerical solution of (8) and the new model for the secondorder longitudinal structure function (23) and are thus hardly visible. A more detailed analysis of these models is thus provided for a few selected Reynolds numbers. In Figs. 3 and 4 it is shown that the semi-empirical model constructed by combining the analytical solution derived for the energy-containing and inertial subranges (21) and the dissipative part of the semi-empirical interpolation formula of Sawford and Hunt 5 is in very good agreement with the numerical quadrature approximation of the second-order longitudinal structure function (8). It is noted that for sufficiently high Reynolds numbers all the profiles are consistent and almost congruent within the inertial subrange. It is further seen that the Kolmogorov structure function valid for the inertial subrange cannot be accurately extrapolated neither to the energy containing nor the dissipation subranges. Moreover, the new approximate algebraic model coincides with acceptable accuracy to the quadrature approximation for all the Reynolds numbers considered. Most important, it is shown that the extend of the inertial subrange increases with increasing Reynolds numbers and for integral scale Reynolds numbers less than about 450 this subrange might even disappear 18 . The Kolmogorov structure function model is surely limited to sufficiently high Reynolds numbers and this model must be avoided for finite integral scale Reynolds numbers less than 450 but it should also be noted that the inertial subrange is very narrow even for integral scale Reynolds numbers around 1000. It is further seen that both the quadrature approximation of Eq. (8) and the combined analytical and semi-empirical formula approach 2v2 for large r in accordance with the turbulent


Page 18 of 38

Page 19 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


kinetic energy constraint given by Eq. (6). This constraint is by far not fulfilled by extrapolation of the Kolmogorov formula valid for the inertial subrange solely. For small r, the structure function should approach the lower limit value ε r2 /15ν . This constraint is certainly not fulfilled by extrapolation of the Kolmogorov formula restricted to the inertial subrange but both the quadrature approximation of Eq. (8) and the combined analytical and semi-empirical formula approach the lower limit value with reasonable accuracy.

Number density of turbulent vortices The consistently extended number density model applicable for the complete range of scales of turbulence is compared with some of the previous models that are designed considering an extended model energy spectrum only 2 (but extrapolating the Kolmogorov structure function formula) (Tab. 1). The comparison of the model predictions is shown in Fig. 2. From Fig. 2 it is seen that, for small r, the results obtained with the vortex number density model designed by inserting the Kolmogorov formulas solely into (24) differ significantly from the extended number density function designed for the entire range of scales by inserting (i) the energy spectrum model for the entire spectrum exclusively (and extrapolating the Kolmogorov structure function formula) as was done by Ghasempour et al. 2 and alternatively (ii) the energy spectrum and second-order longitudinal structure function models for the complete spectrum of turbulence as proposed in this work. There is a small but noticeable difference between the number density predictions obtained by model designs (i) and (ii). These differences in the model results increase slightly with decreasing Reynolds number.

Breakage and coalescence modeling The second-order longitudinal structure function has generally been adopted as a useful approximation for the turbulent fluid particle velocity, the relative fluid particle velocity, the relative vortex–particle velocity, etc without proper validation thus the accuracy of this approach is presently not known. Hence, validations of these assumptions are required but considered outside the scope 19 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

of this article and left for future investigations. In the following the behavior of the breakup model by Coulaloglou and Tavlarides 23 is examined. Fig. 5 holds the breakage frequency predicted by the original 23 model restricted to the inertial subrange and the extended model applicable for the complete spectrum of turbulence. The physical system considered is characterized by ρd = 750 kg/m3 , ε = 8.5 m2 /s3 and σ = 0.053 N/m. Further, the empirical parameters are identical in the original and extended models: c6 = c11 = 0.003 and c7 = c12 = 1. With these conditions, the original and extended models give different predictions of the breakup frequency within the mother particle size range investigated. The extended breakup model is a function of the Reynolds number whereas the original model is not. For ReL = 105 , the mother particles with sizes in the range of 1–13 cm falls in the inertial subrange. Reducing the Reynolds number to ReL = 104 , the mother particles must be larger than 4 cm in order to fall in the inertial subrange. In the case of ReL = 103 , the inertial subrange does not exist in the range of mother particle sizes up to 14 cm. The assumption of adopting the second-order longitudinal structure function which is restricted to the inertial subrange in breakage modeling is thus strictly not valid for many practical applications. For example, in the experimental work on liquid–liquid dispersions by Becker et al. 29 ,30 the dispersed drops were much smaller than 1 cm (i.e., a size comparable to the lower limit vortex scale in the inertial subrange for very high Reynolds numbers, ReL = 105 ) in a stirred tank with a mean energy dissipation rate of the whole tank in the order of unity. On the other hand, holding c11 and c12 fixed in the extended breakup frequency model (i.e. c11 = 0.003 and c12 = 1), similar results of the breakup frequency can be obtained with the original Coulaloglou and Tavlarides 23 model fitting the c6 and c7 parameter values (i.e. ReL = 105 : c6 = 0.0038, and c7 = 1, ReL = 104 : c6 = 0.0035, and c7 = 5, and ReL = 103 : c6 = 0.0028, and c7 = 29). The results are given in Fig. 5. This observation is consistent with the fact that the extended breakup model is Reynolds number dependent whereas the original model is not. Moreover, this trend means that for a fixed Reynolds number, it is possible to refit the original breakup model parameters to comply with the


Page 20 of 38

Page 21 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


extended breakup models. It is noted that the lack of Reynolds number dependence might explain the numerous sets of parameter values for the original breakup models reported in the literature as reviewed by Solsvik et al. 20 . A major aim for the current work was to extend the models and thereby reduce the need for refitting the model parameters for different Reynolds number flows thus less costly experiments are required. However, the extended breakage kernel Reynolds number dependency still requires extensive validation. In the sequel the behavior of the breakage kernel proposed by Luo and Svendsen 24 is examined. For ReL = 105 , Fig. 6 shows the simulated results of the breakage frequency as obtained by use of the original Luo and Svendsen 24 kernel restricted to the inertial subrange solely and the extended kernel applicable for the complete spectrum of turbulence. A significant difference in observed in the simulated breakage frequency results when using the two model approaches. The physical system properties are given by µc = 0.001 Pa s, ρc = 1000 kg/m3 , ε = 8.5 m2 /s3 , σ = 0.053 N/m,

ρd = 750 kg/m3 3 . The simulated breakage frequency results achieved by the extended kernel is also largely sensitive to the Reynolds number. However, this kernel sensitivity to the Reynolds number does not occur for the original kernel which is strictly only valid for infinite Reynolds numbers. Both the energy spectrum and the second order longitudinal structure function models employded defining the extended breakage rate kernel are considered valid for both the dissipative, inertial and energy containing ranges of turbulence. In the integral of Eq. (33) the Kolmogorov microscale, η is the natural lower limit and the upper limit should be infinity when the complete range of turbulence scales are considered. It is noted that to produce these figures from the Luo and Svendsen 24 model the turbulence structure function, spectrum and vortex number density are √ determined by specifying L = 1 and ReL = 104 or 105 , and computing k = [Lε ]2/3 , ν = kL/ReL , and η = [ν 3 /ε ]1/4 . On the other hand, in the present simulation study the lower integral limit has been considered independent of the Reynolds number and computed as η = η (µc = 0.001Pa s). The upper limit is set equal to 2 times the mother fluid particle size. In Fig. 7 the corresponding normalized daughter size distribution function profiles have been



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

plotted. It is seen that the predictions of the original model of Luo and Svendsen 24 deviate considerably compared to those of the extended version of the model. However, the daughter sized distribution function seems to be independent of Reynolds number. This trend is reasonable as the Reynolds number dependency cancels out in the model equation (37). Adjustable parameters in the original breakup models can thus to some extent overcome severe theoretical limitations in the model framework. However, the major aim in this work is that the novel understanding of the statistical turbulence theory might lead to more general closures resolving further physical details and containing more universal parameters having constant values for larger ranges of Reynolds numbers minimizing the requirement for expensive experimental investigations as needed for parameter fitting and model validation. The fluid particle collision rate model by Prince and Blanch 25 is plotted in Fig. 8. The result is obtained by letting k2 = k3 (Eqs. (42) and (43)). Small differences are observed in the predicted collision rate. However, as the Reynolds number decreases the inertial subrange are present only for relatively large fluid particles (as previously discussed for the breakup frequency). Similarly to what was observed for the breakup models discussed earlier, adjusting the parameters k2 in Eq. (42) and k3 in Eq. (43) at a fixed Reynolds number allows the original kernel limited to the inertial subrange and the extended model applicable to the complete range of turbulence to predict the same fluid particle collision rate θi j . However, a major goal of this work is to achieve more general descriptions containing universal parameter values thus having constant values for larger ranges of Reynolds numbers and the need for expensive experimental investigations to enable parameter fitting and model validation might be reduced. To validate the closure model extensions well planned model based experiments are needed. Proper model validation requires both characterizations of the flow determining the local turbulent energy dissipation rate, turbulent energy, integral length scale in order to compute the local integral scale Reynolds number, and measuring single fluid particle breakage characteristics 31–34 and investigating binary coalescence events. Such detailed experimental data is not available yet, even though several research groups are providing local experimental data on the fluid particle size


Page 22 of 38

Page 23 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


characteristics, but normally without a suitable characterization of the local flow properties. In our research group there is an ongoing project in which we intend to design a suitable apparatus for performing the necessary experiments to provide sufficient experimental data for proper model validation. This work is considered to be outside the scope of this article which instead focus on the formal derivation of the population balance closure models based on more rigorous and extended turbulence theory. Pope 1 did show that the intermittency effects have a much smaller influence than the Reynolds number variations upon the energy spectrum model. Hence, because of the uncertainty in the energy spectrum model Reynolds number dependency, the intermittency effects are neglected in this work. The accuracy and validity reflected by the new population balance equation closure models proposed thus rely on the accuracy of the model energy spectrum employed. Any additional limitations of the original closure laws considered are not improved in this work, emphasis is placed exclusively on the extension of the model turbulence spectrum and the derivation of a consistent structure function model.

Conclusion It has been shown that the Kolmogorov second-order longitudinal structure function model is limited to the inertial subrange of turbulence and is thus not applicable for the dissipation and energy containing subranges with sufficient accuracy. Hence, outside the inertial subrange the Kolmogorov stucture function is not consistent with the spectrum model proposed by Pope 1 . It is noted that even in the inertial subrange of turbulence the scaling exponents might deviate from the Kolmogorov scaling of the structure function and the energy spectrum because of insufficiently high Reynolds numbers, flow anisotropy, intermittency or possibly other factors. 35 It is further concluded that the Kolmogorov structure functions and the energy spectrum parameter values are dependent on Reynolds number thus for industrial flows with finite Reynolds numbers these parameter values cannot be considered universal. Likewise, the power laws of Kol-



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

mogorov are not rigorous and might change significantly. For this reason, when modeling fluid particle coalescence and breakage induced by turbulence, the parameter values from turbulence theory cannot be considered universal and are not constants. In this work the turbulent energy spectrum and the structure function Reynolds number dependency is considered and accounted for by the use of exact relations from statistical turbulence theory of isotropic turbulence and the approximate three-dimensional model energy spectrum of Pope 1 . The simulated results achieved by the extended models for the vortex number density, fluid particle breakage frequency, daughter size distribution and collision rate are significantly different compared to those obtained by the original models. The original models were limited to the inertial range of turbulence whereas the extended models are applicable for a wider range of Reynolds numbers. For any given ReL , the deviations in model predictions might in many cases be accounted for by re-fitting the parameters in the original models which are limited to the inertial subrange. However, the integral scale Reynolds number dependency of the model parameters is only taken into account by the extended models proposed in this work. Improved knowledge and understanding of the turbulent phenomena might lead to more general models containing more universal parameters having constant values in larger ranges of Reynolds numbers. Such a modeling approach will reduce the need for gathering costly experimental data as required for parameter fitting and model validation. Among the weakest links in the derivation of the extended closure models applicable for the complete turbulence spectrum, which remain to be elucidated, are whether the second-order longitudinal structure function procures a sufficiently accurate estimate for the relative velocity of fluid particles, the relative vortices–particle velocity, vortex translational velocity, etc. for which this quantity has been used. Moreover, model discrimination of the coalescence and breakage mechanisms and criteria also require further attention. There is a general requirement for model based experimental investigations to validate the Reynolds number dependency of the underlying closure models, the breakage and coalescence criteria, and in order to discriminate the physical mechanisms determining the basis of the different


Page 24 of 38

Page 25 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


model closures.

Acknowledgment The financial support by the Norwegian Research Council through the GASSMAKS program is gratefully appreciated. The discussion on the analytical solution of the structure function with Prof. Stephen B. Pope (Sibley School of Mechanical and Aerospace Engineering, Cornell University) is gratefully acknowledged. We also appreciate the discussion on turbulence theory with Prof. Emeritus Per-Åge Krogstad (Department of Energy and Process Engineering, The Norwegian University of Science and Technology (NTNU)).

References (1) Pope, S. B. Turbulent Flows; Cambridge University Press, Cambridge, UK, 2000. (2) Ghasempour, F.; Andersson, R.; Andersson, B.; Bergstrom, D. J. Number Density of Turbulent Vortices in the Entire Energy Spectrum. AIChE J. 2014, 60, 3989. (3) Han, L.; Gong, S.; Li, Y.; Gao, N.; Fu, J.; Luo, H.; Liu, Z. Influence of Energy Spectrum Distribution on Drop Breakage in Turbulent Flows. Chem. Eng. Sci. 2014, 117, 55. (4) Han, L.; Gong, S.; Ding, Y.; Fu, J.; Gao, N.; Luo, H. Consideration of Low Viscous Droplet Breakage in the Framework of the Wide Energy Spectrum and the Multiple Fragments. AIChE J. 2015, 61, 2147. (5) Sawford, B. L.; Hunt, J. C. R. Effects of Turbulence Structure, Molecular Diffusion and Source Size on Scalar Fluctuations in Homogeneous Turbulence. J. Fluid Mech. 1986, 165, 373. (6) Ghasempour, F. Structures, Properties and Dynamics of Turbulent Vortices. Ph.D. thesis, Chalmers Univerity of Technology, Gothenburg, Sweden, 2015. 25 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(7) Kolmogorov, A. N. Dissipation of Energy in the Locally Isotropic Turbulence. Dokl. Akad. Nauk SSSR 1941, 32, 19, [in Russian]. (8) Kolmogorov, A. N. The Local Structure of Turbulence in an Incompressible Fluid for Very Large Reynolds Numbers. Dokl. Akad. Nauk SSSR 1941, 30, 299, [in Russian]. (9) Kolmogorov, A. N. Dissipation of Energy in the Locally Isotropic Turbulence. In Proceedings: Mathematical and Physical Sciences, Vol. 434, No. 1890, Turbulence and Stochastic Process: Kolmogorov’s Ideas 50 Years On (Jul. 8, 1991), Royal Society Publishing; London, 1991; pp 15–17 . (10) Kolmogorov, A. N. The Local Structure of Turbulence in an Incompressible Fluid for Very Large Reynolds Numbers. In Proceedings: Mathematical and Physical Sciences, Vol. 434, No. 1890, Turbulence and Stochastic Process: Kolmogorov’s Ideas 50 Years On (Jul. 8, 1991), Royal Society Publishing; London, 1991; pp 9–13 . (11) Sreenivasan, K. R. On the Universality of the Kolmogorov Constant. Phys. Fluids 1995, 7, 2778. (12) Ni, R.; Xia, K.-Q. Kolmogorov Constants for the Second-Order Structure Function and the Energy Spectrum. Phys. Rev. E 2013, 87, 023002. (13) Davidson, P. A. Turbulence: An Introduction for Scientists and Engineers; Oxford University Press, New York, USA, 2004. (14) Pope, S. B. Personal Communications. 2015; Professor at Cornell University, Ithaca, NY. (15) Gasper, G. q-Extensions of Barnes’, Cauchy’s, and Euler’s Beta Integrals. In: Topics in Mathematical Analysis; Rassias, T. M., Ed.; Series in Pure Mathematics; World Sci. Publ. Co., Singapore., 1989; Vol II; pp 294–314 .


Page 26 of 38

Page 27 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(16) Solsvik, J.; Jakobsen, H. A. A Review of the Statistical Turbulence Theory Required Extending the Population Balance Closure Models to the Entire Spectrum of Turbulence. AIChE J. 2015, DOI: 10.1002/aic.15128. (17) Ishihara, T.; Gotoh, T.; Y, K. Study of High-Reynolds Number Isotropic Turbulence by Direct Numerical Simulation. Annu. Rev. Fluid Mech. 2009, 41, 165. (18) Solsvik, J.; Jakobsen, H. A. The Foundation of the Population Balance Equation – A Review. J. Dispersion Sci. Technol. 2015, 36, 510. (19) Yeoh, G. H.; Cheung, C. P.; Tu, J. Multiphase Flow Analysis Using Population Balance Modeling: Bubbles, Drops and Particles; Butterworth–Heinemann, Amsterdam, 2014. (20) Solsvik, J.; Tangen, S.; Jakobsen, H. A. On the Constitutive Equations for Fluid Particle Breakage. Rev. Chem. Eng. 2013, 29, 241. (21) Liao, Y.; Lucas, D. A Literature Review of Theoretical Models for Drop and Bubble Breakup in Turbulent Dispersions. Chem. Eng. Sci. 2009, 64, 3389. (22) Liao, Y.; Lucas, D. A Literature Review on Mechanisms and Models for the Coalescence Process of Fluid Particles. Chem. Eng. Sci. 2010, 65, 2851. (23) Coulaloglou, C. A.; Tavlarides, L. L. Description of Interaction Processes in Agitated Liquid– Liquid Dispersions. Chem. Eng. Sci. 1977, 32, 1289. (24) Luo, H.; Svendsen, H. F. Theoretical Model for Drop and Bubble Breakup in Turbulent Dispersions. AIChE J. 1996, 42, 1225. (25) Prince, M. J.; Blanch, H. W. Bubble Coalescence and Break-Up in Air-Sparged Bubble Columns. AIChE J. 1990, 36, 1485. (26) Batchelor, G. K. Pressure Fluctuations in Isotropic Turbulence. Proc. Cambridge Philos. Soc. 1951, 47, 359. 27 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(27) Kolmogorov, A. N. On the Breakage of Drops in a Turbulent Flow. Dokl. Akad. Nauk SSSR 1949, 66, 825, [in Russian]. (28) Xing, C.; Wang, T.; Guo, K.; Wang, J. A Unified Theoretical Model for Breakup of Bubbles and Droplets in Turbulent Flows. AIChE J. 2015, 61, 1391. (29) Becker, P. J.; Puel, F.; Henry, R.; Sheibat-Othman, N. Investigation of Discrete Population Balance Models and Breakage Kernels for Dilute Emulsification Systems. Ind. Eng. Chem. Res. 2011, 50, 11358. (30) Becker, P. J.; Puel, F.; Chevalier, Y.; Sheibat-Othman, N. Monitoring Silicone Oil Droplets During Emulsification in Stirred Vessel: Effect of Dispersed Phase Concentration and Viscosity. Can. J. Chem. Eng. 2013, 92, 296. (31) Maaß, S.; Kraume, M. Determination of Breakage Rates Using Single Drop Experiments. Chem. Eng. Sci. 2012, 70, 146. (32) Maaß, S.; Buscher, S.; Hermann, S.; Kraume, M. Analysis of Particle Strain in Stirred Bioreactors by drop breakage investigations. Biotechnol. J. 2011, 6, 1. (33) Andersson, R.; Andersson, B. On the Breakup of Fluid Particles in Turbulent Flows. AIChE J. 2006, 52, 2020. (34) Solsvik, J.; Jakobsen, H. A. Single Drop Breakup Experiments in Stirred Liquid–Liquid Tank. Chem. Eng. Sci. 2015, 131, 219. (35) Zhang, J.; Xu, M.; Pollard, A.; Mi, J. Effects of External Intermittency and Mean Shear on the Spectral Inertial-Range Exponent in a Turbulent Square Jet. Phys. Rev. E 2013, 87, 053009.


Page 28 of 38

Page 29 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Table 1: Number density models proposed for different energy spectrum subranges. Model

Subrange

n˙ r (r) =

24C r−4 C1 (2π )5/3

Inertial subrange

n˙ r (r) =

24C r−4 fη (2πη /r) fL (2π L/r) C1 (2π )5/3

Entire energy spectrum proposed by Ghasempour et al. 2

{

n˙ r (r) = (2π )5/3

24Cε 2/3 r−10/3 fL (2π L/r) fη (2πη /r) } ( )2/3 r2 4 {1−[T1 (r)+T2 (T3 (r)T4 (r)−T5 (r))]} 2 2 3k

Entire energy spectrum proposed in this study

r +r d



10

0

h(δv)2i

Decreasing ReL

ReL = 20

10

-2

ReL = 100 ReL = 450 ReL = 1000 ReL = 10000 ReL = 100000

10

-4

10

-4

10

-3

10

-2

10

-1

10

0

r

Figure 1: Continuous lines denote the numerical quadrature solutions of the structure function integral formula (8) applicable for the entire spectrum of turbulence. The dashed line represents the exact solution (21) in the energy containing and inertial subranges (extrapolated to r-values in the dissipating subrange).

Entire energy spectrum, Re =100000

Entire energy spectrum, Re =450

L

15

10

10

Number density [#/m3m]

Ghasempour this work inertial subrange

10

5

10

0

10

−5

10

L

10

10

Number density [#/m3m]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 30 of 38

Ghasempour this work inertial subrange

5

10

0

10

−5

10

−10

−4

10

−2

0

10

10

10

2

10

−2

10

r

−1

0

10

10

1

10

r

Figure 2: The vortex number density n˙ r models in Tab 1 are compared (C1 = 2 and C = 1.5).


Page 31 of 38

Re =100000

Re =1000

L

1

L

1

10

10

0

10

0

10

10

〈(δv)2〉

〈(δv)2〉

−1

−2

10

−1

10

−2

−3

−4

10

10

inertial range formula new model numerical quadrature

10

−4

10

−2

0

10

10


−3

10

2

10

−3

10

−1

0

10

10

1

10

L

10

0

10

0

〈(δv)2〉

10 −1

10

−1

10 −2

10


−3

10

−2

10

−1

0

10

10


−2

10

1

0

10

10

r

r

Figure 3: Comparison of the structure function models in log/log-plots.


1

10

r Re =20

L

1

−2

10

r Re =100

〈(δv)2〉

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60



Re =100000

Re =100

L

L

4

4 inertial range formula new model numerical quadrature

3

〈(δv)2〉

3

〈(δv)2〉

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

2

1

0 −4 10

Page 32 of 38


2

1

−2

0

10

10

0 −2 10

2

10

r

−1

0

10

10

r

Figure 4: Comparison of the structure function models in log/linear-plots.


1

10

Re =100000 L

0.25

original model original model − fitted parameters new model

b(d) [1/s]

0.2 0.15 0.1 0.05 0 0

0.05

0.1

d [m] Re =10000 L

0.25


0.2

b(d) [1/s]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


0.15 0.1 0.05 0 0

0.05

0.1

d [m] Re =1000 L

0.25


0.2

b(d) [1/s]

Page 33 of 38

0.15 0.1 0.05 0 0

0.05

0.1

d [m]

Figure 5: Breakup frequency model by Coulaloglou and Tavlarides 23 . Comparison between original model (inertial subrange) and extended model (entire spectrum of turbulence) for different Reynolds numbers. Empirical parameters are (i) held identical in the two model frameworks and (ii) adjusted to give the same prediction of the breakup frequency for the two model approaches. The inertial subrange is indicated by dotted lines in the top figure, in the middle figure the dotted line indicates the border between the dissipation and inertial subranges, whereas in the bottom figure the dissipation subrange prevails over whole size interval considered in the figure.



b(d) [1/s]

102

100

10-2

inertial range formula wide spectrum,Re L=10

5

wide spectrum,Re =10 4 L

10

-4

2

4

6

8

10

d [mm]

Figure 6: Comparison between the predictions achieved by the breakup frequency kernel proposed by Luo and Svendsen 24 valid in the inertial subrange and the extended model valid for the entire spectrum of turbulence. The extended breakup model sensitivity to ReL is illustrated. The size range of the mother particle corresponds to the dissipation subrange.

10 3 Inertial range formula wide spectrum,Re =10 5 L

10 2

wide spectrum,Re =10

4

L

Redistribution

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

10 1

10 0

10 -1

10 -2

0

0.2

0.4

0.6

0.8

1

f v [-]

Figure 7: Normalized volume based daughter size redistribution function by Luo and Svendsen 24 . Comparison between the original model valid for the inertial subrange and the extended model applicable for the complete spectrum of turbulence, and sensitivity to ReL for the extended redistribution function model valid for the complete spectrum of turbulence.


Page 34 of 38

Page 35 of 38

ReL=100000 0.03 original model new model 0.025

3

θij/ninj [m /s]

0.02

0.015

0.01

0.005

0 0

0.02

0.04

0.06

0.08

0.1

0.12

0.14

d [m] ReL=10000 0.03 original model new model 0.025

3

θij/ninj [m /s]

0.02

0.015

0.01

0.005

0 0

0.02

0.04

0.06

0.08

0.1

0.12

0.14

d [m]

ReL=1000 0.1 original model new model

0.09 0.08 0.07

θij/ninj [m3/s]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


0.06 0.05 0.04 0.03 0.02 0.01 0 0

0.02

0.04

0.06

0.08

0.1

0.12

0.14

d [m]

Figure 8: The fluid particle collision rate model by Prince and Blanch 25 . Comparison between the original model valid for the inertial subrange and the extended model to the entire spectrum of turbulence. The inertial subrange is indicated by dotted lines in the top figure, in the middle figure the dotted line indicates the border between the dissipation and inertial subranges, whereas in the bottom figure the dissipation subrange prevails over whole size interval considered in the figure.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure captions: Figure 1 Continuous lines denote the numerical quadrature solutions of the structure function integral formula (8) applicable for the entire spectrum of turbulence. The dashed line represents the exact solution (21) in the energy containing and inertial subranges (extrapolated to r-values in the dissipating subrange). Figure 2 The vortex number density n˙ r models in Tab 1 are compared (C1 = 2 and C = 1.5). Figure 3 Comparison of the structure function models in log/log-plots. Figure 4 Comparison of the structure function models in log/linear-plots. Figure 5 Breakup frequency model by Coulaloglou and Tavlarides 23 . Comparison between original model (inertial subrange) and extended model (entire spectrum of turbulence) for different Reynolds numbers. Empirical parameters are (i) held identical in the two model frameworks and (ii) adjusted to give the same prediction of the breakup frequency for the two model approaches. The inertial subrange is indicated by dotted lines in the top figure, in the middle figure the dotted line indicates the border between the dissipation and inertial subranges, whereas in the bottom figure the dissipation subrange prevails over whole size interval considered in the figure. Figure 6 Comparison between the predictions achieved by the breakup frequency kernel proposed by Luo and Svendsen 24 valid in the inertial subrange and the extended model valid for the entire spectrum of turbulence. The extended breakup model sensitivity to ReL is illustrated. The size range of the mother particle corresponds to the dissipation subrange. Figure 7 Normalized volume based daughter size redistribution function by Luo and Svendsen 24 . Comparison between the original model valid for the inertial subrange and the extended model applicable for the complete spectrum of isotropic turbulence, and sensitivity to ReL for the extended redistribution function model valid for the complete spectrum of turbulence.


Page 36 of 38

Page 37 of 38

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 8 The fluid particle collision rate model by Prince and Blanch 25 . Comparison between the original model valid for the inertial subrange and the extended model to the entire spectrum of turbulence. The inertial subrange is indicated by dotted lines in the top figure, in the middle figure the dotted line indicates the border between the dissipation and inertial subranges, whereas in the bottom figure the dissipation subrange prevails over whole size interval considered in the figure.



graphical illustraction (TOC Graphic) 100

10-2

E(κ)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

10-4

Inertial subrange Energy containing subrange Dissipation subrange

10

-6

10-1

100

101

102

103

104

105

κ


Page 38 of 38

Development of Fluid Particle Breakup and ... - ACS Publications

Recommend Documents