- Email: [email protected]
A continuum model for gas–liquid flow in packed towers
Chemical Engineering Science 56 (2001) 5945–5953
www.elsevier.com/locate/ces
A continuum model for gas–liquid %ow in packed towers Shijie Liu ∗ Department of Chemical and Materials Engineering, University of Alberta, Edmonton AB, Canada T6G 2G6
Abstract Gas–liquid %ow in packed towers is commonly encountered in the chemical and processing industry. A continuum model is developed based on the volume-and-time averaging of multiphase %ows in isotropic rigid porous media=packed columns. Closures are presented for the evaluations of the extra surface=intrinsic phase integral terms. Both inertia and inter-phase interactions are retained in the volume averaged (Navier–Stokes) equations. These governing equations are solved for fully-developed axi-symmetric single and gas–liquid two phase %ows in highly porous packed towers. It is found that the dispersion term is present in the continuity equation as well as the momentum equations. Numerical simulations with the models show that the volume-and-time averaged equations can predict the velocity, phase hold-up and pressure drop quite well for up to the loading point for gas–liquid counter-current %ows. ? 2001 Elsevier Science Ltd. All rights reserved. Keywords: Packed tower; Multiphase %ow; Volume-and-time averaging; Counter-current %ow; Velocity distribution; Phase hold-up
1. Introduction There are three gas–liquid %ow arrangements useful for separation and reaction engineering applications. There is only one possible way for countercurrent %ow arrangement with the heavy %uid (liquid) %owing downward and the light %uid (gas) %owing upward. This countercurrent %ow arrangement is commonly used for distillation, cooling tower, absorption and stripping operations (e.g., Kister, 1992; Billet, 1995). There are two possible arrangements for cocurrent %ows: either upward or downward %ow. While concurrent %ows can be applied for mass and=or heat transfer operations, they are generally used as multiphase reactors (e.g., Shah, 1979; Herskowitz & Smith, 1983). In particular, the downward %ow arrangement is commonly used. In this case, the liquid rate can be small and the unique %ow system is called a trickle bed (e.g., Saroha, Nigam, Saxena, & Kapoor, 1998; Iliuta, Ortiz-Arroyo, Larachi, Grandjean, & Wild, 1999). The cocurrent up %ow requires higher pressure drop and liquid hold-up in tower than the down %ow mode. A continuum model for numerical simulation of %ows in such systems will be developed in this paper. Closure models will be presented suitable for gas–liquid countercurrent %ows. ∗
Tel.: +1-780-492-0981; fax: +1-780-492-2881. E-mail address: [email protected] (S. Liu).
For a laminar saturated single %uid %ow in a straight duct, the pressure drop is proportional to the %ow rate. The same relationship holds also for %ow in curved ducts when the %ow is weak (Liu, Afacan, Nasr-El-Din, & Masliyah, 1994). This unique relationship can be generalized to %ow through porous media as well, dp − gx ˙ −U; dx
(1)
where U is the superEcial %uid velocity, p is pressure, x is the axis of the %ow direction, is the %uid density and gx is the gravity in the direction of %ow. This limiting %ow behavior has received great attention for centuries. The direct expression resulted from the proportionality for %ow in porous media is called the Darcy’s law. Several approaches have been used for the semitheoretical modeling of packed towers. When the %ow velocity is high, Eq. (1) is not strictly adequate. In this case, all the existing models converge to a simple (and experimentally observed) relationship, which is commonly expressed as −
1 dp + G gx = CD G UG2 ; dx 2
(2)
where CD is the so-called packing factor and has a reciprocal length unit. It is inversely proportional to the particle size and is a strong function of the packing
0009-2509/01/$ - see front matter ? 2001 Elsevier Science Ltd. All rights reserved. PII: S 0 0 0 9 - 2 5 0 9 ( 0 1 ) 0 0 2 3 2 - 9
5946
S. Liu / Chemical Engineering Science 56 (2001) 5945–5953
porosity as well. The packing factor is normally supplied by manufactures for each speciEc random packings. Because the expression is usually used for dry gas %ow through packed beds, a subscript G was used in Eq. (2). A collection of data for the packing factor is also available in textbooks, for example, Trebal (1980). Richards (1931) and Muskat (1937) intuitively extended the Darcy’s law to multiphase %ows. In their extension, the equation remains the same for each phase but allowing the %uid properties as well as the permeability to diJer. That is, ˜vi = −
ki ( pi − i˜g); i
(3)
where the subscript i denotes the ith %uid phase, k is the permeability and is the dynamic viscosity. Although widely used for %ow in soils and sandstones, Eq. (3) lacks three key features (Liu, 1999b): namely, inertial eJects; diJusion eJects and the inter-phase interactions (Dullien, 1992). These defects make Eq. (3) inapplicable to the gas–liquid %ows in packed towers. In this paper, we shall apply the volume-and-time averaging principles to derive a governing equation in packed towers.
2. Volume-and-time averaged Navier–Stokes equations for multiphase ows in porous media
(4)
and the ith phase momentum equation is given by @(∗i ˜vi∗ ) + · (∗i ˜vi∗˜vi∗ ) + p∗ − ∗i ˜g @t − · i∗ [˜vi∗ + (˜vi∗ )T ] = 0:
@(i hi ) + · (i hi˜u i ) − · (i˜u di ) = 0: @t
(6)
Eq. (6) is similar to the traditional phase balance (or continuity) equation except for the dispersion %ux ˜u di . For multi-phase %ow with at least one non-continuous phase, there exists a dispersion %ux for each phase. Because of the diJerence in mobility (or resistance to %ow) for each component, the relative distributions of the diJerent components will be diJerent at diJerent spatial locations. The diJerence in the relative distribution (or concentration) will give rise to “hydrodynamic diJusion” or dispersion. However, if all the phases are continuous, the phase dispersion is diminished. Because of the dispersion, one may get a false sense that the mass is not conserved. One should note that the velocities are not volumetric transport %uxes in Eq. (6). They are convection velocities for the respective phases. There are two components in the transport %ux: convection and dispersion. Eq. (6) is a mass balance equation. The mass transport %ux is given by m ˜ i = i (hi˜u i − ˜u di );
Navier–Stokes equations are partial diJerential equations governing the %ow of a %uid. In particular, the continuity equation for the ith phase is given by, @∗i + · (∗i ˜vi∗ ) = 0; @t
In this case, the porous medium will be treated as a continuous medium and thus no internal porous matrix conEguration and %ow patterns need to be speciEed. The volume-and-time averaged continuity equation is given by
(7)
Finally, the volume-and-time averaged momentum equation is given by hi @i hi˜u i + · (i hi˜u i˜u i ) + ( pi − i˜g) @t − · i [ hi˜u i + ( hi˜u i )T ] + i Fi˜u i − · K ij · (hi i˜u i − hj j˜u j ) j
(5)
Although there are no diJerences between the intrinsic phase averaged value and its corresponding local value for (continuous) %ow in free space, the superscript ∗ in the above equations is used to distinguish them from the averaged quantities in a porous medium. For %ow through porous media, Eqs. (4) and (5) are still applicable within each continuous %uid phase. However, solutions of these equations become formidable due to the complex and often unknown conEgurations of the free voids and the %ow straits of each %uid phase. To add to the complexity, at a given location, the occupant phase can be diJerent at diJerent times. Therefore, volume-and-time averaging is necessary in order for a meaningful solution to be reached in most applications.
− · (i K i · hi˜u i ) +
ij (˜u i − ˜u j ) = 0:
(8)
j
The above averaged Navier–Stokes equation is similar to the averaged heat and mass transfer equations.
3. Volume-and-time averaged Navier–Stokes equations and closure models for gas–liquid ows in packed towers For gas–liquid %ows in packed towers, the following assumptions can be made to establish a simple mathematical model: (1) the packing medium is isotropic and random, (2) liquid wets the solid surface and there are negligible contacts between the gas and the solid, (3) there is no surface force eJect on the %ow behavior. In
S. Liu / Chemical Engineering Science 56 (2001) 5945–5953
particular, liquid is not foaming. With this assumption, the applicability of the following analyses will be limited to %ows below loading point, (4) both gas phase and liquid phase are incompressible. This assumption is valid strictly if the system is isothermal and isobaric, (5) the %uids are Newtonian and the viscosity remains constant in the tower. 3.1. The governing equations The continuity equation for both liquid and gas phases are given by Eq. (6) and the momentum equations can be simpliEed from Eq. (8). hL @L hL˜u L + · (L hL˜u L˜u L ) + ( p − L˜g) @t − · L [ hL˜u L + ( hL˜u L )T ] + L FL˜u L − · [K LG · (L hL˜u L − G hG˜u G )] − · (L K L · hL˜u L ) + G FLG (˜u L − ˜u G ) = 0
(9)
and hG @G hG˜u G + · (G hG˜u G˜u G ) + ( p − G˜g) @t − · G [ hG˜u G + ( hG˜u G )T ]
tion can still be recognized. Liu and Long (2000) gave =
1 − b ∗ + b = b + Er [(1 − 0:3pd ) 2 2 2 R−r × cos + 0:3p d ; 1 + 1:6 Er 2 pd d
(11)
where ∗ is the local azimuthal averaged porosity, pd is the period of oscillation normalized by the characteristic (or nominal) particle size d; Er is the exponential decaying function and b is the bulk porosity or the porosity in an unbounded packing. For spheres of uniform size ds ; pd = 0:94 and d = ds . For cylinders of equal height and diameter such as Pall-rings, d is the √ normal size (diameter or height) and pd = 0:94(2 + 2)=3. The exponential decaying function is given by 3=4 R−r : (12) Er = exp −1:2pd d One can note that Eq. (11) is inadequate for the center region as the convenient condition: d=dr |r=0 = 0 may not hold. The modiEed equation that satisfy the convenient condition is given by =
− · [K LG · (G hG˜u G − L hL˜u L )]
+ G FLG (˜u G − ˜u L ) = 0;
5947
(10)
1 − b ∗ + b = b + Er [(1 − 0:3pd ) 2 2 2 R−r × cos + 0:3pd ; a + 1:6 Er 2 pd d
(13)
where the velocities, ˜u i ’s, are “interstitial” convection. The transport velocities can be deEned similar to Eq. (7). The velocities, hold-ups and pressure can be obtained by solving Eqs. (6), (9) and (10) simultaneously together with no-slip conditions for both the gas and liquid phases on the column wall. However, this requires the complete description of all the parameters involved.
where a is a constant depending on the ratio of the particle size to column size. 3=4 2R R − 1:6 exp −2:4pd a = n pd d d
3.2. Closure model parameters
n = int
In total, there are six parameters in the volume-and-time averaged Eqs. (6), (9) and (10). These parameters are porosity, shear factors, interaction shear factor or inter-phase interaction factor, phase dispersion %ux, momentum dispersion coeOcient tensor and interaction dispersion coeOcient tensor. These parameters need to be deEned (modeled or measured) before one can apply the volume-and-time averaged equations to simulate %ows in packed towers. 3.2.1. Porosity As an averaged quantity, short range variations (on or shorter than the particle scale) need not be considered. The averaging approach actively smoothes out short-range local variations. However, long range varia-
and
R 2 : 1 + 1:6 exp[ − 2:4pd (R=d)3=4 ] pd d
(14)
3.2.2. The shear factor for saturated single phase 4ow The shear factor for saturated single phase %ow has been reviewed by Liu and Masliyah (1996). The full implementation with porosity variations is given by Liu and Masliyah (1999). The shear factor, or %uid–solid matrix interaction term, is given by + b F(ds ; ; |˜u|) = F0 ds ; 2 + F1 (ds ; )
Re2 (Re − 3); Re2 + 36
(15)
where subscript 0 denotes for creeping %ow contribution and I denotes for inertial contribution. Re is the %ow
5948
S. Liu / Chemical Engineering Science 56 (2001) 5945–5953
Assuming the liquid phase is %owing along the surface of the porous matrix, one obtains
hL 1 dp FL0 = − − gx ; (25) L |˜u L | dx L∞
Reynolds number which is deEned as, Re = 24 = 24
Kinetic energy of %ow Channel diameter ×Pressure drop (1=2)u2
e
dp de − d x − gx
= 2
ds |˜u| ;
(16)
∞
where ue is the equivalent velocity of %uid in the passage, de is the equivalent hydraulic passage diameter and the subscript ∞ denotes unidirectional negligible inertial %ow in the channel. ue =
|˜u|
A &
=
|˜u|
(17)
1=6
and de =
where the subscript L∞ denotes “single” liquid phase Elled in the channels formed by the porous matrix. Since the liquid phase is covering all the solid surfaces, the void ratio is exactly the porosity, and the pressure drop can be related to the %ow velocity in a passage through
dp 16k1 ueL − − gx = L ; (26) dx de deL L∞ therefore, FL0 =
4 2ds ; = ap 3(1 − )
(18)
where A is the eJective area ratio of the void space, & is the tortuosity of the passage and ap is the speciEc surface area of the bed. Using the model of Liu and Masliyah (1999), the shear factors are deEned as 18(1 − ) 0:363 0:637 + ; (19) F0 (ds ; ) = 10=3 2 ds s' 18(1 − ) 1 + 0:46(e1=2 ; (20) FI (ds ; ) = 0:048 10=3 2 ds s'
18(1 − hL )hL1=3 : 11=3 d2s
Similar to the shear factor for single phase saturated %ow, the liquid phase shear factor is given by
0:363 1 1=2 + 0:048 + 0:46(eL FL = FL0 0:637 + s'L s'L ReL2 × (28) (ReL − 3) ; ReL2 + 36 where the Reynolds number is deEned similar to saturated single phase %ow ReL = 24
where (e is the passage curvature ratio, (e = 1 − 1=2
s' =
passage length : passage diameter
(22)
The shape factor s' can be found from the relationship of dry pressure drop and FI () for the single gas phase %ow. The relationship is given by, 1 CD G UG2 = G FI (ds ; ∗ ) Re|˜u G |: 2
(23)
3.2.3. Shear factor for the liquid phase For multi-phase %ows, only Attou, Boyer, and Ferschneider (1999) and Long and Liu (2000) attempted a shear factor model. These two papers treat the shear factors in a very similar manner except diJerent model bases were used. In this paper, an extension of the model shown in Section 3.2.2 will be presented for its suitability to packed beds of high performance random packings. The eJective velocity in the passage for the liquid %ow is ueL =
|hL˜u L |
AL &L
= hL1=3
|˜u L |
1=2
:
(24)
(1=2) u2
L eL
Qp deL − L − gx
= 5=3 hL1=3
(21)
and s' is the shape factor which is deEned as
(27)
L∞
L ds |˜u L | : L
(29)
Because of the tortuorous voids, the liquid phase tends to move with a straight and uniform pathway when possible. We propose a compromise such that 1=3 1 − hL s' ; (30) s'L = 1 − hL and
1 − hL (eG = 1 − hL
1=3
(1 − 1=2 ):
(31)
This treatment should allow us to model %ows whereby the liquid phase is not only strictly %owing on the solid surface. Because there is negligible solid matrix–gas phase contact, the shear factor of the gas phase is negligible. 3.2.4. The interaction shear factor between gas phase and liquid phase Since we assume that the gas phase is the non-wetting phase, the gas–solid interaction is negligible, that is, FG = 0 and K G = 0. The interface between the two phases is equivalent to the solid surface for saturated single
S. Liu / Chemical Engineering Science 56 (2001) 5945–5953
phase %ow. Similar to the liquid shear factor, FG0 =
18(1 − hG ) 2 h23=6 G dsG
and
(32)
0:363 1 1=2 FGL = FG0 0:637 + + 0:048 + 0:46(eG s'G s'G ReG2 × − 3) ; (33) (Re G ReG2 + 36
where s'G is the shape fact for the gas phase, (eG is the passage curvature ratio for the gas phase and dsG is the eJective porous matrix grain size for the gas phase %ow (1 − hG )1=3 ds : dsG = (1 − )1=3
(34)
Therefore, the low %ow rate limit interface shear factor is given by 18(1 − hG )1=3 (1 − )2=3 : FG0 = 2 h23=6 G ds
(35)
Using the same deEnition as Eq. (16) and noting the porous matrix grain size for the gas phase and interaction with the liquid phase at the boundary, the gas Reynolds number is given by 1=3 √ ds |G hG˜u G − G L hL˜u L | 1 − hG : ReG = hG 1− G (36) Similar to the shape factor and curvature related terms for the liquid phase shear factor, the passage for the gas phase %ow is not rigid and therefore the modiEcations of these terms must allow for the frequent pulsiEed %ow pattern at the constrictions h = 1 − hG G s' (37) s'G = 1 − hG and
1 − hG (eG = 1 − hG
hG =
(1 − 1=2 ):
(38)
The closure model is introduced based on the promise that the liquid phase %ows along the solid surface and forming an internal connected void in the passage proportional to the void in the porous matrix. ModiEcations are subsequently made to account for the %exible liquid Elm and the alternating %ow pattern at the constrictions. 3.2.5. The liquid dispersion tensor The liquid dispersion is the dispersion of liquid phase Eeld variables due to the interaction of the liquid phase
5949
with the solid matrix. The liquid completely wets the solid matrix and therefore the liquid dispersion can be treated in the same manner as the single phase saturated %ow. The liquid phase dispersion tensor for an isotropic medium is known to be of the following form )L 0 0 )L 0 0 K L = KLrr 0 1 0 = ds DT |˜u L | 0 1 0 ; (39) 0 0 1 0 0 1 where DT is the normalized transverse dispersion coefEcient and )L is the normalized longitudinal dispersion factor. For packed beds of spheres DT = 10:711=3 (1 − )2=3 (1 − 0:75 Er):
(40)
Here the exponential decay function is introduced into the dispersion coeOcient because of the insuOcient averaging to be used for the porosity near the column wall. 3.2.6. The interaction dispersion coe5cient tensor The gas phase is completely surrounded by the liquid phase. The dispersion due to the interaction between the gas and the liquid phases can therefore be treated as the gas phase is %owing through a Ectitious porous matrix, the solids coated with a layer of mobile liquid. That is )L 0 0 K GL = KGLrr 0 1 0 0 0 1 )L 0 0 1=3 1 − hG = ds DT |˜u G − ˜u L | 0 1 0 (41) 1− 0 0 1 and 2=3 DT = 10:7h11=3 G (1 − hG ) (1 − 0:75 Er);
(42)
which is the same as that for saturated single phase %ow with porosity replaced by the holdup. 4. Numerical solutions We consider here fully developed two phase %ows in packed columns. The column is circular cylindrical and the %ows are axi-symmetrical. Only one non-zero averaged velocity component exists: the axial %ow velocity. The cylindrical coordinate system is used with the axial (x-) axis pointing upwards (Fig. 1). All the Eeld variables are dependent only on the radial coordinate. The continuity equations for the liquid and the gas phase are obtained from Eq. (6). We treat the case here the bed porosity is high and pores are large. Liquid not only %ow alongside the solid surface, but forms droplets as
5950
S. Liu / Chemical Engineering Science 56 (2001) 5945–5953
Fig. 2. Normalized superEcial velocity distribution for a saturated single phase %ow through a packed bed of spheres.
hG Fig. 1. A schematic of the reference system for a packed tower.
−
dp G d dhG uG r − G g − dx r dr dr
1 d d(G hG uG − L hL uL ) rKGLrr r dr dr
+ G FGL (uG − uL ) = 0; well. Therefore, there exists a dispersion %ux in the radial direction due to the gradient of stress on the liquid=gas phases similar to particulate dispersion in liquids (Liu, 1999a). Thus, we have FI
L hL − G hG = constant; i:e:; not a function of spatial locations; (43)
L=
L R2
G=
G R2
R
0
0
R
2 rhL uL dr;
(44)
2 rhG uG dr;
(45)
hL + hG =
(46)
and the momentum Eqs. (9) and (10), are reduced to hL dp − L g + L FL uL + G FGL (uL − uG ) dx −
L d dhL uL 1 d d(L hL uL − G hG uG ) r rKGLrr − r dr dr r dr dr
−
1 d dhL uL rL KLrr = 0; r dr dr
(47)
(48)
where uL and uG are the phase intrinsic axial %ow velocities or the gas and the liquid phase, respectively. A visual Basic program has been developed to solve the volume-and-time averaged Eqs. (43) – (48) based on central diJerence method. The boundary conditions used are no slip conditions on the packed column wall. The program was running under the Excel environment. Some of the numerical simulations are presented in Figs. 2–5. Fig. 2 shows the computed velocity distribution for a saturated single phase %ow through a paced bed of spheres. One can observe that the current prediction agrees with the experimental data of Kufner and Hofmann (1990) reasonably well. The locations of the peaks and troughs in the velocity proEle agreed quite well because of the treatment imposed by Eqs. (13) and (14). The diJerence in absolute values is due to the fact that the experimental data is taken a few mm above the actual packings—the peak and trough values have been greatly dampened. The oscillation (appearance of peaks and troughs) in velocity proEle is due to the porosity variation in the radial direction. Therefore, the strength of the current model is that it is capable of predicting the %ow behaviors for beds packed with large particles. In this case, only two particles can be packed on the column radius and the predictions are still reliable.
S. Liu / Chemical Engineering Science 56 (2001) 5945–5953
5951
Fig. 4. Pressure drop and liquid hold-up variation with liquid load for air–water %ows in a column of 0:3 m in diameter packed with 50 mm plastic Nor-PacJ rings.
Fig. 3. Pressure drop and liquid hold-up variation with gas load for air–water %ows in a column of 0:3 m in diameter packed with 50 mm plastic Nor-PacJ rings.
Fig. 3 shows the variation of liquid hold-up and pressure drop for air–water %ow in a typical packed tower (for distillation operations). One can observe that the pressure drop calculation agrees quite well with the experimental data of Billet (1995) for %ow below the loading point. Fig. 4 shows the variation of pressure drop and liquid hold-up for air–water %ow variation with the liquid rate. The pressure drop increases with increasing liquid and=or gas loads. Following the traditional presentation, √ the gas capacity factor, Fv = UG G , is also shown in Fig. 3 alongside the average gas superEcial velocity. As shown in Fig. 3, the pressure drop increases with gas load quite smoothly—the proEle is nearly a straight line on a log–log plane. In contrast, the pressure drop increases with increasing liquid load gradually at low liquid load and sharply at high liquid load. The model predictions for the pressure drop agreed quite well with the experimental data (Fig. 3). However, when the gas load exceeds the %ooding point, the model is incapable of predicting the
sharp increase in pressure drop (near the %ooding point) as the gas load (or gas capacity factor) is increased. The liquid hold-up changes with the gas rate slowly below the %ooding point (Fig. 3). Near the loading point, the liquid hold-up increases sharply with increasing gas load. The liquid hold-up increases with increasing liquid load. The predictions for the liquid hold-up are in agreement with the experimental data as depicted in Fig. 4. Fig. 5 shows that the computed radial liquid velocity proEles agree well with the experimental data of Kouri and Sohlo (1996) for various packings and under various operating conditions. In general, the predicted superEcial liquid velocity proEle is oscillatory in the diameter of the column. On the wall, the no-slip condition is imposed and the liquid velocity is zero. The liquid velocity increases sharply to a maximum just away from the wall and then decreases sharply to a minimum further away from the wall. The minimum liquid velocity is located near the minimum porosity. Therefore, the number of peaks in the liquid velocity proEle depends on the number of periods for the porosity variation along the column radius. For example, the nominal size, d, of a Pall-ring (Figs. 5a–c) is the diameter of the ring and is the same as the length of the ring. The saddles, on the other hand, are nearly semicircular on the base and with a height about one third of its base diameter. The nominal size of a saddle is its diameter and is the longest dimension for the particle. Therefore, there are more peaks in the liquid velocity proEle for packed beds of saddles (Fig. 5d) than for packed beds of Pall-rings with similar nominal particle sizes. The experimental data are piecewise averaged liquid velocity, which is a rather good representation of
5952
S. Liu / Chemical Engineering Science 56 (2001) 5945–5953
Fig. 5. Liquid velocity distributions for gas–liquid %ows in a column of 0:5 m in diameter and packed with random packings 3 m deep: (a) and (b) 25 mm plastic Pallrings; (c) 50 mm plastics Pallrings; and (d) 38 mm ceramic Intalox saddles. The dot–dot–dashed lines are piece-wise averaged experimental data from Kouri and Sohlo (1996).
the computed liquid velocity proEle averaged in the respective regions. They also show a high liquid rate near the wall. 5. Conclusions The volume-and-time averaging technique has been applied to multiphase %ow in porous media=packed beds. It is found that the dispersion term is present even in the continuity equation. As such, the %ow velocities are convective velocities, not the volumetric transport %uxes. The volume-and-time averaged Navier–Stokes equations are highly non-linear since the interaction coeOcients are linearly related to the %ow velocities when the %ow is strong. The dispersion coeOcient and=or the eJective viscosity is also linearly related to the %ow velocities. Models of closure parameters are presented in this paper, including the shear factors, porosity variation and dispersion coeOcient tensors. These closure models are specially designed for packed towers of highly porous random packings. Incorporated with the porosity variations in the radial direction, the model can predict the %ow behaviors even for small column to particle diameter ratios. Several numerical simulations are presented. The %ow velocities, liquid hold-up and pressure drop agree well with experimental data.
Notation A ap a CD DT d de ds Er F F0 FI Fv G g h K k k1 L m ˜ n p
aeral porosity solid surface area per unit bed (column) volume constant in porosity variation expression packing factor constant in the dispersion coeOcient expression nominal (or characteristic) particle size hydraulic pore diameter equivalent spherical particle diameter radial exponential decay function for porosity shear factor or phase–phase momentum interaction factor viscous term of the shear factor inertial term of the shear factor √ gas capacity factor, = UG G gas load or gas mass discharge velocity gravitational acceleration phase hold-up dispersion coeOcient tensor (Kxx ; Kxy ; Kxz ; : : :) permeability passage shape factor: e.g., k1 = 2 for a circular duct and k1 = 3 for parallel plates liquid load or liquid mass discharge velocity superEcial mass velocity integer constant in porosity variation function pressure
S. Liu / Chemical Engineering Science 56 (2001) 5945–5953
pd R Re r s' T t U u udi ue ˜v x y z
period of oscillation in multiple of nominal particle diameter column or bed radius Reynolds number radial coordinate particle shape factor temperature time superEcial %uid discharge velocity, i.e., superEcial %uid velocity averaged on the cross-sectional plane intrinsic %uid phase velocity (local, based on the phase hold-up) dispersion %ux of phase i average velocity in a passage superEcial %uid transport velocity x-coordinate y-coordinate z-coordinate
Greek letters ij )L ∗ . (e &
momentum interaction factor between phase i and phase j ratio of longitudinal dispersion coeOcient to the transversal dispersion coeOcient porosity local porosity or azimuthal averaged porosity azimuthal coordinate eJective pore curvature ratio dynamic viscosity density tortuosity, or the ratio of the straight line distance and the length of eJective %uid pathway between two points
Subscripts b G i j L r x ∞
bulk or far away from the wall gas phase phase i phase j liquid phase r-directional x-directional single phase inertia-free %ow
5953
Acknowledgements The author is indebted to the Natural Sciences and Engineering Research Council of Canada for Enancial support. References Attou, A., Boyer, C., & Ferschneider, G. (1999). Modelling of the hydrodynamics of the cocurrent gas–liquid trickle %ow through a trickle-bed reactor. Chemical Engineering Science, 54, 785–802. Billet, R. (1995). Packed towers in processing and environmental technology. New York, USA: VCH publishers. Dullien, F. A. L. (1992). Porous media, 4uid transport and pore structure (2nd ed.). Toronto: Academic press. Herskowitz, M., & Smith, J. (1983). Trickle-bed reactors. A review. A.I.Ch.E. Journal, 29, 118. Iliuta, I., Ortiz-Arroyo, A., Larachi, F., Grandjean, B. P. A., & Wild, G. (1999). Hydrodynamics and mass transfer in trickle-bed reactors: An overview. Chemical Engineering Science, 54(21), 5329–5337. Kister, H. Z. (1992). Distillation design. New York: MacGraw-Hill. Kouri, R. J., & Sohlo, J. (1996). Liquid and gas %ow patterns in random packings. Chemical Engineering Journal, 61, 95–105. Kufner, R., & Hofmann, H. (1990). Implementation of radial porosity and velocity distribution in a reactor model for heterogeneous catalytic gas phase reactions (Torus-Model). Chemical Engineering Science, 45, 2141–2146. Liu, S. (1999a). Particle dispersion for suspension %ow. Chemical Engineering Science, 54, 873–891. Liu, S. (1999b). A continuum approach to multi-phase %ow in porous media. Journal of Porous Media, 2(3), 295–308. Liu, S., Afacan, A., Nasr-El-Din, H. A., & Masliyah, J. H. (1994). An experimental study of pressure drop in helical pipes. Proceedings of the Royal Society London A, 444, 307–316. Liu, S., & Long, J. (2000). Gas–liquid countercurrent %ows through packed towers. Journal of Porous Media, 3(2), 99–113. Liu, S., & Masliyah, J. H. (1996). Single %uid %ow in porous media. Chemical Engineering Communication, 148–150, 653–732. Liu, S., & Masliyah, J. H. (1999). Non-linear %ows in porous media. Journal Non-Newtonian Fluid Mechanics, 86, 229–252. Muskat, M. (1937). The 4ow of homogeneous 4uids through porous media. New york: McGraw-Hill. Richards, L. A. (1931). Capillary conduction of liquids through porous media. Physics, 1, 318–333. Saroha, A. K., Nigam, K. D. P., Saxena, A. K., & Kapoor, V. K. (1998). Liquid distribution in trickle-bed reactors. A.I.Ch.E. Journal, 44(9), 2044–2052. Shah, Y. T. (1979). Gas–liquid–solid reactor design. New York: McGraw-Hill. Trebal, R. E. (1980). Mass transfer operations (3rd ed.). Toronto: McGraw-Hill.
www.elsevier.com/locate/ces
A continuum model for gas–liquid %ow in packed towers Shijie Liu ∗ Department of Chemical and Materials Engineering, University of Alberta, Edmonton AB, Canada T6G 2G6
Abstract Gas–liquid %ow in packed towers is commonly encountered in the chemical and processing industry. A continuum model is developed based on the volume-and-time averaging of multiphase %ows in isotropic rigid porous media=packed columns. Closures are presented for the evaluations of the extra surface=intrinsic phase integral terms. Both inertia and inter-phase interactions are retained in the volume averaged (Navier–Stokes) equations. These governing equations are solved for fully-developed axi-symmetric single and gas–liquid two phase %ows in highly porous packed towers. It is found that the dispersion term is present in the continuity equation as well as the momentum equations. Numerical simulations with the models show that the volume-and-time averaged equations can predict the velocity, phase hold-up and pressure drop quite well for up to the loading point for gas–liquid counter-current %ows. ? 2001 Elsevier Science Ltd. All rights reserved. Keywords: Packed tower; Multiphase %ow; Volume-and-time averaging; Counter-current %ow; Velocity distribution; Phase hold-up
1. Introduction There are three gas–liquid %ow arrangements useful for separation and reaction engineering applications. There is only one possible way for countercurrent %ow arrangement with the heavy %uid (liquid) %owing downward and the light %uid (gas) %owing upward. This countercurrent %ow arrangement is commonly used for distillation, cooling tower, absorption and stripping operations (e.g., Kister, 1992; Billet, 1995). There are two possible arrangements for cocurrent %ows: either upward or downward %ow. While concurrent %ows can be applied for mass and=or heat transfer operations, they are generally used as multiphase reactors (e.g., Shah, 1979; Herskowitz & Smith, 1983). In particular, the downward %ow arrangement is commonly used. In this case, the liquid rate can be small and the unique %ow system is called a trickle bed (e.g., Saroha, Nigam, Saxena, & Kapoor, 1998; Iliuta, Ortiz-Arroyo, Larachi, Grandjean, & Wild, 1999). The cocurrent up %ow requires higher pressure drop and liquid hold-up in tower than the down %ow mode. A continuum model for numerical simulation of %ows in such systems will be developed in this paper. Closure models will be presented suitable for gas–liquid countercurrent %ows. ∗
Tel.: +1-780-492-0981; fax: +1-780-492-2881. E-mail address: [email protected] (S. Liu).
For a laminar saturated single %uid %ow in a straight duct, the pressure drop is proportional to the %ow rate. The same relationship holds also for %ow in curved ducts when the %ow is weak (Liu, Afacan, Nasr-El-Din, & Masliyah, 1994). This unique relationship can be generalized to %ow through porous media as well, dp − gx ˙ −U; dx
(1)
where U is the superEcial %uid velocity, p is pressure, x is the axis of the %ow direction, is the %uid density and gx is the gravity in the direction of %ow. This limiting %ow behavior has received great attention for centuries. The direct expression resulted from the proportionality for %ow in porous media is called the Darcy’s law. Several approaches have been used for the semitheoretical modeling of packed towers. When the %ow velocity is high, Eq. (1) is not strictly adequate. In this case, all the existing models converge to a simple (and experimentally observed) relationship, which is commonly expressed as −
1 dp + G gx = CD G UG2 ; dx 2
(2)
where CD is the so-called packing factor and has a reciprocal length unit. It is inversely proportional to the particle size and is a strong function of the packing
0009-2509/01/$ - see front matter ? 2001 Elsevier Science Ltd. All rights reserved. PII: S 0 0 0 9 - 2 5 0 9 ( 0 1 ) 0 0 2 3 2 - 9
5946
S. Liu / Chemical Engineering Science 56 (2001) 5945–5953
porosity as well. The packing factor is normally supplied by manufactures for each speciEc random packings. Because the expression is usually used for dry gas %ow through packed beds, a subscript G was used in Eq. (2). A collection of data for the packing factor is also available in textbooks, for example, Trebal (1980). Richards (1931) and Muskat (1937) intuitively extended the Darcy’s law to multiphase %ows. In their extension, the equation remains the same for each phase but allowing the %uid properties as well as the permeability to diJer. That is, ˜vi = −
ki ( pi − i˜g); i
(3)
where the subscript i denotes the ith %uid phase, k is the permeability and is the dynamic viscosity. Although widely used for %ow in soils and sandstones, Eq. (3) lacks three key features (Liu, 1999b): namely, inertial eJects; diJusion eJects and the inter-phase interactions (Dullien, 1992). These defects make Eq. (3) inapplicable to the gas–liquid %ows in packed towers. In this paper, we shall apply the volume-and-time averaging principles to derive a governing equation in packed towers.
2. Volume-and-time averaged Navier–Stokes equations for multiphase ows in porous media
(4)
and the ith phase momentum equation is given by @(∗i ˜vi∗ ) + · (∗i ˜vi∗˜vi∗ ) + p∗ − ∗i ˜g @t − · i∗ [˜vi∗ + (˜vi∗ )T ] = 0:
@(i hi ) + · (i hi˜u i ) − · (i˜u di ) = 0: @t
(6)
Eq. (6) is similar to the traditional phase balance (or continuity) equation except for the dispersion %ux ˜u di . For multi-phase %ow with at least one non-continuous phase, there exists a dispersion %ux for each phase. Because of the diJerence in mobility (or resistance to %ow) for each component, the relative distributions of the diJerent components will be diJerent at diJerent spatial locations. The diJerence in the relative distribution (or concentration) will give rise to “hydrodynamic diJusion” or dispersion. However, if all the phases are continuous, the phase dispersion is diminished. Because of the dispersion, one may get a false sense that the mass is not conserved. One should note that the velocities are not volumetric transport %uxes in Eq. (6). They are convection velocities for the respective phases. There are two components in the transport %ux: convection and dispersion. Eq. (6) is a mass balance equation. The mass transport %ux is given by m ˜ i = i (hi˜u i − ˜u di );
Navier–Stokes equations are partial diJerential equations governing the %ow of a %uid. In particular, the continuity equation for the ith phase is given by, @∗i + · (∗i ˜vi∗ ) = 0; @t
In this case, the porous medium will be treated as a continuous medium and thus no internal porous matrix conEguration and %ow patterns need to be speciEed. The volume-and-time averaged continuity equation is given by
(7)
Finally, the volume-and-time averaged momentum equation is given by hi @i hi˜u i + · (i hi˜u i˜u i ) + ( pi − i˜g) @t − · i [ hi˜u i + ( hi˜u i )T ] + i Fi˜u i − · K ij · (hi i˜u i − hj j˜u j ) j
(5)
Although there are no diJerences between the intrinsic phase averaged value and its corresponding local value for (continuous) %ow in free space, the superscript ∗ in the above equations is used to distinguish them from the averaged quantities in a porous medium. For %ow through porous media, Eqs. (4) and (5) are still applicable within each continuous %uid phase. However, solutions of these equations become formidable due to the complex and often unknown conEgurations of the free voids and the %ow straits of each %uid phase. To add to the complexity, at a given location, the occupant phase can be diJerent at diJerent times. Therefore, volume-and-time averaging is necessary in order for a meaningful solution to be reached in most applications.
− · (i K i · hi˜u i ) +
ij (˜u i − ˜u j ) = 0:
(8)
j
The above averaged Navier–Stokes equation is similar to the averaged heat and mass transfer equations.
3. Volume-and-time averaged Navier–Stokes equations and closure models for gas–liquid ows in packed towers For gas–liquid %ows in packed towers, the following assumptions can be made to establish a simple mathematical model: (1) the packing medium is isotropic and random, (2) liquid wets the solid surface and there are negligible contacts between the gas and the solid, (3) there is no surface force eJect on the %ow behavior. In
S. Liu / Chemical Engineering Science 56 (2001) 5945–5953
particular, liquid is not foaming. With this assumption, the applicability of the following analyses will be limited to %ows below loading point, (4) both gas phase and liquid phase are incompressible. This assumption is valid strictly if the system is isothermal and isobaric, (5) the %uids are Newtonian and the viscosity remains constant in the tower. 3.1. The governing equations The continuity equation for both liquid and gas phases are given by Eq. (6) and the momentum equations can be simpliEed from Eq. (8). hL @L hL˜u L + · (L hL˜u L˜u L ) + ( p − L˜g) @t − · L [ hL˜u L + ( hL˜u L )T ] + L FL˜u L − · [K LG · (L hL˜u L − G hG˜u G )] − · (L K L · hL˜u L ) + G FLG (˜u L − ˜u G ) = 0
(9)
and hG @G hG˜u G + · (G hG˜u G˜u G ) + ( p − G˜g) @t − · G [ hG˜u G + ( hG˜u G )T ]
tion can still be recognized. Liu and Long (2000) gave =
1 − b ∗ + b = b + Er [(1 − 0:3pd ) 2 2 2 R−r × cos + 0:3p d ; 1 + 1:6 Er 2 pd d
(11)
where ∗ is the local azimuthal averaged porosity, pd is the period of oscillation normalized by the characteristic (or nominal) particle size d; Er is the exponential decaying function and b is the bulk porosity or the porosity in an unbounded packing. For spheres of uniform size ds ; pd = 0:94 and d = ds . For cylinders of equal height and diameter such as Pall-rings, d is the √ normal size (diameter or height) and pd = 0:94(2 + 2)=3. The exponential decaying function is given by 3=4 R−r : (12) Er = exp −1:2pd d One can note that Eq. (11) is inadequate for the center region as the convenient condition: d=dr |r=0 = 0 may not hold. The modiEed equation that satisfy the convenient condition is given by =
− · [K LG · (G hG˜u G − L hL˜u L )]
+ G FLG (˜u G − ˜u L ) = 0;
5947
(10)
1 − b ∗ + b = b + Er [(1 − 0:3pd ) 2 2 2 R−r × cos + 0:3pd ; a + 1:6 Er 2 pd d
(13)
where the velocities, ˜u i ’s, are “interstitial” convection. The transport velocities can be deEned similar to Eq. (7). The velocities, hold-ups and pressure can be obtained by solving Eqs. (6), (9) and (10) simultaneously together with no-slip conditions for both the gas and liquid phases on the column wall. However, this requires the complete description of all the parameters involved.
where a is a constant depending on the ratio of the particle size to column size. 3=4 2R R − 1:6 exp −2:4pd a = n pd d d
3.2. Closure model parameters
n = int
In total, there are six parameters in the volume-and-time averaged Eqs. (6), (9) and (10). These parameters are porosity, shear factors, interaction shear factor or inter-phase interaction factor, phase dispersion %ux, momentum dispersion coeOcient tensor and interaction dispersion coeOcient tensor. These parameters need to be deEned (modeled or measured) before one can apply the volume-and-time averaged equations to simulate %ows in packed towers. 3.2.1. Porosity As an averaged quantity, short range variations (on or shorter than the particle scale) need not be considered. The averaging approach actively smoothes out short-range local variations. However, long range varia-
and
R 2 : 1 + 1:6 exp[ − 2:4pd (R=d)3=4 ] pd d
(14)
3.2.2. The shear factor for saturated single phase 4ow The shear factor for saturated single phase %ow has been reviewed by Liu and Masliyah (1996). The full implementation with porosity variations is given by Liu and Masliyah (1999). The shear factor, or %uid–solid matrix interaction term, is given by + b F(ds ; ; |˜u|) = F0 ds ; 2 + F1 (ds ; )
Re2 (Re − 3); Re2 + 36
(15)
where subscript 0 denotes for creeping %ow contribution and I denotes for inertial contribution. Re is the %ow
5948
S. Liu / Chemical Engineering Science 56 (2001) 5945–5953
Assuming the liquid phase is %owing along the surface of the porous matrix, one obtains
hL 1 dp FL0 = − − gx ; (25) L |˜u L | dx L∞
Reynolds number which is deEned as, Re = 24 = 24
Kinetic energy of %ow Channel diameter ×Pressure drop (1=2)u2
e
dp de − d x − gx
= 2
ds |˜u| ;
(16)
∞
where ue is the equivalent velocity of %uid in the passage, de is the equivalent hydraulic passage diameter and the subscript ∞ denotes unidirectional negligible inertial %ow in the channel. ue =
|˜u|
A &
=
|˜u|
(17)
1=6
and de =
where the subscript L∞ denotes “single” liquid phase Elled in the channels formed by the porous matrix. Since the liquid phase is covering all the solid surfaces, the void ratio is exactly the porosity, and the pressure drop can be related to the %ow velocity in a passage through
dp 16k1 ueL − − gx = L ; (26) dx de deL L∞ therefore, FL0 =
4 2ds ; = ap 3(1 − )
(18)
where A is the eJective area ratio of the void space, & is the tortuosity of the passage and ap is the speciEc surface area of the bed. Using the model of Liu and Masliyah (1999), the shear factors are deEned as 18(1 − ) 0:363 0:637 + ; (19) F0 (ds ; ) = 10=3 2 ds s' 18(1 − ) 1 + 0:46(e1=2 ; (20) FI (ds ; ) = 0:048 10=3 2 ds s'
18(1 − hL )hL1=3 : 11=3 d2s
Similar to the shear factor for single phase saturated %ow, the liquid phase shear factor is given by
0:363 1 1=2 + 0:048 + 0:46(eL FL = FL0 0:637 + s'L s'L ReL2 × (28) (ReL − 3) ; ReL2 + 36 where the Reynolds number is deEned similar to saturated single phase %ow ReL = 24
where (e is the passage curvature ratio, (e = 1 − 1=2
s' =
passage length : passage diameter
(22)
The shape factor s' can be found from the relationship of dry pressure drop and FI () for the single gas phase %ow. The relationship is given by, 1 CD G UG2 = G FI (ds ; ∗ ) Re|˜u G |: 2
(23)
3.2.3. Shear factor for the liquid phase For multi-phase %ows, only Attou, Boyer, and Ferschneider (1999) and Long and Liu (2000) attempted a shear factor model. These two papers treat the shear factors in a very similar manner except diJerent model bases were used. In this paper, an extension of the model shown in Section 3.2.2 will be presented for its suitability to packed beds of high performance random packings. The eJective velocity in the passage for the liquid %ow is ueL =
|hL˜u L |
AL &L
= hL1=3
|˜u L |
1=2
:
(24)
(1=2) u2
L eL
Qp deL − L − gx
= 5=3 hL1=3
(21)
and s' is the shape factor which is deEned as
(27)
L∞
L ds |˜u L | : L
(29)
Because of the tortuorous voids, the liquid phase tends to move with a straight and uniform pathway when possible. We propose a compromise such that 1=3 1 − hL s' ; (30) s'L = 1 − hL and
1 − hL (eG = 1 − hL
1=3
(1 − 1=2 ):
(31)
This treatment should allow us to model %ows whereby the liquid phase is not only strictly %owing on the solid surface. Because there is negligible solid matrix–gas phase contact, the shear factor of the gas phase is negligible. 3.2.4. The interaction shear factor between gas phase and liquid phase Since we assume that the gas phase is the non-wetting phase, the gas–solid interaction is negligible, that is, FG = 0 and K G = 0. The interface between the two phases is equivalent to the solid surface for saturated single
S. Liu / Chemical Engineering Science 56 (2001) 5945–5953
phase %ow. Similar to the liquid shear factor, FG0 =
18(1 − hG ) 2 h23=6 G dsG
and
(32)
0:363 1 1=2 FGL = FG0 0:637 + + 0:048 + 0:46(eG s'G s'G ReG2 × − 3) ; (33) (Re G ReG2 + 36
where s'G is the shape fact for the gas phase, (eG is the passage curvature ratio for the gas phase and dsG is the eJective porous matrix grain size for the gas phase %ow (1 − hG )1=3 ds : dsG = (1 − )1=3
(34)
Therefore, the low %ow rate limit interface shear factor is given by 18(1 − hG )1=3 (1 − )2=3 : FG0 = 2 h23=6 G ds
(35)
Using the same deEnition as Eq. (16) and noting the porous matrix grain size for the gas phase and interaction with the liquid phase at the boundary, the gas Reynolds number is given by 1=3 √ ds |G hG˜u G − G L hL˜u L | 1 − hG : ReG = hG 1− G (36) Similar to the shape factor and curvature related terms for the liquid phase shear factor, the passage for the gas phase %ow is not rigid and therefore the modiEcations of these terms must allow for the frequent pulsiEed %ow pattern at the constrictions h = 1 − hG G s' (37) s'G = 1 − hG and
1 − hG (eG = 1 − hG
hG =
(1 − 1=2 ):
(38)
The closure model is introduced based on the promise that the liquid phase %ows along the solid surface and forming an internal connected void in the passage proportional to the void in the porous matrix. ModiEcations are subsequently made to account for the %exible liquid Elm and the alternating %ow pattern at the constrictions. 3.2.5. The liquid dispersion tensor The liquid dispersion is the dispersion of liquid phase Eeld variables due to the interaction of the liquid phase
5949
with the solid matrix. The liquid completely wets the solid matrix and therefore the liquid dispersion can be treated in the same manner as the single phase saturated %ow. The liquid phase dispersion tensor for an isotropic medium is known to be of the following form )L 0 0 )L 0 0 K L = KLrr 0 1 0 = ds DT |˜u L | 0 1 0 ; (39) 0 0 1 0 0 1 where DT is the normalized transverse dispersion coefEcient and )L is the normalized longitudinal dispersion factor. For packed beds of spheres DT = 10:711=3 (1 − )2=3 (1 − 0:75 Er):
(40)
Here the exponential decay function is introduced into the dispersion coeOcient because of the insuOcient averaging to be used for the porosity near the column wall. 3.2.6. The interaction dispersion coe5cient tensor The gas phase is completely surrounded by the liquid phase. The dispersion due to the interaction between the gas and the liquid phases can therefore be treated as the gas phase is %owing through a Ectitious porous matrix, the solids coated with a layer of mobile liquid. That is )L 0 0 K GL = KGLrr 0 1 0 0 0 1 )L 0 0 1=3 1 − hG = ds DT |˜u G − ˜u L | 0 1 0 (41) 1− 0 0 1 and 2=3 DT = 10:7h11=3 G (1 − hG ) (1 − 0:75 Er);
(42)
which is the same as that for saturated single phase %ow with porosity replaced by the holdup. 4. Numerical solutions We consider here fully developed two phase %ows in packed columns. The column is circular cylindrical and the %ows are axi-symmetrical. Only one non-zero averaged velocity component exists: the axial %ow velocity. The cylindrical coordinate system is used with the axial (x-) axis pointing upwards (Fig. 1). All the Eeld variables are dependent only on the radial coordinate. The continuity equations for the liquid and the gas phase are obtained from Eq. (6). We treat the case here the bed porosity is high and pores are large. Liquid not only %ow alongside the solid surface, but forms droplets as
5950
S. Liu / Chemical Engineering Science 56 (2001) 5945–5953
Fig. 2. Normalized superEcial velocity distribution for a saturated single phase %ow through a packed bed of spheres.
hG Fig. 1. A schematic of the reference system for a packed tower.
−
dp G d dhG uG r − G g − dx r dr dr
1 d d(G hG uG − L hL uL ) rKGLrr r dr dr
+ G FGL (uG − uL ) = 0; well. Therefore, there exists a dispersion %ux in the radial direction due to the gradient of stress on the liquid=gas phases similar to particulate dispersion in liquids (Liu, 1999a). Thus, we have FI
L hL − G hG = constant; i:e:; not a function of spatial locations; (43)
L=
L R2
G=
G R2
R
0
0
R
2 rhL uL dr;
(44)
2 rhG uG dr;
(45)
hL + hG =
(46)
and the momentum Eqs. (9) and (10), are reduced to hL dp − L g + L FL uL + G FGL (uL − uG ) dx −
L d dhL uL 1 d d(L hL uL − G hG uG ) r rKGLrr − r dr dr r dr dr
−
1 d dhL uL rL KLrr = 0; r dr dr
(47)
(48)
where uL and uG are the phase intrinsic axial %ow velocities or the gas and the liquid phase, respectively. A visual Basic program has been developed to solve the volume-and-time averaged Eqs. (43) – (48) based on central diJerence method. The boundary conditions used are no slip conditions on the packed column wall. The program was running under the Excel environment. Some of the numerical simulations are presented in Figs. 2–5. Fig. 2 shows the computed velocity distribution for a saturated single phase %ow through a paced bed of spheres. One can observe that the current prediction agrees with the experimental data of Kufner and Hofmann (1990) reasonably well. The locations of the peaks and troughs in the velocity proEle agreed quite well because of the treatment imposed by Eqs. (13) and (14). The diJerence in absolute values is due to the fact that the experimental data is taken a few mm above the actual packings—the peak and trough values have been greatly dampened. The oscillation (appearance of peaks and troughs) in velocity proEle is due to the porosity variation in the radial direction. Therefore, the strength of the current model is that it is capable of predicting the %ow behaviors for beds packed with large particles. In this case, only two particles can be packed on the column radius and the predictions are still reliable.
S. Liu / Chemical Engineering Science 56 (2001) 5945–5953
5951
Fig. 4. Pressure drop and liquid hold-up variation with liquid load for air–water %ows in a column of 0:3 m in diameter packed with 50 mm plastic Nor-PacJ rings.
Fig. 3. Pressure drop and liquid hold-up variation with gas load for air–water %ows in a column of 0:3 m in diameter packed with 50 mm plastic Nor-PacJ rings.
Fig. 3 shows the variation of liquid hold-up and pressure drop for air–water %ow in a typical packed tower (for distillation operations). One can observe that the pressure drop calculation agrees quite well with the experimental data of Billet (1995) for %ow below the loading point. Fig. 4 shows the variation of pressure drop and liquid hold-up for air–water %ow variation with the liquid rate. The pressure drop increases with increasing liquid and=or gas loads. Following the traditional presentation, √ the gas capacity factor, Fv = UG G , is also shown in Fig. 3 alongside the average gas superEcial velocity. As shown in Fig. 3, the pressure drop increases with gas load quite smoothly—the proEle is nearly a straight line on a log–log plane. In contrast, the pressure drop increases with increasing liquid load gradually at low liquid load and sharply at high liquid load. The model predictions for the pressure drop agreed quite well with the experimental data (Fig. 3). However, when the gas load exceeds the %ooding point, the model is incapable of predicting the
sharp increase in pressure drop (near the %ooding point) as the gas load (or gas capacity factor) is increased. The liquid hold-up changes with the gas rate slowly below the %ooding point (Fig. 3). Near the loading point, the liquid hold-up increases sharply with increasing gas load. The liquid hold-up increases with increasing liquid load. The predictions for the liquid hold-up are in agreement with the experimental data as depicted in Fig. 4. Fig. 5 shows that the computed radial liquid velocity proEles agree well with the experimental data of Kouri and Sohlo (1996) for various packings and under various operating conditions. In general, the predicted superEcial liquid velocity proEle is oscillatory in the diameter of the column. On the wall, the no-slip condition is imposed and the liquid velocity is zero. The liquid velocity increases sharply to a maximum just away from the wall and then decreases sharply to a minimum further away from the wall. The minimum liquid velocity is located near the minimum porosity. Therefore, the number of peaks in the liquid velocity proEle depends on the number of periods for the porosity variation along the column radius. For example, the nominal size, d, of a Pall-ring (Figs. 5a–c) is the diameter of the ring and is the same as the length of the ring. The saddles, on the other hand, are nearly semicircular on the base and with a height about one third of its base diameter. The nominal size of a saddle is its diameter and is the longest dimension for the particle. Therefore, there are more peaks in the liquid velocity proEle for packed beds of saddles (Fig. 5d) than for packed beds of Pall-rings with similar nominal particle sizes. The experimental data are piecewise averaged liquid velocity, which is a rather good representation of
5952
S. Liu / Chemical Engineering Science 56 (2001) 5945–5953
Fig. 5. Liquid velocity distributions for gas–liquid %ows in a column of 0:5 m in diameter and packed with random packings 3 m deep: (a) and (b) 25 mm plastic Pallrings; (c) 50 mm plastics Pallrings; and (d) 38 mm ceramic Intalox saddles. The dot–dot–dashed lines are piece-wise averaged experimental data from Kouri and Sohlo (1996).
the computed liquid velocity proEle averaged in the respective regions. They also show a high liquid rate near the wall. 5. Conclusions The volume-and-time averaging technique has been applied to multiphase %ow in porous media=packed beds. It is found that the dispersion term is present even in the continuity equation. As such, the %ow velocities are convective velocities, not the volumetric transport %uxes. The volume-and-time averaged Navier–Stokes equations are highly non-linear since the interaction coeOcients are linearly related to the %ow velocities when the %ow is strong. The dispersion coeOcient and=or the eJective viscosity is also linearly related to the %ow velocities. Models of closure parameters are presented in this paper, including the shear factors, porosity variation and dispersion coeOcient tensors. These closure models are specially designed for packed towers of highly porous random packings. Incorporated with the porosity variations in the radial direction, the model can predict the %ow behaviors even for small column to particle diameter ratios. Several numerical simulations are presented. The %ow velocities, liquid hold-up and pressure drop agree well with experimental data.
Notation A ap a CD DT d de ds Er F F0 FI Fv G g h K k k1 L m ˜ n p
aeral porosity solid surface area per unit bed (column) volume constant in porosity variation expression packing factor constant in the dispersion coeOcient expression nominal (or characteristic) particle size hydraulic pore diameter equivalent spherical particle diameter radial exponential decay function for porosity shear factor or phase–phase momentum interaction factor viscous term of the shear factor inertial term of the shear factor √ gas capacity factor, = UG G gas load or gas mass discharge velocity gravitational acceleration phase hold-up dispersion coeOcient tensor (Kxx ; Kxy ; Kxz ; : : :) permeability passage shape factor: e.g., k1 = 2 for a circular duct and k1 = 3 for parallel plates liquid load or liquid mass discharge velocity superEcial mass velocity integer constant in porosity variation function pressure
S. Liu / Chemical Engineering Science 56 (2001) 5945–5953
pd R Re r s' T t U u udi ue ˜v x y z
period of oscillation in multiple of nominal particle diameter column or bed radius Reynolds number radial coordinate particle shape factor temperature time superEcial %uid discharge velocity, i.e., superEcial %uid velocity averaged on the cross-sectional plane intrinsic %uid phase velocity (local, based on the phase hold-up) dispersion %ux of phase i average velocity in a passage superEcial %uid transport velocity x-coordinate y-coordinate z-coordinate
Greek letters ij )L ∗ . (e &
momentum interaction factor between phase i and phase j ratio of longitudinal dispersion coeOcient to the transversal dispersion coeOcient porosity local porosity or azimuthal averaged porosity azimuthal coordinate eJective pore curvature ratio dynamic viscosity density tortuosity, or the ratio of the straight line distance and the length of eJective %uid pathway between two points
Subscripts b G i j L r x ∞
bulk or far away from the wall gas phase phase i phase j liquid phase r-directional x-directional single phase inertia-free %ow
5953
Acknowledgements The author is indebted to the Natural Sciences and Engineering Research Council of Canada for Enancial support. References Attou, A., Boyer, C., & Ferschneider, G. (1999). Modelling of the hydrodynamics of the cocurrent gas–liquid trickle %ow through a trickle-bed reactor. Chemical Engineering Science, 54, 785–802. Billet, R. (1995). Packed towers in processing and environmental technology. New York, USA: VCH publishers. Dullien, F. A. L. (1992). Porous media, 4uid transport and pore structure (2nd ed.). Toronto: Academic press. Herskowitz, M., & Smith, J. (1983). Trickle-bed reactors. A review. A.I.Ch.E. Journal, 29, 118. Iliuta, I., Ortiz-Arroyo, A., Larachi, F., Grandjean, B. P. A., & Wild, G. (1999). Hydrodynamics and mass transfer in trickle-bed reactors: An overview. Chemical Engineering Science, 54(21), 5329–5337. Kister, H. Z. (1992). Distillation design. New York: MacGraw-Hill. Kouri, R. J., & Sohlo, J. (1996). Liquid and gas %ow patterns in random packings. Chemical Engineering Journal, 61, 95–105. Kufner, R., & Hofmann, H. (1990). Implementation of radial porosity and velocity distribution in a reactor model for heterogeneous catalytic gas phase reactions (Torus-Model). Chemical Engineering Science, 45, 2141–2146. Liu, S. (1999a). Particle dispersion for suspension %ow. Chemical Engineering Science, 54, 873–891. Liu, S. (1999b). A continuum approach to multi-phase %ow in porous media. Journal of Porous Media, 2(3), 295–308. Liu, S., Afacan, A., Nasr-El-Din, H. A., & Masliyah, J. H. (1994). An experimental study of pressure drop in helical pipes. Proceedings of the Royal Society London A, 444, 307–316. Liu, S., & Long, J. (2000). Gas–liquid countercurrent %ows through packed towers. Journal of Porous Media, 3(2), 99–113. Liu, S., & Masliyah, J. H. (1996). Single %uid %ow in porous media. Chemical Engineering Communication, 148–150, 653–732. Liu, S., & Masliyah, J. H. (1999). Non-linear %ows in porous media. Journal Non-Newtonian Fluid Mechanics, 86, 229–252. Muskat, M. (1937). The 4ow of homogeneous 4uids through porous media. New york: McGraw-Hill. Richards, L. A. (1931). Capillary conduction of liquids through porous media. Physics, 1, 318–333. Saroha, A. K., Nigam, K. D. P., Saxena, A. K., & Kapoor, V. K. (1998). Liquid distribution in trickle-bed reactors. A.I.Ch.E. Journal, 44(9), 2044–2052. Shah, Y. T. (1979). Gas–liquid–solid reactor design. New York: McGraw-Hill. Trebal, R. E. (1980). Mass transfer operations (3rd ed.). Toronto: McGraw-Hill.