Anal. Chem. 2003, 75, 843-850
A Correlation for the Pressure Drop in Monolithic Silica Columns Nico Vervoort,* Piotr Gzil, Gino V. Baron, and Gert Desmet
Vrije Universiteit Brussel, Department of Chemical Engineering, Pleinlaan 2, 1050 Brussels, Belgium
To gain insight into how the pressure drop in monolithic silica columns is determined by the microscopic details of the pore structure, a series of well-validated computational fluid dynamics simulations has been performed on a simplified model structure, the so-called tetrahedral skeleton column. From these simulations, a direct correlation between the pressure drop and two main structural properties (skeleton thickness and column porosity) of the monolithic skeleton could be established. The correlation shows good agreement with the experimental pressure-drop data available from the literature on silica monoliths, especially when a correction for the flowthrough pore size heterogeneity is made. The established correlation also yields a much more accurate representation of the relation between the flow resistance and the bed porosity than does the Kozeny-Carman model, making it much better suited for porosity optimization calculations. Ever since their first introduction in the field of chromatography,1-3 monolithic columns have been acclaimed for their superior flow characteristics. It was assumed that, because of the possibility to control the size of the flow-through pores independently of the size of the solid-phase skeleton, the use of wider flow-through pores would result in a marked smaller flow resistance and eddy-diffusion than the packed bed of spheres. A vast number of studies4-13 has since then clearly confirmed this * Corresponding author. Phone: (+)0 32 2 629 36 17. Fax.: (+)0 32 2 629 32 48. E-mail:
[email protected]. (1) Hjerte´n, S.; Liao, J.-L.; Zhang, R. J. Chromatogr. 1989, 473, 273-275. (2) Liao, L. J.; Hjerten, S. J. Chromatogr. 1988, 457, 165-174. (3) Svec, F.; Frechet, J. M. Anal. Chem. 1992, 64, 820-822. (4) Minakuchi, H.; Nakanishi, K.; Soga, N.; Ishizuka, N.; Tanaka, N. Anal. Chem. 1996, 68, 3498-3501. (5) Minakuchi, H.; Ishizuka, N.; Nakanishi, K.; Soga, N.; Tanaka, N. J. Chromatogr. 1998, 828, 83-90. (6) Minakuchi, H.; Nakanishi, K.; Soga, N.; Ishizuka, N.; Tanaka, N. J. Chromatogr. 1998, 797, 121-131. (7) Minakuchi, H.; Nakanishi, K.; Soga, N.; Ishizuka, N.; Tanaka, N. J. Chromatogr. 1997, 762, 135-146. (8) Minakuchi, H.; Nakanishi, K.; Soga, N.; Ishizuka, N.; Tanaka, N. Anal. Chem. 1996, 68, 1275-1280. (9) Ishizuka, N.; Minakuchi, H.; Nakanishi, K.; Soga, N.; Nagayama, H.; Hosoya, K.; Tanaka, N. Anal. Chem. 2000, 72, 1275-1280. (10) Ishizuka, N.; Kobayashi, H.; Minakuchi, H.; Nakanishi, K.; Hirao, K.; Hosoya, K.; Tanaka, N.; Ikegami T.; Tanaka N. Anal. Chem. 2002, 960, 85-96. (11) Tanaka, N.; Nagayama, H.; Kobayashi, H.; Ikegami, T.; Hosoya, K.; Ishizuka, N.; Minakuchi, H.; Nakanishi, K.; Cabrera, K.; Lubda, P. J. High Resolut. Chromatogr. 2000, 23, 111-116. (12) Tennikova, T. B.; Belenkii, B. G.; Sˇ vec, F. J. Liq. Chromatogr. 1990, 13, 63. (13) Sˇ trancar, A.; Koselj, P.; Scwinn, H.; Josic´, Dj. Anal. Chem. 1996, 68, 3483. 10.1021/ac0262199 CCC: $25.00 Published on Web 01/15/2003
© 2003 American Chemical Society
intuitive idea. The advantage of monolithic columns is most probably best illustrated by the dramatic reduction of the separation impedance, E.14 Whereas a typical packed sphere column has an impedance on the order of 3000, it has been demonstrated by several groups that the separation impedance of silica monolithic columns can be made as small as a few hundreds.15,16 With the advent of the monolithic silica and polymeric columns,4-13 the monolithic column concept has now matured into a number of successful commercial products. The flow properties of the monoliths are, however, still poorly understood.18 For example, a model linking the observed pressure drop to the specific pore structure of silica monolithic columns is still lacking. Traditionally, the flow resistance of a monolithic column is modeled by comparing the observed pressure drop to the one predicted by the well-known Kozeny-Carman equation (KC equation), using the concept of an equivalent sphere diameter to bring the observed predicted pressure drop into agreement with the one predicted by the KC equation.14,19,20 The KC equation for the superficial velocity usf is based on Darcy’s law for laminar flows (stating that usf ∼ ∆P/L) and on the assumption that the interstitial void space can be modeled as a bundle of parallel capillaries, for which it is well-known that the flow velocity is directly proportional to the second power of the diameter. In its basic form, the KC equation is usually written as22
usf )
1 ∆p Kva2 ηL
(1)
wherein Kv is a shape factor (for a packing of uniform spheres, (14) Gusev, I.; Huang, X.; Horvath, C. J. Chromatogr. 1999, 855, 273-290. (15) Bidlingmaier, B.; Unger, K. K.; von Doehren, N. J. Chromatogr. 1999, 832, 11-16. (16) Bristow, P. A.; Knox, J. H. Chromatographia 1977, 10, 279-289. (17) Tanaka, N.; Kobayashi, H.; Nakanishi, K.; Minakuchi, H.; Ishizuka, N. Anal. Chem. 2001, 420-429. (18) Tanaka, N.; Kobayashi, H.; Ishizuka, N.; Minakuchi, H.; Nakanishi, K.; Hosoya, K.; Ikegami, T. J. Chromatogr. 2002, 965, 35-49. (19) Hahn, R.; Jungbauer, A. Anal. Chem. 2000, 72, 4853-4858. (20) Leinweber, F. C.; Lubda, D.; Cabrera, K.; Tallarek, U. Anal. Chem. 2002, 74, 2470-2477. (21) Coulson, J. M.; Richardson, J. F. Chemical Engineering; Volume 2, Particle Technology and Separation Processes; Pergamon Press: Oxford, 1991. (22) Fourie, J. G.; Du Plessis, J. P. Chem. Eng. Sci. 2002, 57, 2781-2789. (23) Puncochar, M.; Drahos, J. Chem. Eng. Sci. 2000, 55, 3951-3954. (24) Astro¨m, T. B.; Pipes, B. R.; Advani, S. G. J. Compos. Mater. 1992, 26, 13511373. (25) Westerterp, K. R.; Wijngaarden, R. J.; Nijhuis, N. B. G. Chem. Eng. Technol. 1996, 19, 291-298. (26) Thompson, K. E.; Fogler, H. S. AIChE J. 1997, 43, 1337-1389.
Analytical Chemistry, Vol. 75, No. 4, February 15, 2003 843
Kv = 5; for a cylindrical tube, Kv = 2) and a is the specific external surface of the packing, which is, for a packing of spheres given by
a)
1- 6 dp
(2)
Combining eqs 1 and 2, the following well-known expression is obtained,
usf )
∆pdp 3 1 2 Kc (1 - ) ηL
(3)
with Kc)180. The KC equation, including the Kc ) 180 value, is well-validated for the packed bed of spheres, but for obvious reasons, this validation has always been limited to a small range of porosities around ) 0.4. For systems for which a wider range of porosities can be established, such as, for example, fiber beds or metallic foam beds, it is, however, well-known22-24 that the Kozeny constant changes dramatically with the porosity as soon as porosities above a value of ) 0.5-0.6 are considered. For example, for the flow through a bed of silk fibers, it has been found21 that the Kozeny constant increases from a value of 180 around ) 0.4 to a value of 300 around ) 0.9 and to a value of 500 around ) 0.95. This variability of the Kc factor implies that the assumptions underlying Kozeny-Carman’s law fail to represent the flow resistance of porous structures having a much more open structure than the typical packed bed case. Since monolithic silica columns are also intended for use over a broad range of bed porosities up to ) 0.85, it is straightforward to expect that the KC equation is unsuited to represent the variation of the pressure drop in silica monoliths over a broad range of bed porosities. For example, this is illustrated by the fact that the equivalent sphere diameter derived from the KC equation usually appears to be totally unrelated to the equivalent sphere diameter derived from the mass transfer data20 or from the mean size of the skeleton estimated from SEMs.19 One of the shortcomings of macroscopic models such as the KC-model is that they lump the internal bed structure into a number of empirical constants (tortuosity, sphericity factor, etc.), whereas it is a well-established fact from the field of hydrodynamics that the mean flow velocity through a given pore depends heavily on the exact pore topology and that many phenomena occurring in packed beds and porous media require an understanding of the flow at the very smallest scales.25-27 To the best of our knowledge, the single attempt to develop a model relating the observed pressure drop in monolithic columns to the microscopic details of the monolith structure (pore size distribution, pore connectivity number, etc.) has been made by Meyers & Liapis.28-30 They used a network-modeling approach wherein a number of so-called flow nodes is interconnected by nT cylindrically shaped pores with variable diameter. The pressure (27) Sederman, A. J.; Johns, M. L.; Alexander, P.; Gladen, L. F. Chem. Eng. Sci. 1998, 53, 2117-2128. (28) Meyers, J. J.; Liapis, A. J. Chromatogr. 1998, 827, 197-213. (29) Liapis, A. I.; Meyers, J. J.; Crosser, O. K. J. Chromatogr. 1999, 865, 13-25. (30) Meyers, J. J.; Liapis, A. I. J. Chromatogr. 1999, 852, 3-23.
844
Analytical Chemistry, Vol. 75, No. 4, February 15, 2003
drop is then calculated from the pressure drop in the cylindrical pores. The drawback of this model is that it considers only tangentially oriented flows (i.e., flowing in parallel with the bed structure) and does not account for the fact that a large part of the flow field in a silica monolith is directed perpendicularly to the skeleton (so-called frontal flow). Because it can be inferred that the latter will contribute significantly to the observed flow resistance, it should be obvious that a more faithful representation of the internal pore structure might lead to an improved pressuredrop correlation, representing the link between the bed porosity and the observed pressure drop in a more accurate way than the KC model. Such an improved correlation would allow interpretation of the experimental pressure drop data of silica monoliths in a much more specific way, that is, directly related to the main structural properties of the bed (skeleton thickness, bed porosity, flowthrough pore size, specific skeleton surface). The availability of an accurate correlation between the pressure drop and the bed porosity would also enable a theoretical determination of the optimal bed porosity of silica monoliths, for example, by setting off the need for large theoretical plate numbers (requiring a small flow resistance and, hence, a large bed porosity) against the need for a large retention capacity and column loadability (requiring small bed porosities). In the present study, we have considered a geometric model that mimics the internal structure of a silica monolith as closely as possible. To obtain the desired information on the velocity field and the corresponding pressure drop with the largest possible accuracy, we have resorted to so-called computational fluid dynamics (CFD) simulations. CFD simulations represent nothing else but a numerical solution of the Navier-Stokes equations determining the flow field, using an integrated software package combining a sophisticated grid generator and an equation solver. Since the introduction of commercial CFD software and with the dramatic increase in computational power of commercial PCs and workstations, CFD simulations have over the past decade become an important and widely accepted tool in the field of chemical engineering.31 They are ideally suited to calculate the velocity field in complex geometries and provide a unique opportunity to yield detailed velocity data for flows through spaces that are inaccessible for experimental measurement methods. CFD simulations have, for example, also become an increasingly used design tool in the lab-on-a-chip field.32-34 Since the flow under investigation is without any doubt purely laminar, the Navier-Stokes equation are, unlike for turbulent flows, 100% deterministic. This implies that the CFD results can be made fully reliable, provided care is taken that the obtained results are independent of the grid shape and size. CONSIDERED STRUCTURE AND CFD SOLUTION METHOD Conceiving a simplified structure mimicking the geometric features of the spongelike geometry7,18 of monolithic silica beds as closely as possible, one of the straightforward options is the (31) Bakker, A.; Haidari, A. H.; Marshall, E. M. Chem. Eng. Prog. December 2001, 30-39. (32) Patankar, N. A.; Hu, H. H. Anal. Chem. 1998, 70, 1870-1881. (33) Ermakov, S. V.; Jacobson, S. C.; Ramsey, J. M. Anal. Chem. 2000, 72, 35123517. (34) Molho, J. I.; Herr, A. E.; Mosier, B. P.; Santiago, J. G.; Kenny, T. W.; Brennen, R. A.; Gordon, G. B.; Mohammadi, B. Anal. Chem. 2001, 73, 1350-60.
Table 1. Overview of Considered Combinations of ds and ls and Corresponding E Values
Figure 1. (a) Detailed view of skeleton structure and definition of the geometric parameters ds and ls. (b) Global view of unit cell of tetrahedral skeleton column (TSC) for the case of ds ) 2 µm and ls ) 5 µm (corresponding to ) 0.863). The (x-y)- and (x-z)-side planes are treated as a “periodic wall” boundary; the (y-z)-side plane is treated as a “symmetry wall” boundary.
crystal lattice structure of diamond in which each carbon atom is bonded tetrahedrally to its four neighbors. An appealing advantage of this so-called tetrahedral skeleton column (TSC) structure is that all of its structural properties (porosity, flow-through pore size, specific skeleton surface, etc.) can be represented using only two single parameters: the skeleton diameter, ds, and the skeleton unit length, ls. Both parameters are defined in Figure 1a. Varying the skeleton unit length ls between ds and a given large multitude of ds’s, the TSC model allows, for each possible value of the skeleton thickness, establishment of a complete range of possible bed porosities going from the minimal porosity min ) 0.378 (for ds ) ls) all the way up to ) 1. It is, of course, obvious that the TSC model is only a simplification of the reality, because it lacks, for example, the stochastic nature of the distribution of the flow-through pore size of a true monolithic column, but it has been considered that the best way to gain insight into the complex relation between the internal structure of a monolithic column and its flow resistance is to start off with a system exhibiting an intermediate degree of complexity before moving on to more sophisticated models. The CFD simulations were carried out with a commercial CFD package, FLUENT (v.6). The software was installed on a PC with an Intel Pentium IV processor running at 1500 MHz and equipped with 1 Gb RAM. The grids were generated with GAMBIT (v.2) software run on the same hardware configuration. To drastically
ds (µm)
ls (µm)
τ
1 1 1 1 1 1 1 1 1 2 3 4 5
1 1.05 1.15 1.35 1.5 2 3 4 5 5 5 5 5
0.378 0.419 0.492 0.603 0.667 0.797 0.902 0.943 0.963 0.863 0.720 0.553 0.378
1.200 1.185 1.158 1.121 1.101 1.064 1.037 1.026 1.021 1.048 1.086 1.135 1.200
reduce the computational effort, all simulations were carried out on a unit cell (Figure 1b), representing the minimal building block from which an infinite tetrahedal skeleton column can be reconstructed. Since the tetrahedral skeleton has one symmetry plane, two walls of the unit cell could be set as “symmetry wall” (implying the use of a slip boundary condition for the calculation of the velocity field); the other boundaries were set as “periodic wall” (see Figure 1b). With these symmetry and periodic boundary conditions, the unit cell behaves as if it were embedded in an infinitely replicated structure. The periodic condition requires the specification of the total mass flow rate across the periodic inand outlet planes. For each combination of ls and ds, the simulation was repeated for at least five different mass flow rates, chosen such that the imposed average interstitial fluid velocities were in the range of 10-3 and 2 × 10-2 m/s. Water, with a density of 1000 kg/m3 and a viscosity of 1.003 10-3 kg/(m.s), was chosen as the working fluid. The corresponding Reynolds numbers (based upon the skeleton diameter as the characteristic length) were typically on the order of 0.001- 0.1, such that it can safely be concluded that the flow conditions were always strictly laminar. The pore volume space was discretized with an unstructured tetrahedral grid, because the setup of a structured hexahedral grid on the rather complicated void-space geometries would have taken too long. The size of the computational domain ranged from a value of 40 000 for the smallest cylinder length up to a value of 400 000 cells for the largest cylinder length. The total volume (Vs) of the skeleton structure could be directly obtained from the volume statistics report function of the Fluent software. From this Vs value, the effective bed porosity () of the TSC could be directly calculated using
)1-
Vs Vuc
(4)
wherein Vuc is the total volume () sum of solid and fluid portions) of the unit cell depicted in Figure 1b. Table 1 gives an overview of all the different combinations of the ds- and ls values considered in the present study, together with the Vs and values obtained from the Fluent software. As can be noted, the ds and ls values were varied between 1 and 5 µm, allowing us to cover nearly the entire range of published skeleton sizes and porosities.4-11 It can also be noted from Table 1 that the ds- and ls combinations have Analytical Chemistry, Vol. 75, No. 4, February 15, 2003
845
Figure 2. Example of calculated flow field in the top (y-z) plane of the unit cell shown in Figure 1b (color scale in m/s).
been selected such that we could compare the influence of the skeleton diameter at constant bed porosity, as well as the influence of the bed porosity at constant skeleton diameter. Figure 2 shows an example of the calculated velocity field in the top (y, z)-plane of the unit cell. The velocity map clearly shows that in the regions far away from the skeleton structure, the velocity is much larger than in the regions immediately adjacent to the skeleton, where the fluid is nearly stagnant. Fully similar plots were obtained for the other considered flow conditions and geometries. All of the calculated velocity fields were subsequently characterized by calculating the volume-averaged velocity magnitude 〈u〉, the volume-averaged x component of the velocity 〈ux〉, the total volumetric flow rate Q passing through each cross-sectional plane of the unit cell, and the total pressure gradient ∆P/L. It should be noted that the thus-obtained 〈ux〉 velocity values correspond to the linear flow velocity traditionally used in pressure drop correlations and that the ratio of the volumetric flow rate Q and the total cross-sectional area of the unit cell Auc corresponds to the superficial velocity usf.
usf )
Q Auc
(5)
The volume-averaged velocity magnitude 〈u〉 is a measure for the true interstitial velocity. It is one of the interesting features of the CFD approach that this value can actually be determined. From the theory of hydrodynamics, it is a well-known fact33 that the 〈u〉 value is related to the 〈ux〉 value via the bed tortuosity τ,
τ)
〈u〉 〈ux〉
(6)
but the superficial velocity usf and the volume-averaged x component of the velocity 〈ux〉 are in turn related by the bed porosity.
)
usf 〈ux〉
(7)
RESULTS AND DISCUSSION Geometrical Properties of the Tetrahedron Model Structure. Using a number of geometrical considerations and an Excell846
Analytical Chemistry, Vol. 75, No. 4, February 15, 2003
Figure 3. Variation of with ds and ls for the tetrahedral skeleton column structure. The circular symbols relate to the data given in Table 1. The full lines are calculated on the basis of eq 8.
based least-squares-of-error parameter fitting procedure, the following relation between the through-pore porosity, , and the basic parameters ds and ls of the tetrahedral skeleton column could be established.
)1-
ds2(ls/2 - ds/5.27) 4(ls/2)3
(8)
Figure 3 shows the excellent agreement between the actual values, as obtained from the Fluent software, and the values predicted by the correlation given in eq 8. Validation and Discussion of the Calculated Velocity Fields. The first prerequisite of a reliable set of CFD data is that they are grid-independent. This has been ensured by performing, after the initial 100-300 iteration steps, a velocity gradient adaptation of the grid before running another 200 iterations. The change in the velocity field values and pressure drop was always within 1-2%, a small enough difference to conclude that the solutions are grid-independent and sufficiently converged. A second point of validation can be obtained from the observation that the relation between the calculated pressure drop gradient (∆P/L) and the volume-averaged linear flow velocity 〈ux〉 is perfectly linear (Figure 4), as predicted by Darcy’s law for laminar flows through porous media.36 Whereas Figure 4 shows only the data for the ds ) 1 µm and ls ) 2 µm and the ds ) 1 µm and ls ) 3 µm cases, it should be noted that the linearity of the (〈ux〉,∆P/L) data for the other considered ds and ls combinations was equally good. A third validation of the CFD calculations follows from the fact that, when calculating the ratio of the obtained 〈ux〉 values to the superficial velocities (calculated via eq 5) and plotting this ratio as a function of the actual porosity () listed in Table 1, a nearly perfect one-to-one relation is obtained (r2-value of >0.999), in full agreement with eq 7. (35) Seguin, D.; Montillet, A.; Comiti, J. Chem. Eng. Sci. 1998, 53, 3751-376. (36) Kaviany, M. Principles of Heat Transfer in Porous Media, 1st ed.; SpringerVerlag: New York, 1995.
Figure 4. Plot of ∆P/L versus 〈ux〉 for two different combinations of ds and ls, showing the excellent linearity of the (∆P/L, 〈ux〉)-relationship of the CFD data.
Figure 5. 〈u〉/〈ux〉 versus 〈ux〉 for five different combinations of ds and ls. The ratio of 〈u〉 to 〈ux〉 corresponds to the bed tortuosity τ (cf. eq 6).
A final validation originates from the fact that, for a given combination of ds and ls, the ratio of the 〈u〉 and 〈ux〉 values remains constant over the entire range of imposed flow velocities (Figure 5). Because this ratio corresponds to the bed tortuosity (cf. eq 6), which should be a velocity-independent geometric feature of the bed structure, it is, indeed, obvious to expect that the 〈u〉/ 〈ux〉 ratio is independent of the flow velocity. Figure 5 shows the 〈u〉/〈ux〉 ratio for 5 different ds and ls combinations. Comparison of the Obtained Flow Resistance Values with Experimental Data. To generalize our results, the obtained ∆P/L values have been interpreted in terms of the flow resistance parameter φ, defined as
φ)
ds2 K
(9)
with the column permeability, K, in turn being defined as
K)
〈ux〉ηL ∆P
(10)
Figure 6 shows a plot of the column permeability, K, as a function of the skeleton diameter, ds, for a given porosity. Because
Figure 6. Variation of the column permeability K (KTSC) of the tetrahedral skeleton column with the skeleton diameter ds (constant porosity ) 0.378).
Figure 7. Variation of the flow resistance (φTSC) of the tetrahedral skeleton column with (red solid circle) for all considered combinations of ds and ls and comparison with data taken from literature (see Table 2 for origin of experimental data points: 2, entry 1; b, entries 2-4; ), entries 5-6; 4, entries 7-10; 0, entries 11-14; 9, entry 15; O, entry 16).
of the log-log representation, it can clearly be concluded that, for systems with the same porosity, the permeability of the tetrahedal skeleton column varies exactly according to the second power of the skeleton diameter ds (K ∼ ds2). This finding, in fact, constitutes an additional validation of our calculations, because the K ∼ ds2 relation is known to be a characteristic feature of all laminar flows through porous media.36 Moving then to a representation in terms of flow resistance values (φTSC) and considering the definition of φ given in eq 9, it is not surprising to find that all of the data points of Figure 6 reduce to a single point in the (, φ) representation (Figure 7). As can be noted, the φ values for all of the other considered ls and ds combinations, all fall on a single curve. Exploring several mathematical expressions, we found that the dependency of the φTSC data on the bed porosity can be very accurately described by the following simple correlation.
φ ) 55
(1 - )
1.55
Analytical Chemistry, Vol. 75, No. 4, February 15, 2003
(11) 847
Q h )
Table 2. List of Flow Resistance Data Taken from Literature
1.(15) 2.(10) 3.(10) 4.(10) 5.(16) 6.(16) 7.(6) 8.(6) 9.(6) 10.(6) 11.(7) 12.(7) 13.(7) 14.(7) 15.(5) 16.(13)
ds (µm)
K (10-14 m2)
φa
2.20 2.0 2.0 2.0 1.60 2.0 2.31 1.58 1.16 1.06 1.59 1.34 1.16 1.00 0.7 0.8e
0.85 0.83 0.83 0.83 0.62 0.86 0.62 0.62 0.62 0.61 0.46b 0.52b 0.58b 0.65 0.65 0.75e
100 120 100 130 40 100 40 17 13 5.6 4.22c 4.37c 5.93c 6.23c 2.49d 8.50e
4.8 3.3 4.0 3.1 6.4 4.0 13.3 14.7 10.3 20.1 59.9 41.1 22.7 16.1 19.7 7.5
a Calculated from the K values on the basis of eqs 9-10. b Estimated from the value of entry 14 (also taken from ref 7) and the variation of with ds and ls of the TSC model. c Value estimated from Figure 8 of ref 7, assuming a value of η ) 6 × 10-4 kg/(m‚s) for the viscosity. d Value estimated from Figure 3 of ref 5, assuming a value of η ) 6 × 10-4 kg/(m‚s) for the viscosity. e Values for rod 178 of ref 13, ds estimated from the given porosity and the given flow-through pore diameter, and K taken as the mean of the values given in Figure 1 of ref 13. (f) Unless otherwise stated, the ds data in column 2 are taken directly from the corresponding literature reference.
Figure 7 shows that the only small deviation between the actual φTSC values and the ones predicted by eq 11 occurs in the > 0.95 range, but this range is nevertheless considered to be irrelevant for practical purposes. In Figure 7, we have also added a large number of experimental data found in the literature. The origin of these data is given in Table 2. It should be noted that some of the parameters needed to establish these data points were lacking or had to be estimated from the figures in the papers. This undoubtedly explains the relatively large scatter of the data. Despite this scatter, and although there is no perfect agreement, the tetrahedral skeleton column model clearly yields a good approximation for the experimentally observed (φ, ) relation. Being confronted with the relatively large scatter, and especially with the fact that most of the data points fall below the φTSC line, we postulated that part of the scatter should be due to the flow-through pore size heterogeneity of the columns. When a bed is composed of pores with different sizes, and considering that the flow rate through the larger pores will be more than proportionally larger than the flow through the smaller pores (flow rate ∼ (pore diameter)4) it is obvious to expect that the observed flow resistance in the case of a distribution of pore sizes will always be smaller than the flow resistance calculated on the basis of the average pore size. To quantify this effect, we have adapted a procedure outlined by Schisla et al.37 for the calculation of the effective pressure drop in a bundle of capillaries with a polydisperse diameter. In the present calculation, we represent the variation of the flow-trough pore diameter as a variation of the local value in a small region of the cross-sectional area of the bed. The observed flow rate is then defined as (37) Schisla, D. K.; Ding, H.; Carr, P. W.; Cussler, E. L. AIChE J. 1993, 39, 946.
848
Analytical Chemistry, Vol. 75, No. 4, February 15, 2003
∫u
x,locloc
A
dA,
(12)
wherein A represents the entire cross-sectional area of the bed and ux,loc is the local x velocity component measured at a region of the cross-sectional plane of the bed where the local porosity has a value loc. From the value of Q h , it is straightforward to define the mean x component of the interstitial velocity 〈ux〉 as
Q h
〈ux〉 )
∫ A
loc
(13) dA
To represent the variation of the local porosity in terms of the parameters of the TSC model, we have assumed that the skeleton thickness, ds, remains invariable and that the changes in local porosity are exclusively due to changes in the segment length ls. Although this is only an approximate approach, it has the advantage that the variation of the local bed porosity can be described by a single parameter. We believe, furthermore, that the outcome of this approach is not fundamentally different from what would have been obtained when the skeleton thickness would also have been treated as a variable. Defining now a dimensionless parameter, z,
ls - lsm lsm
z)
(14)
to represent the variation of the unit skelet length around its mean value lsm, eq 8 can be rewritten as
(z) ) 1 -
ds2[lsm(1 + z)/2 - ds/5.27] , 4[lsm(1 + z)/2]3
(15)
Assuming, then, a Gaussian distribution (with standard deviation σz) for the value of the parameter z, eqs (12-13) can be reformulated to yield
〈ux〉 )
∫
+∞
-∞
∫
ux(z)A(z)e-z /2.σz2 dz 2
+∞
-∞
A(z)e-z /2.σz2 dz 2
(16)
with (z) given by eq 15, and with ux(z) given by (cf. eqs 9-10)
ux(z) )
2 1 ds ∆P φ(z) ηL
(17)
wherein φ(z) is calculated using eq 11 and 15. From eq 16 and from the definition of φ (cf. eq 9), we can now calculate φcorr, representing the flow resistance of the TSC corrected for the existence of a certain degree of flow-through pore heterogeneity, as
φcorr )
2 1 ds ∆P 〈u 〉 ηL x
(18)
Figure 8. Reconsideration of the flow resistance data of Figure 8, but now compared to the theoretical φ values corrected for the flowthrough-pore heterogeneity. The dashed lines have been obtained by plotting the φcorr values given by eq 18 versus the corrected values given by eq 19 using, respectively, σz ) 0.1, σz ) 0.2, σz ) 0.3, and σz ) 0.4 as the value for the standard deviation of the skeleton unit length variation.
Plotting this value as a function of the corrected value (corr),
∫ (z)e ∫ e +∞
corr )
-∞ +∞
-z2/2.σz2
-z2/2.σz2
-∞
dz (19)
dz
for different values of σz, the dashed lines shown in Figure 8 are obtained. As can be noted, the experimental data now perfectly fit into the bundle of φcorr lines. The slight “wavy” nature of the φcorr lines is due to the complex interplay between z and (cf. eq 15). The trend exhibited by the φcorr lines is also in agreement with our intuitive supposition that the influence of a given degree of heterogeneity is maximal for the cases with the smallest average porosity and should vanish nearly completely when the average porosity tends to unity. Comparison with Kozeny-Carman’s Law. From eq 3, together with eqs 9-10, it can easily be found that the KozenyCarman model leads to the following expression for the flow resistance.
(1 - )
φKC ) Kc
2
(20)
In eq 20, the Kc constant (Kc ) 180 for the packed bed of spheres) has been deliberately left as a freely selectable parameter, because Kc is precisely the constant that is modified when the KC equation is used12,17,18 to define an equivalent sphere diameter on the basis of the experimental flow resistance values. To compare the TSC model with this traditional KC-modeling approach, Figure 9 shows a plot of the ratio of φKC to φTSC as a function of for four different values of Kc. One value is the traditional Kc ) 180 value, whereas the three other Kc values relate to the value that is obtained when using the KC model to fit the flow resistance of the tetrahedal skeleton column for a given value of . As can be noted, the KC model can be brought into agreement with the TSC-model for
Figure 9. Comparison of flow resistances obtained for the tetrahedral skeleton column and those predicted with the KC equation. Ratio of φKC to φTSC as a function of for four different values of Kc.
only one given value of . Doing so, the φKC values for the other values are inevitably strongly offset. CONCLUSIONS The TSC model allows linking of the observed flow resistance of a given silica monolith to two easily determinable structural parameters of the bed (skeleton thickness and bed porosity) and, hence, omits the need to resort to an effective particle diameter, as is the case for the traditionally employed Kozeny-Carman equation. The TSC model also allows for a much more accurate representation of the variation of the flow resistance with the bed porosity, hence, making it much better suited for porosity optimization calculations. The TSC model can also easily be corrected for the throughflow pore heterogeneity. In doing so, an excellent agreement with the experimental data is obtained. The proposed correction method even seems to open the road to the development of a procedure that permits rapid estimation of the degree of flowthrough pore heterogeneity from the difference between the experimentally observed flow resistance and the resistance predicted by eq 11. This is, however, still a speculative idea and requires further validation. Future work will therefore include a more in-depth investigation of the influence of the flow-through pore size heterogeneity. Other types of model structures (with both uniform and nonuniform flow-through pores) will be considered as well. As a final remark, it can be concluded that the currently proposed approach to use CFD simulations to make detailed calculations of the internal velocity field in an idealized model structure closely mimicking the internal structure of a chromatographic column is promising enough to continue to explore the use of this method for the determination of other important chromatographic parameters (eddy-diffusion, local film transfer coefficients, etc.). Symbols a
external specific surface (m2/m3)
Auc
cross-sectional area of the unit cell (m2)
ls
TSC skeleton length (m) Analytical Chemistry, Vol. 75, No. 4, February 15, 2003
849
dp
particle diameter (m)
τ
ds
TSC skeleton diameter (m)
Subscripts
K
permeability (m2)
Kc
Kozeny-Carman constant
Kv
shape factor
L
column length (m)
lsm
mean skeleton length (m)
Abbreviations
Q
volumetric flow rate (m3/s)
CFD
computational fluid dynamics
usf
superficial velocity (m/s)
KC
Kozeny-Carman
〈u〉
volume-averaged velocity (m/s)
TSC
tetrahedal skeleton column
〈ux〉
volume-averaged x component of velocity (m/s)
Vs
skeleton volume (m3)
Vuc
total volume (solid + fluid fraction) of the considered unit cell (m3)
Greek Symbols
tortuosity
corr
corrected
loc
local
uc
unit cell
ACKNOWLEDGMENT The authors greatly acknowledge a Research Grant (FWO KNO 81/00) of the Fund for Scientific Research-Flanders (Belgium). P.G. is supported through a specialization grant from the Instituut voor Wetenschap en Technologie (IWT) of the Flanders Region (Grant no. SB/11419).
∆p
pressure drop (Pa/m)
porosity
φ
flow resistance
η
viscosity (kg/ms)
Received for review October 8, 2002. Accepted December 4, 2002.
σz
standard deviation
AC0262199
850
Analytical Chemistry, Vol. 75, No. 4, February 15, 2003
