Article pubs.acs.org/IECR
Application of Population Balance Model in the Simulation of Slurry Bubble Column Lijia Xu, Zihong Xia, Xiaofeng Guo, and Caixia Chen* Key Laboratory of Coal Gasification and Energy Chemical Engineering of Ministry of Education, East China University of Science and Technology, Shanghai, 200237, China ABSTRACT: Numerical simulations of a slurry bubble column with particle loadings up to 40 vol % and superficial gas velocities up to 0.26 m/s have been performed using the Euler−Euler approach with the renormalized group (RNG) k−ε turbulence model. Special attention was paid to the bubble size distributions, interphase closure models, and liquid−solid drag forces. A modified Luo−Lehr population balance model was used to simulate the changes of mean bubble size and the overall gas holdup as functions of particle loadings. A pseudo gas−slurry closure model was proposed to overcome the drawbacks of existing interphase momentum exchange closure models. The newly constructed model was different from the traditional gas−slurry closure model in that the hydrodynamics of three phases are solved by three sets of momentum equations describing gas, liquid, and solid phases, respectively. Three different drag force models, i.e., the Schiller−Naumann model, the Wen−Yu model, and the energy-minimization multiscale (EMMS) model, were used for the computation of the liquid−solid interaction force and the simulated particle settlings were compared with experiments. A series of numerical simulations were performed. The following conclusions were drawn from the simulation results: (1) The Luo−Lehr population balance model could simulate the mean bubble size changes as a function of particle loadings. The gas holdup and its variation trend predicted by the modified Luo−Lehr population balance model were in agreement with the experimental data by setting a coalescence coefficient C0 = 0.4. (2) A pseudo gas−slurry closure model that took the influence of solid concentration into account was recommended for gas−liquid− solid three-phase simulation when the particle size was much smaller than the bubble size. The simulated average gas holdups showed good agreement with the experimental data. (3) The axial distribution of solid concentration simulated by the EMMS drag model was much closer to experimental results than other drag models for the simulation of slurry bubble columns.
1. INTRODUCTION Slurry bubble column reactors are widely used in petrochemical and coal liquefaction processes due to their advantages of high capacity, easy operation, and good heat and mass transfer characteristics. Traditionally, the design and scale-up of these reactors rely on a large number of empirical correlations obtained from different scale experiments. Those correlations include estimations of the gas holdup, the interphase heat transfer coefficient, the axial solid concentration, and other hydrodynamic parameters.1 However, applications of the empirical correlations are limited by the experimental conditions, such as the superficial gas velocities, liquid properties, particle sizes and concentrations, among others. Computational fluid dynamics (CFD) has provided the state-of-the-art method for the description of fluid dynamics in gas−liquid flows (e.g., Chen and Fan2 and Buwa et al.3). Nevertheless, the application of the CFD approach to the predictions of gas−liquid−solid fluid dynamics is still immature due to the inherent complexities including the description of bubble size distributions, the uncertainty of the interphase closure models, and the choice of liquid−solid drag model. Over the years, the Eulerian multifluid model incorporated with the population balance model (PBM) has been applied in the simulations of air−water bubble columns operated in the churn-turbulence regime, and the predicted distributions of the axial liquid velocity and the gas holdup agree with the experiments.4,5 The application of the PBM technique to the simulation of a gas−liquid−solid three-phase system needs considerable elaboration of modeling the influence of solids on © 2014 American Chemical Society
hydrodynamics, including the breakup and coalescence of bubbles, and the interaction forces between liquid−solids and gas−solids. Chen and Fan6 suggested that the effect of bubble− particle collisions on the bubble breakup could be neglected when the particle size is much smaller than the bubbles. Hooshyar et al.7 studied the dynamics of a single rising bubble in liquid−solid suspensions, and found that the small particles (78 and 587 μm) remained on the liquid flow stream around the bubble. Their experimental findings demonstrated that the liquid−particle mixture behaved more like a single pseudoclear liquid than two discrete media. In this regard, the gas−liquid− solid three-phase problem can be simplified as a gas−slurry two-phase problem, and the influence of small solid particles on the bubble breakup and coalescence are taken into account through a change of physical properties of the slurry instead of direct bubble−particle collisions. The viscosity of slurry was reported to increase with particle concentrations, and a higher particle concentration would usually result in an increase of bubble size and a decrease of gas holdup.8 Hence, the increases of bubble size can be interpreted as a result of increases in the viscosity of slurry in a gas−liquid−solid bubble column. Following this scenario, the gas−liquid−solid three-phase problem could be simplified as a gas−slurry two-phase problem, as long Received: Revised: Accepted: Published: 4922
October 14, 2013 February 22, 2014 February 27, 2014 February 27, 2014 dx.doi.org/10.1021/ie403453h | Ind. Eng. Chem. Res. 2014, 53, 4922−4930
Industrial & Engineering Chemistry Research
Article
Figure 1. Closure law for the gas−liquid−solid interphase momentum exchange.
as the particles mix with the liquid perfectly or the particle size is small.9 The gas−slurry two-phase model is referred to as “closure A” in the present paper, and is depicted in Figure 1a. Such a simplification avoided the explicit closures of the gas− solid and liquid−solid interphase momentum exchanges and provided reasonable results in certain cases.9,10 However, when the particle settling is obvious in some gas− liquid−solid flows, the liquid and solids need to be treated in two distinct phases, and closure models which describe the interaphase momentum exchanges between gas−liquid−solid three phases are required. In most published literature, only the liquid−gas and liquid−solid drag forces were considered,11,12 and the gas−solid interaction force was neglected, referred to as “closure B” as shown in Figure 1b. In recent years, an increasing number of researchers suggested that the gas−solid interaction did exist in the gas−liquid−solid flows, based on the fact that the particles in the vicinity of a bubble tend to follow the bubble. Padial et al.,13 Schallenberg et al.,14 Panneerselvam et al.,15 and others directly used a general form of the drag force model (such as the Wen−Yu and Schiller−Naumann drag models) to calculate the gas−solid interaction forces, and treated the gas phase as a continuous phase and the particles as a dispersed one, respectively. This treatment made a closed gas−liquid−solid−gas loop of interface exchanges, and is referred to as “closure C” as shown in Figure 1c. Here, the closure C model involves three major uncertainties: first, in the majority of slurry bubble column reactors, the gas phase is presented as bubble foam and the particles are more likely to be presented in the liquid. That is, the particles can barely penetrate the gas−liquid interface and dance around inside a gas bubble. A direct, explicit gas−solid interfacial drag force would overpredict the gas−solid interactions. Second, when bubbles rise in a slurry, a negative pressure generated in the bubble wake will induce the liquid to move, and the movement of the liquid drives the particles to rise faster than those that are
not in the wake. Rather than a direct action, movements of gas bubbles influence the solids motion in an indirect way. Third, the liquid is simultaneously susceptible to the upward force from bubbles and the downward force from solids, and the solids are driven by the upward forces from bubbles and liquid, simultaneously. This effect may cause a higher velocity of particles than the liquid phase, or a decrease of the relative velocity between liquid and solids. Therefore, the closure C model may underestimate the particle settling. For the simulations of particle settling, a reasonable closure model of liquid−solid interaction forces plays a critical role. The Wen−Yu and Gidaspow drag models were recommended for dilute and dense gas−solid or liquid−solid flow systems, respectively. Visuri et al.16 simulated the liquid−solid (particle size 1.5 mm) flows with different drag models, and found that the bed height obtained with the Wen−Yu model was higher than that obtained in the experiments. On the other hand, the Gibilaro model was reported to be appropriate for flows with high particle Reynolds numbers, and it predicted the particle concentrations in good agreement with experiments.16 Yang et al.17 simulated a gas−solid flow in a circulating fluidized bed, and found that the Wen−Yu/Ergun drag model produced a more homogeneous structure. In addition, they found that the Wen−Yu and Gidaspow drag models overestimated the gas− solid drag force in a dilute gas−solid system. In the present paper, the gas−liquid−solid flows were simulated by the Euler−Euler approach implemented in ANSYS Fluent software. The experimental setup of Gandhi et al.18 was used as a benchmark case. The PBM equations were solved in combination with a mixture manifold of the renormalized group (RNG) k−ε model. A pseudo gas−slurry closure model referred to as the “closure D” model as indicated in Figure 1d was proposed. The newly constructed model was different from the existing closure A model in that the hydrodynamics of three phases were solved by three sets of momentum equations 4923
dx.doi.org/10.1021/ie403453h | Ind. Eng. Chem. Res. 2014, 53, 4922−4930
Industrial & Engineering Chemistry Research
Article
describing gas, liquid, and solid phases, respectively. Using the closure D model, the particle settlings with various particle loadings were simulated using different liquid−solid drag models. The simulated results were compared with experimental data of literature.
Cd,distorted
fα = αslu1.5
2. DESCRIPTION OF COMPUTATIONAL METHODS 2.1. Governing Equations of the Pseudo Gas−Slurry Closure Model. For the Eulerian multifluid model, bubbles and particles are treated as dispersed phases. The continuity and momentum of gas, liquid, and solids are solved by using the phase-coupled SIMPLE algorithm (PC-SIMPLE). The conservation equations of each phase are presented as follows: ∂ρi αi ∂t
+ ∇·(αiui) = 0
∂(αiρi ui) ∂t
Cd,cap =
(10)
The liquid−solid drag force used in closures B, C, and D are computed using the Schiller−Naumann correlation22 (eqs 11 and 12), which is also used to calculate the gas−solid interaction force in closure C.
(1)
K sl,Schiller − Naumann =
(2)
3 αlαsρl |ul − us| Cd 4 ds
(11)
⎧ 24 0.687 ) Res ≤ 1000 ⎪ (1 + 0.15Res Re Cd = ⎨ s ⎪ Res > 1000 ⎩ 0.44 23
(12) 17
For comparison, the Wen−Yu and the EMMS liquid−solid drag models are also used in the simulations of the closure D model. 3 αlαsρl |ul − us| −2.65 αl Cd 4 ds
(13)
3 αlαsρl |ul − us| −2.65 Cd HD αl 4 ds
(14)
K sl,Wen − Yu =
(4)
K sl,EMMS =
HD = a(Res + b)c
(15)
Here, the definitions of a, b, and c are referred to Hong’s work.24 It should be noted that the existing EMMS model based on superficial gas velocities and solid flux is only used to simulate a gas−solid fluidized bed. It is not easy to build a EMMS model for liquid−solid interaction in a slurry bubble column. In the present paper, the Hong-EMMS model is directly used to calculate the liquid−solid interaction force. The drag force coefficient Cd of the Wen−Yu (eq 13) and EMMS (eq 14) models can be calculated as ⎧ 24 (1 + 0.15(αlRes)0.687 ) Res < 1000 ⎪ α Re ⎨ Cd = l s ⎪ 0.44 Res > 1000 ⎩
(16)
4
μslu = μ l (1 + 2.5αs + 10.05αs 2 + 0.00273e16.6αs)
(5)
ρslu = ρl αl + ρα s s
(6)
In our previous research work, it was found that reasonable gas holdup distributions could be obtained only by taking the lift force into account, and the Tomiyama lift force25 was used: Flift = −C liftρslu αg(ul − ug) × (∇ × ul)
When the population balance model is incorporated in the simulation, the Sauter mean diameters of bubbles change in a wide range. The Ishii−Zuber model21 (nonspherical bubble drag model) is used to calculate the liquid−gas drag coefficient: 24 (1 + 0.1Re 0.75) Re
(9)
= min(Cd,distorted , Cd,cap)
In the closure D model, the density and viscosity of the liquid−solid slurry are used, which poses the major difference compared to the closure B and C models, where the pure liquid properties are used in the computations. The exchange coefficients Kl,s of liquid−solid drag forces are presented in eqs 11−15. Note that the closure D model is different from closure A in that the hydrodynamics of three phases are solved by three sets of momentum equations describing gas, liquid, and solid phases, respectively. The liquid and solids are treated as a single slurry phase in the closure A model, and the viscosity and density of the slurry phase are functions of solid concentrations; the surface tension of slurry remains unchanged when the solid weight loading increases.19 The viscosity of the liquid−solid slurry is calculated using the Thomas semitheoretical correlation,20 and the density is the volume-weighted average density:
Cd,spherical =
8 2 αslu 3
= Cd,spherical ; else Cd
where αi is the volume fraction of the ith phase. The sum of the volume fractions of all phases equals 1. ρi, τi, ui, and Fi,j represent the density, viscous stress tensor, velocity, and interface momentum exchange term, respectively. The liquid−gas (Fl,g) and liquid−solid (Fl,s) drag models used in the pseudo gas−slurry closure model (closure D) are expressed as 3 αlαg Fl,g = Cslu,g ρslu |ul − ug|(ul − ug) 4 dg (3)
Fl,s = Kl,s(ul − us)
(8)
If Cd,spherical > Cd,distorted , Cd
+ ∇·(αiρi uiui)
= −αi·∇p − ∇(αiτi) + Fi , j + αiρi g
2 g (ρslu − ρg ) ⎛ 1 + 17.67f 6/7 ⎞ 2 α ⎜ ⎟ = d ⎜ ⎟ , 18.67fα 3 σ ⎝ ⎠
(17)
where C lift
(7) 4924
⎧ min[0.288 tanh(0.121Re), f (Eo)] Eo < 4 ⎪ = ⎨ f (Eo) 4 < Eo < 10 ⎪ ⎩− 0.27 10 < Eo
(18)
dx.doi.org/10.1021/ie403453h | Ind. Eng. Chem. Res. 2014, 53, 4922−4930
Industrial & Engineering Chemistry Research
Article
Table 1. Models for Bubble Breakup Kernelsa authors
breakage rate kernel 2/3
Martinez-Bazan et al.28
8.2(εdi) Ω(di) = 0.25
Ω(di) = Luo and Svendsen
29
Lehr et al.30
∫0
c i
di
0.5
Ω(di , dj) dfbv
⎛ ε ⎞1/3 Ω(di , dj) = 0.9238αslu⎜ 2 ⎟ ⎝ di ⎠
Ω(di) =
∫0
Ω(di , dj) =
Ω(di) = Zhao and Ge31
σ
− 12 ρ d
∫0
∫ξ
(1 + ξ)2
1
11/3
ξ
min
⎛ 12σc ⎞ f exp⎜⎜ ξ−11/3⎟⎟ dξ 2/3 5/3 ⎝ ρslu ε di ⎠
0.5
Ω(di , dj) dfbv 1.19σ
∫ξ
ρslu ε1/3di 7/3fbv1/3
(1 + ξ)2
1
13/3
ξ
min
⎞ ⎛ 2σWecrit exp⎜⎜ ξ−2/3⎟⎟ dξ 2/3 5/3 1/3 ⎠ ⎝ ρslu ε di fbv
0.5
Ω(di , dj) dfbv
Ω(di) = 0.923(1 − φ)ε1/3
∫ξ
1
min
(1 + ξ)2 di
5/3
11/3
(ξ)
⎛ e (d , λ) ⎞ exp⎜− c i ⎟ dξ ei̅ (λ) ⎠ ⎝
⎞ ⎛ c πd 2σ πσλ ⎟ ec(di , λ) = max⎜⎜ f i , 3di min[(fbv ), 1 − (fbv )]1/3 ⎟⎠ ⎝ Ced
The volume fraction of daughter bubble f bv = (dj/di)3, the ratio of eddy size to the parent bubble diameter ξ = λ/di, the mean eddy energy ei̅ (λ) = πβLρslurryε2/3λ11/3/12, and the coefficient of surface increment during bubble breakup cf = f bv2/3 + (1 − f bv)2/3 −1. a
f (Eo) = 0.00105Eo3 − 0.0159Eo2 − 0.0204Eo + 0.474
where CB, CD, BB, and BD represent the birth rate and death rate due to coalescence and the birth rate and death rate due to breakup, respectively. By far, the population balance model has been mainly used in the simulation of air−water systems. In the present work, an attempt has been made to test the existing population balance models in the simulation of gas−liquid−solid flows with different particle loadings. As the first step of the effort, the effects of particle concentrations on the breakup rates calculated by different bubble breakup models (Table 1) have been investigated and are compared in Figure 2. The breakup rates of
(19)
where Eo is the Eötvös number. 2.2. Numerical Scheme and Benchmark Cases. The simulations were performed based on the experimental setup of Gandhi et al.18 The slurry bubble column was operated with particle loadings up to 40 vol % and gas velocities up to 0.26 m/s (ds = 35 μm, ρs = 2452 kg/m3). The diameter and height of the slurry bubble column were 0.15 m and 2.50 m, respectively. The static liquid height of the column was 1.50 m. Hamidipour et al.26 performed three-phase fluidized bed simulations with different solid wall boundary conditions (no slip, partial slip, and free slip), and discovered the simulation results in excellent agreement with experimental results when a near-zero value of the specularity coefficient was used. The same free-slip boundary condition was used by Schallenberg et al.14 and Panneerselvam et al.15 Therefore, in all simulations of the present paper, the free-slip wall boundary condition is used for the bubbles and particles. For bubbly flow, 3-D time dependent simulation is recommended.27 The mixture RNG k−ε turbulence model is used for the modeling of bubble induced turbulence. The QUICK discretization scheme is used for momentum and volume fraction equations, and the first order upwind scheme is employed to discretize other equations. More information about simulation numerical scheme and result processing can refer to our previous work.4
Figure 2. Comparison of different bubble breakup kernels with the change of solid concentration.
3. RESULTS AND DISCUSSION 3.1. Selection of Population Balance Model. The bubble number density equation of population balance model is expressed as
Martinez-Bazan, Luo and Svendsen, and Zhao and Ge increase with the increase of solid concentrations; only that of the Lehr et al. model decreases with the increase of solid concentrations. Therefore, the Lehr breakup kernel is supposed to be able to model the effect of solid concentration on the overall gas holdup. Subsequently, the Lehr breakup model and the Luo coalescence model are used to simulate the gas−slurry flows with 4925
dx.doi.org/10.1021/ie403453h | Ind. Eng. Chem. Res. 2014, 53, 4922−4930
Industrial & Engineering Chemistry Research
Article
different particle concentrations (0, 10, 20, 30 vol %) and a gas velocity of 0.25 m/s. The coalescence rate cij is defined as
In eq 21, tdrainage and tcontact represent film drainage time and bubble contact time, respectively. It is worth noting that C0 is an adjustable parameter and takes very different values in the existing collision frequency models. For example, in the models of Lee et al.,32 Prince and Blanch,33 and Luo,34 the values of C0 are 1.4142π/4, 0.28, and 1.43π/4, respectively. Theoretically, the bubble collision frequency is contributed by three different mechanisms: turbulence, buoyancy-driven, and laminar shear stress collisions. However, in the Luo coalescence model, only the turbulence collision is implemented and the other two mechanisms are ignored; thus a much larger value of C0 (∼1.1) was used. In our previous work,4 we found this value overestimated the bubble collision frequency, and 0.5 was used for C0 to better fit the experimental data. For small particles, the influence of particle concentration on bubble breakup and coalescence is the property change of the slurry instead of the particle−bubble collision. The effect of liquid viscosity on the bubble coalescence rates can be modeled by introducing a film drainage time; see, e.g., Chesters35 and Martiń et al.36 Since the Luo coalescence model does not take into account the influence of liquid viscosity on the bubble coalescence, a slightly smaller C0 of 0.4 is used to simulate a lower bubble coalescence rate due to an increase in slurry viscosity. Figure 3 shows the simulated gas
Figure 4. Comparison of Sauter mean diameters obtained by Wilkinson correlation and simulation method with different population balance models.
Figure 3. Comparison of gas holdups with particle loadings up to 30 vol % using different population balance models.
Figure 5. Gas holdup and liquid streamline of gas−slurry flow with 30 vol % particle concentration.
holdups using original and modified Luo−Lehr models and the experimental data as well. The Sauter mean bubble sizes are also simulated and compared with an empirical correlation as shown in Figure 4. The simulated bubble sizes as a function of particle loading by using the modified Luo−Lehr model agree with the Wilkinson model, which reads37 dg,Wilkinson = 3g −0.44 σ 0.34μ l 0.22 ρl −0.45 ρg −0.11U −0.02
shown in Figure 6. An increase of slurry viscosity leads to increases of bubble size and axial liquid velocities with the increases of the particle loading. Bubble size distribution of gas−liquid−solid flows will be simulated by the above modified Luo−Lehr model in the following simulations. 3.2. Comparison of Interphase Closure Models. As described in section 1, closure A is a gas−slurry two-phase model, and closures B, C, and D are three-phase models; i.e., bubbles, liquid, and solid particle three-phase momentum equations are solved, separately. While the particle settling cannot be simulated in the closure A model, the other three models do treat the liquid and solids as distinct phases, and the particle settling can be modeled. The newly constructed closure D model is similar to the closure A model in treating the liquid− solid as a slurry while computing the gas−liquid drag force; the
(22)
The gas holdup and liquid streamline of gas−slurry flow with 30 vol % particle concentration are shown in Figure 5. A high velocity of liquid is located in the center of the column, and complicated fine structures of reverse flow are formed in the near-wall zone. The radial distributions of bubble size and axial liquid velocity profiles with different particle loadings are simulated by the modified Luo−Lehr model, and the results are 4926
dx.doi.org/10.1021/ie403453h | Ind. Eng. Chem. Res. 2014, 53, 4922−4930
Industrial & Engineering Chemistry Research
Article
Figure 7. Time-averaged radial distributions of gas holdup and axial liquid velocities simulated with different interphase closure models.
Figure 6. Time-averaged radial distributions of Sauter diameters and axial liquid velocities with different particle loadings.
Table 2. Comparison of the Whole Column Gas Holdups Simulated by Four Closure Models with Experimental Value
closure D model treats the solids as a distinct phase and models the liquid−solid drag force in a similar way to the closure B and C models. The simulation results of the four closure schemes are compared and discussed as follows. The distributions of gas volume fractions and axial liquid velocities of a slurry bubble column operated at a superficial gas velocity of 0.25 m/s and a slurry concentration of 0.2 were simulated using the closure A, B, C, and D models, respectively. The simulated radial distributions of axial liquid velocities and gas holdups are shown in Figure 7. As can be found in Figure 7, since closure B underestimates the interaction between bubbles and other phases, the simulated liquid velocity of closure B is relatively low, and the gas holdup is underestimated. Taking the gas−solid interaction into account, the bubble lifting resistance increases, and the gas holdup is overestimated in closure C. The axial liquid velocity increases with the increase of gas holdup. The overall gas holdups simulated by the four closure models are compared in Table 2. The gas holdups of closure A and closure D are in good agreement with the experimental data. The closure B model underestimates and the closure C model overestimates the gas holdup. Figure 8 compares the axial velocities of liquid and solid simulated by closure C. It is interesting to find that the solid velocity is higher than the liquid’s. Therefore, closure C is incapable of predicting the particle settling phenomenon when the particle size in the slurry reactor is much smaller than the bubble size. The simulated gas volume fraction and axial liquid velocity of closure D are similar to those of closure A due to the contributions of solids to the increase of gas−slurry interactions. 3.3. Comparison of Drag Force Models. The axial distribution of the solid phase is largely dependent on the
gas holdup
expt
closure A
closure B
closure C
closure D
0.185
0.173
0.154
0.279
0.169
Figure 8. Comparison of axial velocities of liquid and solid simulated by closure C.
superficial gas velocity and the properties of the liquid and particles. Smith and Ruether38 studied the particle dynamics in a slurry bubble column and proposed a one-dimensional sedimentation dispersion model. The model is represented as follows:
⎛ z⎞ Cs = Cs0 exp⎜B ⎟ ⎝ H⎠
(23)
Cs0 = CsB(exp(B) − 1)
(24)
B=− 4927
ψlu pH Ds
(25) dx.doi.org/10.1021/ie403453h | Ind. Eng. Chem. Res. 2014, 53, 4922−4930
Industrial & Engineering Chemistry Research
Article
u p = 1.10ug 0.026ut 0.80ψl 3.5 ugDT Ds
= 9.6(Frg 6/Reg)0.1114 + 0.019Rep1.1
(26)
(27)
In eqs 26 and 27, Frg, ut, Reg, and Rep represent the Froude number, particle terminal settling velocity, gas Reynolds number, and particle Reynolds number, respectively. For a cylindrical slurry bubble column, an accurate solid concentration distribution can be obtained using the Smith− Ruether empirical correlation (Figure 9a). When the geometry shape is complex or internals are presented in a column, the empirical model is no longer applicable. In this situation, CFD provides an alternative approach. The Schiller−Naumann model, Wen−Yu model, and EMMS model are used to simulate the axial solid concentrations in slurry bubble columns with different particle concentrations. The axial solid concentration distributions simulated by the Schiller−Naumann, Wen−Yu, and EMMS liquid−solid drag models are shown in Figure 9, parts b, c, and d, respectively. Figure 9 shows very little difference between the axial solid concentration distributions of the slurry bubble columns with different particle loadings simulated with the Schiller− Naumann model, which does not consider the effect of the volume fraction of the continuous phase as indicated in eqs 11 and 12. The simulated results obtained with the Wen−Yu and EMMS models reveal that the particle settling tendency becomes weaker with the decrease of solid holdup. The simulation results are in good agreement with the experimental results. The distributions of solid holdup computed with the Wen−Yu model are more homogeneous than that of the EMMS model. The reason for a large variation obtained with the EMMS model is that the model takes the local heterogeneous structure into account for the drag force computation. Similar results can be found in the simulations of gas−solid flow of Yang et al.17 and Chalermsinsuwan et al.39
4. CONCLUSION A slurry bubble column was simulated using the Euler−Euler approach implemented with a population balance model. A pseudo gas−slurry closure model was proposed and applied in the predictions of gas holdup and axial liquid velocity, and the results were compared with the existing closure models. Three different drag force models, i.e., the Schiller−Naumann model, the Wen−Yu model, and the EMMS model, were used for the computation of the liquid−solid interaction force, and the simulated particle settlings were compared with experiments. The following conclusions were drawn from the simulation results: (1) The Luo−Lehr population balance model could simulate the mean bubble size changes as a function of particle loading. The gas holdup and its variation trend predicted by the modified Luo−Lehr population balance model were in agreement with the experimental data by setting the Luo coalescence coefficient C0 = 0.4. (2) A pseudo gas−slurry closure model which took the influence of solid concentration into account was recommended for gas−liquid−solid three-phase simulation when the particle size was much smaller than the bubble size. The simulated average gas holdups showed good agreement with the experimental data. (3) The axial distribution of solid concentration simulated by using the EMMS drag model was much closer to the experiments than those of other drag models for the simulation of slurry bubble columns.
Figure 9. Comparison of radial distributions of solid concentration calculated with (a) one-dimensional sedimentation dispersion model and CFD method using (b) Schiller−Naumann, (c) Wen−Yu, and (d) EMMS drag models.
where up and Ds represent the particle hindered settling velocity and axial dispersion coefficient which are calculated by the following equations, respectively: 4928
dx.doi.org/10.1021/ie403453h | Ind. Eng. Chem. Res. 2014, 53, 4922−4930
Industrial & Engineering Chemistry Research
■
Article
(2) Chen, C.; Fan, L. S. Discrete simulation of gas-liquid bubble columns and gas-liquid-solid fluidized beds. AIChE J. 2004, 50, 288. (3) Buwa, V. V.; Deo, D. S.; Ranade, V. V. Eulerian-Lagrangian simulations of unsteady gas-liquid flows in bubble columns. Int. J. Multiphase Flow 2006, 32, 864. (4) Xu, L.; Yuan, B.; Ni, H.; Chen, C. Numerical simulation of bubble column flows in churn-turbulent regime: comparison of bubble size models. Ind. Eng. Chem. Res. 2013, 52, 6794. (5) Chen, P.; Sanyal, J.; Dudukovic, M. P. CFD modeling of bubble columns flows: implementation of population balance. Chem. Eng. Sci. 2004, 59, 5201. (6) Chen, Y.; Fan, L. Bubble breakage mechanisms due to collision with a particle in liquid medium. Chem. Eng. Sci. 1989, 44, 117. (7) Hooshyar, N.; van Ommen, J. R.; Hamersma, P. J.; et al. Dynamics of Single Rising Bubbles in Neutrally Buoyant Liquid-Solid Suspensions. Phys. Rev. Lett. 2013, 110, 244501. (8) Luo, X.; Lee, D. J.; Lau, R.; Fan, L. Maximum stable bubble size and gas holdup in high-pressure slurry bubble columns. AIChE J. 1999, 45, 665. (9) Maretto, C.; Krishna, R. Modelling of a bubble column slurry reactor for Fischer−Tropsch synthesis. Catal. Today 1999, 52, 279. (10) Troshko, A. A.; Zdravistch, F. CFD modeling of slurry bubble column reactors for Fischer−Tropsch synthesis. Chem. Eng. Sci. 2009, 64, 892. (11) Mitra-Majumdar, D.; Farouk, B.; Shah, Y. T. Hydrodynamic modeling of three-phase flows through a vertical column. Chem. Eng. Sci. 1997, 52, 4485. (12) Rampure, M. R.; Buwa, V. V.; Ranade, V. V. Modelling of GasLiquid/Gas-Liquid-Solid Flows in Bubble Columns: Experiments and CFD Simulations. Can. J. Chem. Eng. 2003, 81, 692. (13) Padial, N. T.; VanderHeyden, W. B.; Rauenzahn, R. M.; et al. Three-dimensional simulation of a three-phase draft-tube bubble column. Chem. Eng. Sci. 2000, 55, 3261. (14) Schallenberg, J.; Enß, J. H.; Hempel, D. C. The important role of local dispersed phase hold-ups for the calculation of three-phase bubble columns. Chem. Eng. Sci. 2005, 60, 6027. (15) Panneerselvam, R.; Savithri, S.; Surender, G. D. CFD simulation of hydrodynamics of gas-liquid-solid fluidised bed reactor. Chem. Eng. Sci. 2009, 64, 1119. (16) Visuri, O.; Wierink, G. A.; Alopaeus, V. Investigation of drag models in CFD modeling and comparison to experiments of liquid− solid fluidized systems. Int. J. Miner. Process. 2012, 104, 58. (17) Yang, N.; Wang, W.; Ge, W.; et al. CFD simulation of concurrent-up gas−solid flow in circulating fluidized beds with structure-dependent drag coefficient. Chem. Eng. J. 2003, 96, 71. (18) Gandhi, B.; Prakash, A.; Bergougnou, M. A. Hydrodynamic behavior of slurry bubble column at high solids concentrations. Powder Technol. 1999, 103, 80. (19) Brian, B. W.; Chen, J. C. Surface tension of solid-liquid slurries. AIChE J. 1987, 33, 316. (20) Thomas, D. G. Transport characteristics of suspension: VIII. A note on the viscosity of Newtonian suspensions of uniform spherical particles. J. Colloid Science 1965, 20, 267. (21) Ishii, M.; Zuber, N. Drag coefficient and relative velocity in bubbly, droplet or particulate flows. AIChE J. 1979, 25, 843. (22) Schiller, L.; Naumann, Z. Ü ber die grundlegenden Berechnungen bei der Schwerkraft-aufbereitung. Z. Ver. Dtsch. Ing. 1935, 77, 318. (23) Wen, C. Y.; Yu, Y. H. Mechanics of Fluidization. Chem. Eng. Prog. Symp. Ser. 1966, 62, 100. (24) Hong, K.; Wang, W.; Zhou, Q.; et al. An EMMS-based multifluid model (EFM) for heterogeneous gas−solid riser flows: Part I. Formulation of structure-dependent conservation equations. Chem. Eng. Sci. 2012, 75, 376. (25) Tomiyama, A.; Sou, I.; Zun, I.; Kanami, N.; Sakaguchi, T. Effects of Eötvös number and dimensionless liquid volumetric flux on lateral motion of a bubble in a laminar duct flow. Adv. Multiphase Flow 1995, 3. (26) Hamidipour, M.; Chen, J.; Larachi, F. CFD study on hydrodynamics in three-phase fluidized bedsApplication of
AUTHOR INFORMATION
Corresponding Author
*E-mail:
[email protected]. Notes
The authors declare no competing financial interest.
■
ACKNOWLEDGMENTS This work is supported by China’s National Science Foundation (NSFC) under Grants 21276085 and 20876049.
■
NOTATION B = dimensionless parameter defined in eq 23 BB = production of bubble by breakup, m−3·kg−1·s−1 BD = destruction of bubble by breakup, m−3·kg−1·s−1 c = coalescence rate, m−3·s−1 CB = production of bubble by coalescence, m−3·kg−1·s−1 CD = destruction of bubble by coalescence, m−3·kg−1·s−1 Cd = drag coefficient, dimensionless Clift = lift coefficient, dimensionless Cs = particle concentration, kg·m−3 d = bubble or solid particle diameter, m Ds = axial dispersion coefficient, dimensionless Eo = Eötvös number, dimensionless f = volume fraction, Vi/Vj F = interface momentum exchange term Flift = lift force, N·m−3 Fr = Froude number, dimensionless g = gravity acceleration, m·s−2 g(V′) = breakup frequency, s−1 H = liquid height, m n = number of bubbles per unit volume, m−3 p = pressure, Pa Re = Reynolds number, dimensionless t = time, s tcontact = bubble contact time, s tdrainage = film drainage time, s T = temperature, K u = velocity, m·s−1 up = particle hindered settling velocity, m·s−1 ut = particle terminal settling velocity, m·s−1 V = bubble volume, m3 z = axial position from bottom of column, m
Greek Symbols
α = void fraction, dimensionless β(V,V′) = probability density function (pdf), m−3·kg−1 ξ = eddy size divided by parent bubble size ρ = density, kg·m−3 μ = viscosity, kg·m−1·s−1 σ = surface tension, N·m−1 ψ = volume fraction in slurry ε = turbulent dissipation rate, m2·s−3 ΩB(vi,vj) = breakup frequency, m−3·s−1 Subscripts
g = gas phase l = liquid phase s = solid phase slu = slurry
■
REFERENCES
(1) Kantarci, N.; Borak, F.; Ulgen, K. O. Bubble column reactors. Process Biochem. 2005, 40, 2263. 4929
dx.doi.org/10.1021/ie403453h | Ind. Eng. Chem. Res. 2014, 53, 4922−4930
Industrial & Engineering Chemistry Research
Article
turbulence models and experimental validation. Chem. Eng. Sci. 2012, 78, 167. (27) Ekambara, K.; Dhotre, M. T.; Joshi, J. B. CFD simulations of bubble column reactors: 1D, 2D and 3D approach[J]. Chem. Eng. Sci. 2005, 60, 6733. (28) Martinez-Bazan, C.; et al. On the breakup of an air bubble injected into a fully developed turbulent flow. Part 1. Breakup frequency. J. Fluid Mech. 1999, 401, 157. (29) Luo, H.; Svendsen, H. F. Theoretical model for drop and bubble breakup in turbulent dispersions. AIChE J. 1996, 42, 1225. (30) Lehr, F.; Millies, M.; Mewes, D. Bubble-size distributions and flow fields in bubble columns. AIChE J. 2002, 48, 2426. (31) Zhao, H.; Ge, W. A theoretical bubble breakup model for slurry beds or three-phase fluidized beds under high pressure. Chem. Eng. Sci. 2007, 62, 109. (32) Lee, C. H.; Erickson, L.; Glasgow, L. Bubble breakup and coalescence in turbulent gas-liquid dispersions. Chem. Eng. Commun. 1987, 59, 65. (33) Prince, M. J.; Blanch, H. W. Bubble coalescence and break-up in air-sparged bubble columns. AIChE J. 1990, 36, 1485. (34) Luo, H. Coalescence, breakup and liquid circulation in bubble column reactors. D.Sc. Thesis, Norwegian Institute of Technology, Trondheim, 1993. (35) Chesters, A. K. The modelling of coalescence processes in fluid liquid dispersion: a review of current understanding. Trans. Inst. Chem. Eng. 1991, 69, 259. (36) Martín, M.; Montes, F. J.; Galán, M. A. Mass Transfer Rates from Oscillating Bubbles in Bubble Columns Operating with Viscous Fluids. Ind. Eng. Chem. Res. 2008, 47, 9527. (37) Wilkinson, P. M. Physical aspects and scale-up of high pressure bubble columns. Ph.D. Thesis, University of Groningen, Netherlands, 1991. (38) Smith, D. N.; Ruether, J. A. Dispersed solid dynamics in a slurry bubble column. Chem. Eng. Sci. 1985, 40, 741. (39) Chalermsinsuwan, B.; Piumsomboon, P.; Gidaspow, D. Kinetic theory based computation of PSRI riser: Part IEstimate of mass transfer coefficient. Chem. Eng. Sci. 2009, 64, 1195.
4930
dx.doi.org/10.1021/ie403453h | Ind. Eng. Chem. Res. 2014, 53, 4922−4930
