3652
Ind. Eng. Chem. Res. 2008, 47, 3652-3663
Dispersion of Nanoparticle Clusters in a Rotor-Stator Mixer Jerzy Bałdyga,* Wojciech Orciuch, Łukasz Makowski, and Katarzyna Malik Faculty of Chemical and Process Engineering, Warsaw UniVersity of Technology, Warsaw, Poland
Gu1 l O 2 zcan-Tas¸ kin, Warren Eagles,† and Gustavo Padron
Downloaded via DURHAM UNIV on July 3, 2018 at 13:03:33 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.
BHR Group Ltd, The Fluid Engineering Centre, Cranfield, Bedfordshire MK43 0AJ, UK
The unique properties of nanoparticles and nanoparticle clusters show high potential for nanomaterials to be formulated into numerous products. In this paper, nanosuspensions are formulated by breaking up nanoparticle clusters (called agglomerates) in high-shear flows. A new breakage model is introduced to interpret erosive dispersion of agglomerates, and the population balance modeling is applied to account for effects of breakage on agglomerate size distribution. Effects of suspension structure on its rheology and flow are included in modeling. The population balance equations are solved using the quadrature method of moments (QMOM) that is linked directly to the k- model of the computational fluid dynamics (CFD) code FLUENT. In dispersion experiments, the aqueous suspensions of fumed silica particles, Aerosil 200V, are used. The test rig consists of an in-line Silverson rotor-stator mixer and a stirred tank. The head is a two-stage rotor-stator design with the inner stator consisting of round holes and the outer stator consisting of smaller square holes. Experimental results are compared with model predictions. 1. Introduction The unique properties of nanoparticles and nanoparticle clusters and their suspensions as well as recent advances allowing the manipulation of these have shown huge potential for nanomaterials to be formulated into numerous products. Examples of commercialized products with improved performance that could not be achieved previously include scratch/abrasionresistant transparent coatings, sunscreen lotions providing UV protection, polishing slurries providing pristine surfaces for optics, and environmental catalysts to reduce pollution. It is expected that within 5 to 10 years, applications will be much broader, spanning from healthcare and medicine to electronics. The manufacture of products that consist of nanoparticle suspensions requires the incorporation of nanoparticles in the liquid phase, breakup of nanoparticle clusters, and stabilization.1 Often the process is carried out in several devices run in series to achieve the desired final product quality. In order for nanoparticles to be highly functional in the final product, it is critical that a fine and stable dispersion is produced. In this paper, the findings of both an experimental and numerical study on dispersing nanoparticle clusters using an in-line rotor-stator are presented. The objectives of the study have been to determine the effect of operating parameters and identify the mechanisms of breakup using an in-line rotor-stator mixer. On the basis of experimental evidence that breakup of Aerosil 200V particles occurs through erosion, a new model of erosive breakup is introduced in this work, and numerical investigations of Aerosil 200V deagglomeration are performed. Breakup that controls agglomerate size distribution depends on the agglomerate structure and is caused by stresses generated by fluid deformation and inertia. Of course, the stresses depend on the structure of the flow field. On the other hand, rheological properties of the suspension depend on suspension structure, which affects the flow field. These effects are included in the model. * To whom correspondence should be addressed. Telephone: +48 22 234 6376. Telefax: +48 22 825 6037. E-mail: baldyga@ ichip.pw.edu.pl. † Present address: GlaxoSmithKline Research and Development Limited, Gunnels Wood Road, Stevenage, Herfordshire SG1 2NY, UK.
Figure 1. Schematic of rotor-stator rig.
To describe these phenomena and to simulate the process, a complex model based on the population balance equation is used in this paper. The population balance with adequate breakage and restructuring terms is solved using the quadrature method of moments (QMOM)2 that is linked to the commercial computational fluid dynamics (CFD) code Fluent supplemented with user defined functions (UDFs) for suspension rheology. 2. Experimental Details The rotor-stator used in the experimental and numerical work is shown in Figure 1. The test rig consisting of an in-line Silverson rotor-stator mixer, Figure 2, and a stirred tank (T ) 0.61 m) was arranged as shown in Figure 1. The slurry was kept suspended off the tank base with a down-pumping pitched blade turbine (D ) 0.42T, blades inclined at 45°, C ) T/4); rotated at 155 rpm. The corresponding average power input was 0.3 to 3% of the power input by the rotor-stator depending on the operating condition. The head used in this work (Figure 2) was a two-stage rotor-stator design with the inner stator consisting of round holes (general purpose disintegrating head) and the outer stator consisting of smaller square holes (square hole high shear screen). The test material used in this work was an aqueous suspension of fumed silica particles, Aerosil 200V, manufactured by
10.1021/ie070899u CCC: $40.75 © 2008 American Chemical Society Published on Web 01/08/2008
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008 3653
Figure 2. Geometry of the rotor-stator mixer. Table 1. Physicochemical Properties of Aerosis 200V Silica Particles (Provided by Degussa AG) property
Aerosil 200V m2
g-1
specific surface area (BET), tapped density, g L-1 moisture (2 h at 105 °C), wt % pH (in 4% dispersion) SiO2 content (based on ignited material) wt %
200 ( 25 approximately 120 e1.5 3.7-4.7 g99.8
Degussa. Fumed silica is produced by the continuous flame hydrolysis of silicon tetrachloride (Degussa technical bulletin 11).3 Aerosil 200V is a densified hydrophilic fumed silica with a specific surface area of 200 m2 g-1 (average primary particle size 12 nm). The densification process consists of applying mechanical pressure under vacuum in order to compact the material and increase its tapped density from 50-60 g dm-3 to approximately 120 g dm-3. Aerosil 200V has many applications in paints and coatings, adhesives, sealants, and printing inks, among others. It is used as an antisettling and antisagging agent or rheology modifier. Table 1 shows typical values for several physicochemical properties of Aerosil 200V. Experiments were performed by varying the rotor speed (3000-9000 rpm) at a fixed flow rate for different Aerosil 200V concentrations (1-20% by weight) (Figure 3). Samples were taken periodically during the experiment. The particle size distribution of the samples was determined using a Malvern Mastersizer S. The morphology of a limited number of samples were investigated under an electron microscope, SFEG SEM (Schottky field-emission gun scanning electron microscope). The rheology of Aerosil 200V suspensions were reported previously:4 at low concentrations, i.e., 1 and 5%, they are almost waterlike; at 10, 15, and 20 wt %, the viscosity is 4.7 × 10-3, 1.9 × 10-2, and 1.0 × 10-1 Pa s, respectively. This is in line with the findings of Chen et al.5 who measured the viscosity of Aerosil 200 suspensions of up to 6 vol % (∼13% by weight). A typical evolution of the particle size distribution is shown in Figure 4. The bimodal particle size distribution suggests that small size aggregates are eroded from large agglomerates. Similar findings were reported earlier for different systems.6,7 This is supported by images obtained using the electron microscopy. Figure 5 shows both large agglomerates and also submicrometer size aggregates of Aerosil 200V, and in Figure 6, aggregates of about 50-200 nm can be seen which correspond to the lowest size range measured by the particle sizing equipment.
Experimental results on the effects of process parameters, such as suspension concentration, rotor speed, and number of tank turnovers, on agglomerate size will be presented in Section 4 of this paper together with model predictions. 3. Model Presentation In this section, the model of the structure of Aerosil 200V agglomerates is presented. On the basis of the agglomerate structure, the strength of agglomerates is estimated and a fragment distribution function is derived. This allows to formulate a population balance equation for fractal agglomerates and express the effective suspension viscosity as a function of the resulting agglomerate size distribution. Results of application of the complete model will be presented and tested by performing simulation for a simple one-dimensional plug flow system. Afterward, the model will be implemented into a CFD environment to model Aerosil agglomerate dispersion in the rotor-stator mixer. 3.1. Structure, Strength, Breakage, and Fragment Distribution of Agglomerates. Aerosil 200V has a complex structure that can be represented by a cluster-fractal model proposed by Logan.8 This is related to observations of Gunko et al.9 that there is a structural hierarchy in fumed silica. Large clusters of size Li consist of smaller primary silica aggregates of a size La (100-500 nm) that are composed of primary silica particles of a size L0 (5-50 nm). The primary particles forming primary aggregates are connected by strong chemical bonds. Primary aggregates cannot be disintegrated by mechanical forces; they are characterized by average fractal dimension Df0 and contain N0 primary particles. The large agglomerates of size Li are relatively unstable;9 they are connected by adhesion forces and can be disintegrated in the high shear flows. A schematic of the aggregate structure is presented in Figure 7. Application of fractal geometry gives a direct relation between the mass of agglomerate mi and its size Li.
( )( )
mi ) ma
Df
δa θa δp
Df/3
(1)
where θa and δa are the packing factor and shape factor for agglomerates, respectively, δp is the shape factor of aggregates of size La, and Df represents the fractal dimension of agglomerates. Similarly, one can calculate the number of primary particles, N0, in the aggregate of size La
N0 )
Figure 3. Experimental procedure.
Li La
()( ) La L0
Df0
δp θp δp0
Df0/3
(2)
where θp is the packing factor of aggregates, and δp0 is the shape factor of primary particles. The specific surface area of fumed silica is mainly determined by the size of primary particles
3654
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008
Figure 4. Evolution of particle size distribution for 5 wt % Aerosil 200V.
Figure 5. Large agglomerates and submicrometer aggregates of Aerosil 200V.
Figure 6. Submicrometer Aerosil 200V aggregates.
(∼ AL0-1), and only slightly by the type of particle packing9 (A = 4-6). Hence, in what follows, we assume that the packing and shape factors take constant values, independent of agglomerate size and its fractal dimension. On the basis of eq 1, one gets the average porosity of agglomerates, a
()
θaDf/3δa(Df/3)-1 Li a ) 1 δ (Df/3)-1 La p
Df-3
(3) Figure 7. Schematic of a structure of a large agglomerate (cluster).
and the average volume fraction of aggregates in an agglomerate is equal to φ ) 1 - a
()
θaDf/3δa(Df/3)-1 Li φ) δ (Df/3)-1 La p
Df-3
and therefore, it is the local value of volume fraction on the external layer, φ(L), not the average one, that would be of interest. This is derived here using fractal geometry:
(4) φ(L) ) al.10
Equations 3 or 4 were introduced by Tang et to the model equations by Rumpf11 to calculate the tensile strength of agglomerates, σT. However, if erosion is the dominant breakup mechanism, it will occur on the periphery of the agglomerate,
()
Df DfθaDf/3δa(Df/3)-1 Li φ) 3 3 δ (Df/3)-1 La p
Df-3
(5)
and gives as well “external porosity”, which determines the number of bonds on the agglomerate periphery, a(L) ) 1 φ(L).
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008 3655
This value can be now used to calculate the tensile strength for the external layer of particles
1 - a(L) F a(L) La2
σT ) 1.1
(6)
where F represents here the bonding force between aggregates in the agglomerate that is calculated using the classical Derjaguin-Landau-Verwey-Overbeek (DLVO) theory. The bonding forces are attributed to the van der Waals attractive force (FA) and the electrostatic repulsion force (FR):
F ) F A + FR
(7)
The van der Waals forces between primary particles are calculated here assuming that they have equal radii, R0:
FA ) -
HaRo
(8)
24H2
where Ha is the effective Hamaker constant and H represents the surface to surface distance. The electrostatic repulsive forces are given by
FR )
{
16πχRokBT2 exp[zeψd/(2kBT)] - 1 e2z2H
exp[zeψd/(2kBT)] + 1
}
2
exp(-κH) (9)
χ is the static permittivity, kΒ, the Boltzmann constant, e, electron charge, and κ, the Debye-Hu¨ckel parameter
κ)
x
2e2NAVIs χkBT
(10)
Finally using fractal geometry, we estimate the number of aggregates on the agglomerate periphery.
() ( )
Li N ) 2Df La
Df-1
δaθa δp
fitted constants, if any, should depend on this presumed value; in the new model presented in this paper, N is calculated directly from eq 11 that is based on mechanistic interpretation of erosion. In turbulent flows, breakup occurs when the hydrodynamic stresses that are expressed in what follows based on the Kolmogorov theory are larger than the tensile strength σT. The frequency of breakup events is described then by the breakage kernel, Γ, given by
b(L/λ) ) Nδ(L - La) + δ[L - (λDf - NLaDf)1/Df] (12) b(L/λ) represents the daughter or fragment distribution, which defines the probability that a fragment of size L is formed from the breakage of an agglomerate of size λ; then, b(L/λ) dλ represents the number of particles of size L formed by the breakage of a particle of a size between λ and λ + dλ. On the other hand ∫∞0 b(L/λ) dL gives an average number of fragments per one breakage event; in the case of eq 12, this number is equal to N + 1. It should be pointed out that to apply eqs 11 and 12 with the population balance equation, one does not need to introduce new physical constants and model parameters; all of them are present in eq 5 that together with eq 6 is used to model the agglomerate resistance against erosion. This reduces the number of model constants and parameters that need to be fitted. Let us compare this model with the most popular in the subject literature model of binary erosion. In the case of binary erosion, one assumes by definition that N + 1 ) 2 and then all
() ()
} }
(13)
(14)
where λK is the Kolmogorov microscale. Equations 13 and 14 show that the frequency of successful breakage events is proportional to the characteristic frequency of turbulence at a scale that is connected with the agglomerate size. It is interesting that when the proportionality constant Cb is fitted to experimental data it takes small values (Cb , 1; see Section 3.3 of this paper), which means that the breakage frequency is smaller than the characteristic frequency of turbulence. This shows that not all turbulent events are efficient enough to break up agglomerates, which can be explained, at least in part, by the fine scale intermittency of turbulence.12 Another interpretation of eq 11 can be offered by combining it with eqs 13 and 14 and assuming binary breakage. Then, the number N can be attributed to increased frequency of binary events with a well-defined size-dependent breakage kernel Γ′(L,) ) Γ(L,)‚N(L), where Γ(L,) is given by eqs 13 and 14. This leads to the following expressions:
(11)
The number of eroded particles should be proportional to this value; what is most important, eq 11 gives effects of agglomerate size and fractal dimension of the average number of fragments. This results in the following form of fragment distribution:
{
()
1/3 2/3 2/3 C > σT b 2/3 if F Li ν Li for Li > λK ) 1/4 : Γ ) if F2/3Li2/3 e σT 0 3/4
[{
for Li e λK:
Df/3
{
1/2 1/2 Cb if µ > σT ν3/4 ν ν for Li e λK ) 1/4 : Γ ) 1/2 e σT if µ 0 ν
2DfCb
Γ′ )
( ) δaθa δp
]
Df /3
ν1/2LaDf-1
1/2LiDf-1
if µ
0
{[
for Li > λK )
0
1/2
(ν)
1/2
> σT e σT
}
(15)
ν3/4 : 1/4
2DfCb
Γ′ )
(ν)
if µ
( ) δa θa δp
LaDf-1
]
Df /3
1/3LiDf -5/3
if F2/3Li2/3 > σT if F2/3Li2/3 e σT
}
(16)
This means that the considered mechanism yields a popular power-law breakage kernel when binary breakage is assumed, showing the physical meaning of the proportionality constant and the rate exponent; the power-law is often presented13,14 in a form Γ(V) ) Γ0Vη, where V represents particle volume. In such a case, the following values of the rate exponent are predicted for erosion of fractal agglomerates: η ) (Df - 1)/3 for Li < λK (so 0.27 < η < 0.67 for fractal dimensions between 1.8 and 3) and η ) (3Df - 5)/9 for Li > λK (so 0.044 < η < 0.44).
3656
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008
3.2. Suspension Viscosity. In the subject literature, one can find many empirical or semiempirical relations that allow the prediction of effective viscosity of suspensions. In this work, we focus on dense suspensions of nanoaggregates; hence, we consider and apply in computations the model15 that is valid for both high and low volume fraction of particles and which covers the whole range of shear rate values in the flow. For the concentrated suspension, the Pe´clet number can be expressed as
D∞B γ˘ a2 kT Pe ) , where DB(φ) ) , and D∞B ) 6πµliquida DB(φ) χ(φ) (17) where χ(φ) represents effect of particle concentration on Brownian self-diffusion. Buyevich and Kapbsov15 considered both the flows with strong effects of Brownian diffusion and negligible shear effects (Pe f 0) and the flows with negligible effects of Brownian diffusion and high effects of shear (Pe f ∞), to formulate finally the relation valid also for intermediate Pe values. Applying the following formula for relative viscosity, M(φ)
M(φ) )
µ µliquid
(18)
Buyevich and Kapbsov15 derived the following relation for the narrow concentration range that adjoins φ*, φ* - φ , φ,
(
M(φ) ) -C ln 1 -
)
φ1/3 (φ*)1/3
(19)
where φ* represents the close-packed limit for φ,
C)
zφ* 2
(20)
and z is an effective “coordination number” of a sphere in the close-packed state. Using eq 19 together with the empirical formula by Frankel and Akrivos16 that is valid for low and moderate concentrations,
µ ) µliquid(1 - φ)-5/2 for φ e φ* - φ
(21)
Buyevich and Kapbsov15 developed expressions covering a whole range of concentrations: for high shear rates, Pe f ∞,
{ [ ( ) ] ∑( ) }
M∞(φ) ) (1 - φ)-5/2 - C ln 1 -
φ
6
1/3
+
φ*
1 φ
j/3
with C ≈ 2 and, for low shear rates, Pe f 0
M0(φ) ) (1 - φ)
[( )
+ 1.3 1 -
φ
2
-2
φ*
-
( )] φ
(1 + j) ∑ φ* j)0
j
(23) For intermediate Pe values, the following equation can be used for the whole φ range:
M)
M0 + M∞Pe 1 + Pe
n
φeff ) φeff,a
(24)
Chen et al.5 have shown that the rheology of Aerosil suspension
wmiNi(3-D )/D ∑ i)1 f
f
(25)
where φeff,a ) φ0N0(3-Df0)/Df0 is the effective volume fraction of primary aggregates and wmi represents the mass fractions of each of n classes applied in QMOM. For just suspended and not deagglomerated Aerosil 200V, the model predicts for 20 wt % M0 ) 25 000 and M∞ ) 25.9, which represents the range of variation of relative suspension viscosity. One can see that, at a very low shear rate, when Pe f 0, the viscosity can be very high, but it decreases significantly when the suspension flows. For 15 wt %, the relative suspension viscosity is between 6.6 and 18, and for 5 wt %, it is between 1.56 and 1.63 (so, waterlike and Newtonian in this last case). In calculations, the close packing concentration φ* ) 0.695 was applied for fractal agglomerates that is higher than the value 0.67 used for compact, spherical particles. 3.3. Population Balance for Agglomerate Erosion. The size ratio between agglomerates to be disintegrated and primary aggregates is about 3 orders of magnitude, which means that the mass ratio for a high value of the fractal dimensions Df ≈ 2.8 that is observed in the considered process is about 108. As a result, when the population balance equations are solved with the QMOM,2 there are very large differences in values of abscissas representing particle mass, which creates problems with convergence during numerical computations. Hence, in this work, in order to follow the evolution of agglomerate population, the population balance equation is used in a form based on a number density function defined in terms of particle characteristic length L.
∂n(L,x b,t) ∂t
3
+
∑i
∂[upi(L,x b,t)n(x b,t)]
+
∂[G(L,x b,t)n(L,x b,t)]
∂xi
)
∂L b,t) - Db(L,x b,t) (26) Bb(L,x
j)1 j φ*
(22)
-5/2
can be described using the model of Krieger and Dougherty,17 similar to eq 21. The model applied in the present paper can be considered as generalization of models by Frankel and Akrivos16 and Krieger and Dougherty;17 it is worth noticing that the variation in time and space of effective volume fraction φ of agglomerates that result from variation of the agglomerate structure and size distribution, allows the time and space dependent rheology of agglomerated suspensions to be represented. In this paper, the effective volume fraction of Aerosil agglomerates is recalculated from the agglomerate size distribution based on the QMOM
where the birth and death rate functions for particle breakage are given by
b,t) ) Bb(L,x
∫L∞Γ(λ)b(L/λ)n(λ,xb,t) dλ
Db(L,x b,t) ) Γ(L)n(L,t)
(27) (28)
The growth term G(L,x b,t) reflects the variation of agglomerate size resulting from variation of its fractal dimension. Consider an agglomerate of mass m and size L; when the agglomerate mass is conserved and the agglomerate size varies only due to variation of its fractal dimension Df, then differentiation of an expression L = La(m/ma)1/Df (see eq 1) gives the following:
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008 3657
G)
()
dL L 1 dDf L )ln dt Df dt La
(29)
It should be noted that a complete form of eq 29 reads
G)-
[ ( ) ( )]
1 δaθa L 1 dDf + ln L ln Df dt La 3 δp
but because
δa θa δaθa L >>> and =1 La δp δp a reduced form given by (29) is used in this paper. It should be also pointed out that when agglomerate mass is used as an internal coordinate, then such additional growth term will not be necessary because for a given agglomerate mass its size can be recalculated directly from eq 1. Due to action of stresses, as well as aggregation and breakage events, agglomerates change their structure which is expressed here by the aggregate fractal dimension. This phenomenon is known in the subject literature as restructuring. Variation of the fractal dimension can be followed using eq 30, where evolution of Df is interpreted as a superposition of relaxation effects related to the three limiting values, characteristic for each involved phenomenon
Figure 8. Evolution of cumulative number distribution of agglomerate sizesmodel predictions for erosion of agglomerates using direct integration of the population balance; tend ) 0.15 s.
dDf ) Cshearγ˘ (Df,shear(max) - Df) + Caggβm0(Df,agg(min) dt Df ) + CbreakΓ(Df,break(min) - Df) (30) where m0 denotes particle concentration (equal to the zero moment of the distribution) and β and Γ are calculated by summing up participations from all QMOM nodes using related weights and abscissas. The constants Cshear, Cagg, and Cbreak at the present stage of model development should be fitted to experimental data. Results of application of the complete eq 30 are presented in ref 18. In the present paper, we only consider effects of breakage (Cshear ) 0, Cagg ) 0). When agglomerates break up and the internal porosity is conserved, then the fractal dimension of smaller fragments has to be smaller than that of the primary agglomerate. To solve the population balance equation, the moment transformation is applied to eq 26, which leads to the system of equations with unclosed birth, death, and growth terms. To close this system of equations, the quadrature method of moments (QMOM) is applied. The method is based on a following quadrature approximation: n
wiδ(L - Li) ∑ i)1
n(L) =
(31)
which results in the following expressions for moments
mk )
∑i wiLik
(32)
and the following closed form of the moment balance equation:
b,t) ∂mk(x ∂t
+ ui
∂mk(x b,t)
)
∂xi
∂ ∂xi
n
[
]
∂mk(x b,t)
DT
∂xi
n
1
n
+
Γib(k) ∑ i wi i)1
() Li
Li Γiwi - GD k∑Likwi ln ∑ D L i)1 i)1 k
f
f
a
(33)
Figure 9. Comparison of model predictions for direct integration of the population balance equation (solid lines) and QMOM (broken lines). (a) Evolution of the average size of large agglomerates. (b) Evolution of the average size of the complete population. (c) Evolution of effective viscosity; tend ) 0.15 s.
where b(k) is given by b(k) ) NLak + (LiDf - NLaDf)k/Df for i i erosion and GDf ) dDf/dt ) CbreakΓ(Df,break(min) - Df). Results of simulations presented in this and the next section of this paper are obtained for following parameters characterizing Aerosil structure and the process parameters: La ) 250 nm, L0 ) 22 nm, Df0 ) 2.5, initial value of Df equal to 2.85, θp ) 1, θa ) 1, δp ) π/6, δa ) π/6, ma ) 5.689 × 10-19 kg, F0 ) 2200 kg m-3, and with agglomerate strength characterized by Ha ) 5 × 10-21 J, H ) 10-9 m for breakage constant Cb ) 2 × 10-3. A unimodal distribution of agglomerate size at the system inlet was applied. Predictions of agglomerate dispersion in a simple onedimensional plug flow system for 15 wt % of Aerosil with the
3658
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008
population balance equation. The average sizes are defined by Lij ) mi/mj and plotted for both: large agglomerates only (Figure 9a) and complete population (Figure 9b). One can see that the difference between differently defined average sizes for large agglomerates initially increases and afterward decreases with the progress of deagglomeration; this reflects the variation of the size distribution of large agglomerates that is also seen in Figure 4. In the case of complete population, Figure 9b, due to a very large number of small fragments, the average size L10 shows in fact only the size of small aggregates, which agrees with results presented in Figure 8. One can also see that deagglomeration results in a decrease of suspension viscosity (Figure 9c). 4. Results and Discussion
Figure 10. Geometry of the rotor-stator mixer: the jacket of the rotorstator device and details of the rotor-stator geometry represented by the 3D grid.
rate of energy dissipation ) 105 m2 s-3 are shown in Figures 8 and 9. Figure 8 shows evolution of the cumulative number distribution of the agglomerate size obtained by direct integration of the population balance eq 26; the bimodality of distribution is very clearly observed. One can also see that the number of small aggregates is so large that it dominates the number distribution after a short deagglomeration time. Figure 9 shows that application of QMOM gives predictions that are almost identical with results obtained by direct integration of the
To interpret and predict experimental results, the model was combined with the k- model of Fluent, and the multiple reference frame model was used to simulate the processes in the high-shear rotor-stator mixer. The jacket of the rotor-stator device and the 3D computational grid illustrating the rotorstator geometry are shown in Figure 10. As it is not possible to follow the process during 14 h for about 300 tank turnovers, the simulations for model-experiment comparison are carried out using the two-dimensional grid. Typical CFD results characterizing details of flow and deagglomeration, as well as comparison of 3D and 2D predictions are shown in Figures 11-13. The figures show that the results for 3D and 2D simulations do not differ much and that the zone of high shear is localized in the region where both rotors operate and in jets emerging from the holes of the disintegration screen. What is important in both cases is that the effect of deagglomeration is
Figure 11. Velocity distribution in a rotor-stator mixer for Q ) 0.60 × 10-3 kg s-1; comparison of 3D and 2D simulations.
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008 3659
Figure 12. Distribution of the rate of energy dissipation, , in a rotor-stator mixer for Q ) 0.60 × 10-3 kg s-1; comparison of 3D and 2D simulations.
Figure 13. Distribution of agglomerate concentration, Na, in a rotor-stator mixer for Q ) 0.60 × 10-3 kg s-1; comparison of 3D and 2D simulations.
3660
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008
Figure 15. Effect of the rotor speed and flow rate on mean particle size L32 for the first passage through the rotor-stator mixer.
Figure 16. Changes of effective viscosity during suspension circulation through the stirred tank and the rotor-stator mixer for different values of Aerosil concentration.
Figure 14. Distribution of agglomerate size, L30, first moment, m1, and effective viscosity, µeff, in a rotor-stator mixer for Q ) 0.60 × 10-3 kg s-1, 20 wt %, 9000 rpm, after the 239th tank turnover; 2D simulations.
similar, which can be seen by comparing the concentration of agglomerates in both cases at the mixer outlet as presented in Figure 13. The maximum difference between 2D and 3D simulations is in increase of the moment of order zero for the first passage through the rotor-stator mixer; this increase is 2.7% higher in the case of 2D simulation. Figure 14 shows typical predictions of deagglomeration effects for one passage through the rotor-stator system. Clearly, one passage through the mixer results in rather small effect of deagglomeration, which explains why in practical applications the suspension circulates many times through the mixing head. Moreover, increasing the flow rate through the rotor-stator mixer, Q, decreases deagglomeration effects,19 which means that reducing the process time simply by increasing the flow rate through the mixer is hardly possible. This effect is presented in Figure 15 for the first passage through the rotor-stator device. Figure 16 shows variation of the effective viscosity at the outlet from the
rotor-stator mixer during the process of deagglomeration; decrease of viscosity is clearly seen at higher concentration (15 and 20 wt %) and is hardly observed for small suspension concentration (5 wt %). Figures 17-19 show comparisons of experimental results and model predictions for 5, 15, and 20 wt %, respectively. The progress of deagglomeration is rather well-predicted, and the agreement between model predictions and experimental results is good up to the moment of order 3. The L30 average size is better predicted than the L32 size, which results most probably from the fact that simulated deagglomeration of very large agglomerates is slower and that of large agglomerates is faster than observed in experiments; see Figure 20, where this effect is clearly observed. Most probably in reality the forces bonding large fragments of very large agglomerates are weaker than predicted by the model. Then, erosion of large agglomerates would be predated by rupture of very large agglomerates. Figure 20 shows evolution of the agglomerate size distribution during suspension circulation through the stirred tank and the rotor-stator mixer, and one can see that the evolution of distribution is well predicted including distribution bimodality. The fractal dimensions Df decrease from 2.85 to about 2.75 when the agglomerate size reaches its final value. Experimental data for direct comparison of this prediction with experimental data obtained using the Silverson rotor-stator are not available but similar decrease of Df was observed in other equipment applied in Proform project.
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008 3661
Figure 17. Changes of moments of agglomerate size distribution (in relation to initial value mk0) and average agglomerate sizes during suspension circulation through the stirred tank and the rotor-stator mixer for an Aerosil concentration of 5 wt %. Experimental results are marked by points, and lines represent results of simulations.
Figure 18. Changes of moments of agglomerate size distribution (in relation to initial value mk0) and average agglomerate sizes during suspension circulation through the stirred tank and the rotor-stator mixer for an Aerosil concentration of 15 wt %. Experimental results are marked by points, and lines represent results of simulations.
In section 3.1, it has been shown that the model of erosive deagglomeration proposed in this work predicts as well a powerlaw size-dependent breakup kernel when binary erosion is assumed. Let us compare now predictions of this new erosion model with predictions of modeling based on a rupture mechanism. Rupture refers to the breakage of an agglomerate into fragments of comparable size. It can be interpreted using two simple binary rupture mechanisms: the first one results in two equisized fragments, and the second one assumes a uniform mass fragment distribution with minimum fragment size La. Because the first one does not lead to bimodal distribution, let us consider in what follows to be the second one. The fragment distribution reads then:
breakup similarly as observed in experiments, whereas its value drifts very slowly toward La in the case of uniform breakage. Also, the experimentally determined behavior of L30 (Figures 17-19) is better predicted by the erosion model. Finally, it is known20-22 that the main qualitative difference between the mechanisms of rupture and erosion is energy input. It is low for erosion and high for rupture. Figures 15 and 17-20 show that the effect of breakage per one passage is small, which indicates relatively small hydrodynamic stresses and thus small fragmentation numbers. For compact agglomerates,20-22 the rate of erosion was shown to be proportional to the surface area of the agglomerate, which agrees with assumptions leading to eqs 11 and 12. The kinetics of erosion presented in the subject literature (see eq 11 in ref 20) valid for short times leads to linear decrease of parent agglomerates,
b(L/λ) )
2DfLDf-1 λDf - LaDf
(34)
for La < L < λ and b(L/λ) ) 0 otherwise. Computations are performed for a constant value of the rate of energy dissipation for the process carried out in the one-dimensional system. The results are presented in Figure 21 using dimensionless time θ ) t(/ν)1/2 with θend ) 6 and tend ) 6(ν/)1/2. Comparison of Figure 9b with Figure 21 shows that in both cases the resulting distribution has a bimodal character; however, the uniform breakage predicts different character of bimodality than presented in Figure 4. Note that the L10 size that is dominated by smallest fragments is almost constant in the case of erosion
L(t) ) Linitial(1 - const‚γ˘ t)
(35)
which agrees well with predictions of the new model presented in Figure 9a. 5. Conclusions The breakage of silica nanoparticle agglomerates can be described by an erosion mechanism by which small fragments, approximately of the size of the primary aggregates, are chipped off from larger agglomerates. This is evidenced by the formation
3662
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008
Figure 19. Changes of moments of agglomerate size distribution (in relation to initial value mk0) and average agglomerate sizes during suspension circulation through the stirred tank and the rotor-stator mixer for an Aerosil concentration of 20 wt %. Experimental results are marked by points, and lines represent results of simulations.
of a bimodal particle size distribution during the deagglomeration process and the images obtained using electron microscopy. A new model for erosive deagglomeration of Aerosil 200V agglomerates has been proposed to interpret the agglomerate disintegration processes. Also, an effective method of application of population balances to describe the evolution of fractal aggregates and effective suspension viscosity has been presented. It is shown that a model can be used to interpret processes in the rotor-stator mixer, a device commonly applied in industry. Comparison of model predictions with experimental results shows that the model predicts well the trends observed in experiments. Acknowledgment This study was carried out within the project PROFORM (“Transforming Nanoparticles into Sustainable Consumer Products Through Advanced Product and Process Formulation” EC Reference NMP4-CT-2004-505645) which was partially funded by the 6th Framework Programme of EC. The contents of this paper reflect only the authors’ view. The authors gratefully acknowledge the useful discussions held with other partners of the Consortium: BHR Group Limited; Karlsruhe University, Institute of Food Process Engineering; Bayer Technology Services GmbH; University of Loughborough; Unilever UK Central Resources Limited; Birmingham University School of Engineering; Poznan University of Technology, Institute of
Figure 20. Changes of the agglomerate size distribution during suspension circulation through the stirred tank and the rotor-stator mixer for an Aerosil concentration of 20 wt %. Experimental results are marked by solid lines, and results of simulations are given by broken lines.
Figure 21. Variation of the mean particle size Lij for the uniform breakage process.
Chemical Technology and Engineering; Rockfield Software Limited; C3M d.o.o. Ljubljana. Nomenclature Bb(L) ) birth function for breakage, m-4 s-1 b(L/λ) fragment distribution, m-1 Cb ) proportionality constant in breakage kernel C ) proportionality constant in eq 19 defined in eq 20 Db(L) ) death function for breakage, m-4 s-1 DT ) turbulent diffusivity of particles, m2 s-1 Df ) fractal dimension of agglomerates Df0 ) fractal dimension of primary aggregates
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008 3663
e ) electron charge, C; e ) 1.6022 × 10-19C F ) particle interaction force, N FA ) attractive force, N FR ) repulsive force, N G ) agglomerate growth rate due to variation of fractal dimension, m s-1 GDf ) rate of change of fractal dimension, s-1 H ) surface to surface distance, m Ha ) Hamaker constant, J Is ) ionic strength, kmol m-3 kB Boltzmann constant, J K-1; kB ) 1.38 × 10-23 J K-1 L ) agglomerate size, m La ) primary aggregate size, m Li ) agglomerate size, abscissas of the quadrature approximation, m Lij ) average particle size defined with moments of order i and j, m L0 ) primary particle size, m M ) relative viscosity m, mi ) particle mass, kg ma ) mass of aggregate, kg N ) number of eroded fragments NAV ) Avogadro number, mol-1; NAV ) 6.0221 × 1023 mol-1 Ni ) number of primary aggregates in an agglomerate N0 ) number of primary particles in an aggregate n ) order of quadrature approximation Pe ) Pe´clet number T ) temperature, K upi ) particle velocity component, m s-1 V ) agglomerate volume, m3 wi ) weights of the quadrature approximation, m-3 wmi ) mass fraction of class i, kg kg-1 z ) valence of ions, effective coordination number Greek Letters Γ(L) ) breakage kernel, s-1 γ˘ ) rate of shear, s-1 δ(x) ) Dirac delta function δa ) shape factor for agglomerates δp ) shape factor for aggregates δp0 ) shape factor for primary particles ) rate of energy dissipation per unit mass, m2 s-3 a ) porosity of agglomerates 0 ) porosity of primary aggregates θp ) packing factor for primary particles in aggregates θa ) packing factor for aggregates in the agglomerate κ ) Debye-Hu¨ckel constant, m-1 λ ) agglomerate size, m λK ) Kolmogorov microscale, m µ ) dynamic viscosity, kg m-1 s-1 ν ) kinematic viscosity, m2 s-1 F0 ) density of primary particles, kg m-3 σT ) tensile strength of agglomerates, Pa φ ) volume fraction φ0 ) volume fraction of primary particles φeff,a ) effective volume fraction of primary aggregates φeff ) effective volume fraction
χ ) dielectrict constant, C V-1 m-1 ψd ) surface potential, V Literature Cited (1) Parfitt, G. D.; Barnes, H. A. The dispersion of fine powders in liquid media. In Mixing in the Process Industries, second ed.; ButterworthHeinemann: Oxford, 1992. (2) McGraw, R. Description of aerosol dynamics by the quadrature method of moments. Aerosol Sci. Technol. 1997, 27, 255. (3) Fine Particles Technical Bulletin Number 11: Basic Characteristics of Aerosil Fumed Silica; Degussa AG: Essen, Germany, 2006. (4) Padron, G.; Eagles, W.P.; O ¨ zcan-Tas¸ kin N. G.; McLeod, G.; Xie, L. Effect of particle properties on the break up of nanoparticle clusters using an in-line rotor-stator. 1st International Conference on Industrial Processes for Nano and Micro Products, London, UK April 2-3, 2007. (5) Chen, S.; Øye, G.; Sjo¨blom, J. Rheological Properties of Aqueous Silica Particle Suspensions. J. Dispersion Sci. Technol. 2005, 26, 495. (6) Ottino, J. M.; DeRoussel, P.; Hansel, S.; Khakar, D. V. Mixing and dispersion of viscous liquids and powdered solids. AdV. Chem. Eng. 2000, 25, 105. (7) Wengeler, R.; Nirschl, H. Turbulent stress induced dispersion and fragmentation of nanoscale agglomerates. J. Colloid Interface Sci. 2006, 306, 262. (8) Logan, B. E. EnViromental Transport Processes, Wiley: New York, 1999. (9) Gunko, V. M.; Zarko, V. I.; Leboda, R.; Chibowski, E. Aqueous suspensions of fumed oxides: particle size distribution and zeta potential, AdV. Colloidal Interface Sci. 2001, 91, 1. (10) Tang, S.; Ma, Y.; Shiu, C. Modeling of mechanical strength of fractal aggregates. Colloids Surf. A: Phys. Eng. Aspects 2001, 180, 7. (11) Rumpf, H. Agglomeration; Knepper, W. A., Ed.; Interscience: New York, 1962. (12) Baldyga, J.; Podgorska, W. Drop break-up in intermittent turbulence. Maximum stable and transient sizes of drops. Can. J. of Chem. Eng. 1998, 76, 456. (13) Ziff, R. M. New solutions for the fragmentation equation J. Phys. A: Math. Gen. 1991, 24, 2821. (14) Diemer, R. B.; Olson, J. H. A moment methodology for for coagulation and breakage problems: Part I -analytical solution of the steady-state population balance. Chem. Eng. Sci. 2002, 57, 2193. (15) Buyevich Yu, A.; Kapbsov, S. K. Segregation of a fine suspension in channel flow. J. Non-Newt. Fluid. Mech. 1999, 86, 157. (16) Frankel, N. A.; Akrivos, A. On the viscosity of a concentrated suspension of solid spheres. Chem. Eng. Sci. 1967, 22, 847. (17) Krieger, I. M.; Dougherty, T. J. A mechanism for non-Newtonian flow in suspensions of rigid spheres. Trans. Soc. Rheol. 1959, 3, 137. (18) Bałdyga, J.; Krasinski, A.; Orciuch, W. Rheological effects in concentrated aggregated suspensions - application of population balance modeling. In Proceedings of the 12th European Conference on Mixing, Bologna, Italy, June 25-27, 2006. (19) Bałdyga, J.; Orciuch, W.; Makowski, L -.; Malski-Brodzicki, M.; Malik, K. Break up of nano-particle cluster in high-shear devices. Chem. Eng. Process. 2007, 46 (9), 851. (20) Hansen, S.; Khakhar, D. V.; Ottino, J. M. Dispersion of solids in inhomogeneous viscous flows. Chem. Eng. Sci. 1998, 53, 1803. (21) Rwei, S. P.; Manas-Zloczower, I.; Feke, D. L. Observation of carbon black agglomerate dispersion in simple shear flows. Polym. Eng. Sci. 1990, 30, 701. (22) Rwei, S. P.; Manas-Zloczower, I.; Feke, D. L. Characterization of agglomerate dispersion by erosion in simple shear flows. Polym. Eng. Sci. 1991, 31, 558.
ReceiVed for reView July 1, 2007 ReVised manuscript receiVed November 2, 2007 Accepted November 2, 2007 IE070899U