Krist V. Gernaey, Jakob K. Huusom and Rafiqul Gani (Eds.), 12th International Symposium on Process Systems Engineering and 25th European Symposium on Computer Aided Process Engineering. 31 May – 4 June 2015, Copenhagen, Denmark © 2015 Elsevier B.V. All rights reserved.
Modelling the Hydrodynamics of Bubble Columns using Coupled OPOSPM-Maximum Entropy Method Menwer Attarakiha,b*, Ferdaous Al-Slaihata, Mark W. Hlawitschkab, Hans-Jörg Bartb a
Department of Chemical Engineering, University of Jordan, Amman 11942, Jordan Chair of Separation Science and Technology, TU Kaiserslautern, 67653, Kaiserslautern, Germany
[email protected] b
Abstract In this contribution, we used a reduced population balance model to describe the hydrodynamics of bubble columns, which play a major role in in determining the bubble size distribution and hence the interfacial area concentration. This model consists of a set of transport equations to track the total number and total interfacial area concentrations and the gas phase volume fraction of bubbles. The model is essentially derived using the One Primary and One Secondary Particle Method (OPOSPM) and its higher extension using an Implicit Two-Equal Weight Quadrature (TwoEqWQ). This is coupled with the Shannon Maximum Entropy Method to predict the bubble size distribution along the bubble column axial direction. The model predictions show good agreement with the published experimental data in the bubbly flow regime. Keywords: Bubble column, Population balances, Maximum Entropy Method.
1. Introduction Bubble columns are multiphase equipment which is widely used in the chemical, petrochemical and biochemical industries (Kantarci et al., 2005). Estimation of the specific interfacial area is necessary in bubbly flows (e.g. the two-fluid model) to provide a closure for the momentum, mass and energy transport equations which are solved for each phase. The design and operation of these two-phase equipment are dependent on the adequate modelling of these transport phenomena (Jakobsen, 2008). In such equipment, the gas phase is dispersed in a continuous environment with a population of bubbles that may undergo breakage, coalescence, expansion and growth. Modelling such a process using the mixture models, which could not take into account these instantaneous discrete events, is not sufficiently accurate to predict the hydrodynamics and mass transport. On the other hand, using the multi-fluid models which could include bubble-bubble interactions were proved to be computationally expensive when detailed column geometry is concerned. For example, typical number of groups is a 15-20 as used by the MUSIG model (Jakobsen, 2008). As a result, the need for a transport equation to describe the transport of the interfacial area concentration has been reflected by the interest of many researchers during the last two decades (Ishii and Hibiki, 2010). Among these is the work of Hibiki and Ishii (2000) which is of great importance. They used the population balance equation (PBE) to derive an Interfacial Area Transport Equation (IATE) for bubbly flow in vertical tubes using section averaging. As the PBE is concerned, many issues arise with respect to the
204
M. Attarakih et al.
accuracy and consistency of the reduced interfacial transport equation (Drum et al., 2010). These authors developed a two-equation model that is consistent w.r.t the continuous PBE and hence its accuracy can easily be extended. This is unlike the IATE of Ishii and Hibiki (2000) that could not be further improved in terms of accuracy and consistency w.r.t. the continuous PBE. In this work, we used a reduced population balance model to describe bubble transport, expansion, breakage and coalescence . The reduced population balance model consists of a set of transport equations to track the total number and volume concentrations, which are coupled through mean mass diameter (d30), while the interfacial area concentration is calculated from the number and volume concentrations of bubbles. This model is essentially derived using the One Primary and One Secondary Particle Method (OPOSPM), which is coupled with the Shannon Maximum Entropy Method (MaxEntM) to predict the bubble size distribution along the bubble column axial direction. To increase the quadrature order for approximating source terms in the continuous PBE, an implicit TwoEqWQ is proposed, which require an additional transport equation for interfacial area concentration.
2. Bubble hydrodynamics The present model for the bubble column hydrodynamics is derived from the population balance equation which is a special type of the Boltzmann transport equation. This equation evolves the number concentration function in space, time and the particle property space. The transported quantity, which is used to describe the population of bubbles, is the number concentration function f(x,d,t). This is assumed to be continuous, satisfies regulatory conditions and describe the number density concentration of bubbles moving with velocity (¢vg²) at a given time t, in a spatial space (x) and along particle size (d). The rth moment of this bubble number concentration function is written as: f w\ wmr + ( ¢vg ² mr ) = ³ w] ]& f ( x, ] , t )d ] + Srb + Src (1) 0 wt
In Eq.(1), ¢vg² is the mean bubble velocity which is function of a pre-defined mean bubble diameter, \ = dr (where r = 0, 1,2, …), ] = v (where v is the bubble volume) and the dot over ] denotes derivative with respect to time. The last two source terms on the right hand side Sb and Sc describes the the rth moment of particle appearance by breakage and disappearance by coalescence respectively(Attarakih and Bart, 2014).These sources terms contain unclosed integrals due to the general dependency of bubble breakage and coalescence kernels on bubble size (Hibiki and Ishii, 2000). To close these source terms, Hibiki and Ishii (2000) and Ishii and Hibiki (2010) used two sets of bubble mean diameters: The mean volume (d30) and the Sauter mean diameters (d32) and derived what is called the Interfacial Area Transport Equation (IATE) followed by the simplifications of the source terms by assuming an equal-size binary bubble breakage and the bubble coalescence was allowed by only two bubbles of the same size. By following Coulaglou and Tavlarides (1977), they derived bubble breakage and coalescence kernels which were included in the source terms (Sb and Sc). Unfortunately, these kernels were merged with the numerical parameters which were used to close the integrals in the IATE. This approach has two major disadvantages: Firstly, it complicates the addition of new bubble kernels to the IACE, and secondly, it prevents increasing the numerical accuracy of the model by further increasing the order of quadrature approximation of the integrals in the source term (Drumm et al., 2010). Due to the mixing between two types of mean droplet diameters: Sauter mean diameter (d32) and volume mean diameter (d30), conflicting source term for the IATE with that
205
Modelling the Hydrodynamics of Bubble Columns uisng Coupled OPOSPM-Maximum Entropy Method
derived from the PBE is observed as shown in Table (1). This is because d30 is a quadrature node based on OPOSPM (Drumm et al., 2010) while d32 is only a physical definition which is not related to a quadrature approximation of the source terms in the PBE. Therefore, using this mathematically consistent mean droplet diameter (d30), the bubbly flow hydrodynamics is derived from Eq.(1) using the One Primary and One Secondary Particle Method (OPOSPM) as a reduced population balance model. OPOSPM for the bubble flow is obtained by using the zero, 2nd and 3rd moments (r =0, 2, 3) from Eq.(1) with proper transformations to obtain the total bubble number concentration N = m0, interfacial area concentration a = Sm2 and volume fraction D = (S/6)m3: wN + ( ¢ vg ² N ) = wt
wa + ( ¢ vg ² a ) = wt w (Į) wt
2 § a · ª wD º §D · + ( ¢vg ²D ) » + ¨ ¸ ¨ ¸« 3 © D ¹ ¬ wt ¼ ©a¹
2
¦kS
j =b,c
j
¦S
j =b,c
j 0
1 ª wĮ º + . ( ¢ v g ² Į ) = « + ( ¢ v g ² Į ) » 3 ¬ wt ¼
j 0
(2)
(3)
(4)
1 b Where S0 = *(d30 ) N and S0c = Z (d30 , d30 ) N 2 . The mean gas velocity for counter 2 current flow is defined as ¢vg ² = j f / (1 D ) + ur , where jf is the superficial velocity of the liquid phase and ur is the relative velocity of the bubbles which is given by (Ishii and Hibiki, 2010): ur =
d 30 g 'U 3C D U f
(5)
In the above equation, CD is the drag coefficient that is function of the mixture Reynolds number. Note that in Eqs.(3 and 4) the term in the square brackets represent bubble expansion due to pressure drop in the physical space. It is worthwhile to mention that the interfacial area concentration could be predicted from a = 6D/d32 using the N and D transport equations. However, we proposed here Eq.(3) for consistent prediction of (a) since d32 could not be estimated from the N and Dtransport equations (only d30 can be expressed in terms of N and D). Table (1): Constants for the source term in the interfacial area concentration equation. OPOSPM
kb kc
36S ( 2
1/ 3
36S ( 2
2/3
Ishii and Hibiki (2010)
1)
2)
36S / 3 36S / 3
According to the model numerical constants (kb and kc) in Table (1), the bubble breakage and coalescence kernels (* and Z respectively) of Hibiki and Ishii (2000) with their adjustable parameters need to be modified. This is to recover the internal consistency with respect to the continuous population balance equation, which allows the increase of the quadrature order approximation to achieve the desired accuracy in closing the integral source term.
206
M. Attarakih et al.
3. The Maximum Entropy Method The Quadrature Method of Moments (QMOM) is one of the most used solver for the PBE because of the low number of transport equations and its high accuracy. This is required to overcome the long CPU time when coupling discrete population balances with CFD codes (Drum et al.,, 2010). The QMOM solves the PBE, where a finite set of population moments is transported. Due to its inherit averaging principle, the QMOM cannot reproduce the particle size distribution, which is replaced by a few set of non-physical moving particles along the particle property space. However, in industrial particulate system applications, the particle size distribution has a significant impact on the physical, chemical and mechanical product properties (Attarakih and Bart, 2012). To preserve the advantages of the QMOM and concentrating on distribution reconstruction away from the discrete methods, the MaxEntM seems to be natural candidate for the problem at hand. The MaxEntM is used to estimate the least biased probability density function subject to the available moment information about the sought distribution. The consistency of the recovered solution is achieved by maximizing the Shannon entropy functional with the constraint that the available moments of the averaged distribution are reproduced from the reconstructed one. Following Attarakih and Bart (2014), this constrained nonlinear program has a unique positive solution provided that the sequence of the targeted moments is completely monotone. This maximum entropy solution is function of Lagrange multipliers, which can be found by minimizing the a convex potential function or using the MaxEntM: § N · f ( x, d , t ) = d G exp ¨ ¦ O j d j ¸ j = 0 © ¹
(6)
In Eq.(6) Oj’s (Oj corresponds to mj) are the Lagrange multipliers which guarantee the reproduction of the targeted moments and G is a positive real number that is introduced here to satisfy the regulatory conditions (f(0) = f(dof) = 0). In this work, we selectively chose {m0, m2, m3} as a set of moments which are satisfied by the transport equations (Eq. (1)). In addition to this, an Implicit Two-Equal Weight Quadrature (TwoEqWQ) was derived which reproduces exactly the above targeted set of moments. This TwoEqWQ is used to close the integral source term in Eq.(1) which is given by the following set of weights and nodes w1 = w2 = 0.5P0 , d1 = z , d 2 = 2mˆ 2 z andz is a solution of {(2 mˆ 2 z ) (2 mˆ 3 z
{
3/ 2 2 / 3
)
}
= 0} . The hat denotes normalization w.r.t. m0.
ˆ 3 )2/3 ]. This equation has two roots, where the first onelies in the interval [0, (m
4. Numerical results and discussion 4.1. Analytical validation In this section, the proposed two-equal weight implicit quadrature to close the source term in Eq.(1) is validated analytically along with the MaxEntM to reconstruct the averaged distribution. The Lagrange multipliers (Oj, j = 0, 2,3) were found by solving a convex potential function (Attarakih and Bart, 2014) using the Levenberg-Marquardt algorithm. The case study is particle aggregation with constant kernel and an exponential (w.r.t. particle volume) initial condition. The left panel of Figure (1) shows the evolution of the targeted moments {m0, m2, m3} as function of dimensionless time (W). Due to the remarkable accuracy of of the QMOM, not only the targeted moments were produced with high accuracy, but also the interpolated moment (m1).
Modelling the Hydrodynamics of Bubble Columns uisng Coupled OPOSPM-Maximum Entropy Method 1.2 0.8
analytical MaxEntM
m2
0.6 0.4 0.2 0
0.04 0.02
m0 0
analytical MaxEntM
0.06
f(x)
moments
0.08
m3
1
207
4
8
12
16
0
20
0
2
W (-)
4
6
x
Figure (1): Validation of the MaxEntM with two implicit equal-weight quadrature for constant kernel aggregation in batch system against analytical solution (Gelbard and Seinfeld, 1978).
This moment was predicted using the implicit two-equal weight quadrature and the reconstructed MaxEnt solution (Eq.(6)) with Lagrange multipliers (O0 = -3.830, O1 = 0.090, O3 = -0.115) at W = 20. The mean percentage relative error in the prediction of m1 using the two-weight quadrature is 0.311 and that predicted from the MaxEntM is 0.519. The recovered distribution at W = 20 is shown on the right panel of Figure (1) and compared to the analytical solution (Gelbard and Seinfeld, 1978). It is clear that the coupled implicit TwoEqWQ and the MaxEntM achieved two targets: Accurate tracking of distribution moments and accurate reproduction of the particle number concentration.
volume fraction
0.8
0.4 OPOSPM TwoEqWQ IATE Symbols: Experiment
0.6
d32 = 2.65 mm D = 0.049
Case 3
0.4
0.2
Case 2
0.2 0
-1 0.3 a = 111 m
0.1
Case 1 0
20
Z/D
40
z/D = 53.5 Case I
60
0 2.2
2.4
2.6
2.8
bubble diameter (mm)
3
Figure (2): Experimental and numerical validation of OPOSPM, the implicit TwoEqWQ using air-water bubbly flows in a 50.8 mm diameter column (Hibiki and Ishii, 2000) . In Case 1: jg = 0.0275 m/s , jf =0.491 m/s, d32in = 2.49 mm, Case 2: jg = 0.19 m/s, jf = 0.491 m/s, d32in= 3.31 mm, and Case 3: jg = 3.90 m/s, jf = 5.00m/s, d32in = 2.83.
4.2. Experimental validation As shown in the previous section, the results show that the present model can be extended to higher-order approximations of the PBE using more secondary particles(quadrature nodes) with distribution reconstruction. This is an obvious advantage over the existing models such as that of Hibiki and Ishii (2000). As a first step of model validation, the pilot plant bubble column of Hibiki and Ishii (2000) was modeled and simulated after the necessary modifications of the adjustable parameters in the breakage and coalescence kernels to recover the mathematical consistency with the PBE. The column diameter is 50.8 mm with 3060 mm height which is operated adiabatically in the bubbly flow regime using air-water mixture in a cocurrent flow. The
208
M. Attarakih et al.
investigated inlet superficial gas velocity ranges from 0.0275 to 3.90 m/s and that for the liquid phase is from 0.491 to 5.00 m/s with a void fraction in the range 0.01270.468. The axial measuring points (z/D) are at 6, 30.3 and 53.5 with detailed boundary conditions as described by Hibiki and Ishii (2000). Figure (2-left) shows the variation of the gas void fraction along column height and the recovered bubble size distribution (right) as compared to the experimental data Hibiki and Ishii (2000) using three different models (OPOSPM, Implicit TwoEqWQ and IATE). The results are almost identical after adjusting the IATE model parameters to reflect the changes imposed by the new consistent formulation (see Table 1). As discussed by Hibiki and Ishii (2000),the simulated bubbly flow is dominated by bubble expansion for Case 1 due to the axial pressure drop along the column height. This is reflected by the increase in the mean size of the bubbles when compared to the mean inlet size (2.49 mm) as shown in Figure(2-right).For this case, the Lagrange multipliers used in Eq.(6) are found to be ((O0 = -369.825, O1 = 156.147, O3 = -39.335)).
5. Summary and conclusions In this work, we derived the transport equations of bubble column hydrodynamics using OPOSPM and an implicit TwoEqWQ as reduced population balance models. The equations from the latter method were found consistent with the bubble number, interfacial area and volume concentrations and the Sauter mean bubble diameter. The coupling of these reduced bubble hydrodynamic models with the constrained MaxEntM allowed the reconstruction of a unique bubble number concentration by minimizing a convex potential function. A very good agreement between the predictions of the present model and those of the IATE and the Experimental data was obtained.
References M. Attarakih and H.-J. Bart, 2014, Solution of the population balance equation using the differential maximum entropy method (DMaxEntM): An application to liquid extraction columns. Chem. Eng. Sci, 108, 123-133. M. Attarakih and H.-J. Bart, 2012, Integral formulation of the smoluchowski coagulation equation using the cumulative quadrature method of moments (CQMOM), Computer - Aided Chem. Eng., 31, 1130-1134. C. A. Coulaloglou and L. L. Tavlarides,1977, Description of interaction processes in agitated liquid-liquid dispersions, Chem. Eng. Sci., 32, 1289-1297. C. Drumm, M. Attarakih, M. W. Hlawitschka and H.-J. Bart, 2010, One-group reduced population balance model for cfd simulation of a pilot-plant extraction column. Ind. Eng. Chem. Res., 49, 3442–3451. F. Gelbard and J. H. Seinfeld, 1978, Numerical solution of the dynamic equation for particulate systems. J. Comput. Phys., 28, 357-375. T. Hibiki and M. Ishii, 2000, One-group interfacial area transport of bubbly flows in vertical round tubes. Int. J. Heat Mass Tran., 43, 2711 - 2726. M. Ishii and T. Hibiki, 2010, Thermo-fluid dynamics of two-phase flow, Springer New York. H. A. Jakobsen, 2008, Chemical reactor modeling: Multiphase reactive flows, Springer-Verlag, Heidelberg. N. Kantarci, F. Borak and K. O. Ulgen, 2005, Bubble column reactors: Review. Process Biochemistry, 40, 2263–2283.