Anal. Chem. 2003, 75, 901-909
Electroosmotic Dispersion in Microchannels with a Thin Double Layer Emilij K. Zholkovskij
Institute of Bio-Colloid Chemistry of Ukrainian Academy of Sciences, Vernadskogo, 42, 03142 Kiev, Ukraine Jacob H. Masliyah*
Department of Chemical and Materials Engineering, University of Alberta, 536 Chemical-Mineral Engineering Building, Edmonton, Alberta, Canada T6N 2G6 Jan Czarnecki
Edmonton Research Centre, Syncrude Canada Ltd., 9421-17 Avenue, Edmonton, Alberta, Canada T6N 1H4
Dispersion of a nonelectrolyte solute due to the electroosmotic flow in long straight microchannels was analyzed theoretically. A version of the Aris-Taylor procedure was employed to predict the dispersion coefficient for arbitrary geometry of the microchannel cross section. The analysis was conducted using a thin double-layer approximation, which is valid when the Debye length is much smaller than the characteristic dimensions of the cross section. For thin double layers, the obtained results describe the electroosmotic dispersion for arbitrary surface potential, electrolyte type, and cross-section geometry. Dispersion for several cases of the cross-section geometries was discussed. It was shown that, for given values of the surface potential and the Debye length, both the cross-section geometry and the electrolyte content of the driven solution substantially affect the dispersion of a nonelectrolyte solute. In the relevant particular cases, the obtained results agree with predictions of the previous theories. Microchannel networks are increasingly used in chemistry, biochemistry, and medicine for the separation, analyses, and synthesis of substances.1-5 In such systems, electroosmosis is often employed as a method for sample transport and separation. Longitudinal spreading (dispersion) of a solute normally takes place in a convective transport through a long straight channel. At zero solute adsorption, the spreading occurs due to the solute molecular diffusion and due to the nonuniformity of the velocity distribution within the channel cross section. Such a solute dispersion is fairly profound in the case of the pressure-driven * Corresponding author. E-mail:
[email protected]. (1) Harrison, D. J.; Fluri, K.; Seiler, K.; Fan, Z.; Effenhauser, C. S.; Manz, A. Science 1993, 261, 895-897. (2) Jacobson, S. C.; Hergenroder, R.; Kounty, L. B.; Ramsey, J. M. Anal. Chem. 1994, 66, 1114-1118. (3) Manz, A.; Effenhauser, C. S.; Burgraf, N.; Harrison, D. J.; Seiler, K.; Fluri, K. J. Micromech. Microeng. 1994, 4, 257-265. (4) Jacobson, S. C.; Ramsey, J. M. In High Performance Capillary Electrophoresis; Khaledi, M. G., Ed.; John Wiley & Sons: New York, 1998; Vol. 146, pp 613-633. (5) Gulbertson, C. T.; Jacobson, S. C.; Ramsey, J. M. Anal. Chem. 2000, 72, 5814-5819. 10.1021/ac0203591 CCC: $25.00 Published on Web 01/22/2003
© 2003 American Chemical Society
(Poiseuille) flow when the flow velocity profile is of a parabolic type. On the other hand, the electroosmotic flows are characterized by a velocity profile that is nearly uniform within the whole cross section except the thin double-layer regions adjacent to the wall. Therefore, at equal mean velocities, the electroosmotic dispersion is much weaker than that due to the Poiseuille flow. Nevertheless, dispersion produced by purely electroosmotic flow can play an important role. At sufficiently low electrolyte concentrations (lower than 10-3 mol/L) and strong applied electric fields (on the order 1 kV/cm), electroosmotic dispersion can dominate as compared with that due to molecular diffusion.6-15 Electroosmotic dispersion was theoretically studied by different authors.6-15 Analysis of the dispersion phenomenon is often based on the approach developed by Taylor16 (for the Poiseuille flow in a circular tube) and generalized by Aris17,18 (for arbitrary crosssection shape). According to the Taylor-Aris theories, for sufficiently long times, the mean cross-sectional concentration of a solute satisfies an equation having the same structure as the onedimensional convective-diffusion equation. This equation contains a longitudinal “diffusion” term that is proportional to the apparent diffusion coefficient, D/, which is referred to as the dispersion coefficient and is given by the relation
D/ ) D[1 + (〈u〉d/D)2]
(1)
where D is the molecular diffusion coefficient of the nonelectrolyte solute, u is the local liquid velocity, and 〈...〉 ) (1/S)∫S ... dS (6) Martin, M.; Guiochon, G. Anal. Chem. 1984, 56, 614-620. (7) Martin, M.; Guiochon, G.; Walbroehl, Y.; Jorgenson, J. W. Anal. Chem. 1985, 57, 559-561. (8) Datta, R.; Kotamarthi, V. R. AIChE J. 1990, 36, 916-926. (9) Datta, R. Biotechnol. Prog. 1990, 6, 485-493. (10) McEldoon, J. P.; Datta, R. Anal. Chem. 1992, 64, 227-230. (11) Andreev, V. P.; Lizin, E. E. Electrophoresis 1992, 13, 832-837. (12) Andreev, V. P.; Lizin, E. E. Chromatographia 1993, 37, 202-210. (13) Gas, B.; Stedry, M.; Kenndler, E. J. Chromatogr., A 1995, 709, 63-68. (14) Griffiths S. K.; Nilson, R. H. Anal. Chem. 1999, 71, 5522-5529. (15) Griffiths S. K.; Nilson, R. H. Anal. Chem. 2000, 72, 4767-4777. (16) Taylor, G. Proc. R. Soc. London 1953, 219A, 186-203. (17) Aris, G. Proc. R. Soc. London 1956, 235A, 67-77. (18) Aris, G. Proc. R. Soc. London 1959, 252A, 538-550.
Analytical Chemistry, Vol. 75, No. 4, February 15, 2003 901
signifies an average over the cross-sectional area, S. Here, d is a length-scale parameter that can be obtained using the TaylorAris procedure. At zero solute adsorption, d is completely defined by the velocity distribution within the cross section. In previous studies,6-15 dispersion is often described in equivalent terms of a plate height, H, which is interrelated with D/ in a unique manner (H ) 2D//〈u〉). Within the framework of the Taylor-Aris approach, Martin and Guiochon,6 Martin et al.,7 Datta and Kotamarthi,8 Datta,9 and McEldoon and Datta,10 evaluated D/, H, or both by assuming low electric potentials at the wall/solution interface. Consequently, in these studies, the electroosmotic velocity distribution was obtained using the linearized version of the Poisson-Boltzmann equation (Debye equation). The Debye equation was also employed by Griffiths and Nilson,14 who suggested an original perturbation technique for describing the hydrodynamic and concentration profiles inside the channel. To evaluate the electroosmotic dispersion for high interfacial potentials, some authors (Andreev and Lizin,11,12 Gas et al.,13 Griffiths and Nilson15) used the nonlinear Poisson-Boltzmann equation. In the paper of Gas et al., dispersion for a circular tube was studied on the basis of a simplified version of the method proposed by Martin and Guiochon6 and Martin et al.7 The result of Gas et al. is equivalent to eq 1 provided that the parameter d ) β/21/2, where β is the effective double-layer thickness. Generally, the latter parameter differs from the Debye length, κ-1, and depends on the interfacial potential. The parameter β is evaluated by Gas et al. by making use of a numerical solution of the nonlinearized Poisson - Boltzmann problem. The result of Gas et al. holds when β , a, where a is the circular cross-section radius, i.e., for thin double layers. The predictions of Gas et al. agree with the results of Griffiths and Nilson,15 who obtained the velocity and concentration distributions for arbitrary surface potentials. For a thin double layer, Griffiths and Nilson15 suggested an asymptotic expression which corresponds to (1) when the parameter d is given by
d 2 ) g(ζ)K/κ2
Analytical Chemistry, Vol. 75, No. 4, February 15, 2003
ELECTROOSMOTIC DISPERSION FOR ARBITRARY CROSS-SECTION GEOMETRY Let us formulate the Taylor-Aris procedure in the form that is convenient for describing an arbitrary cross-section geometry. The parameter d from eq 1 can be represented as
d 2 ) (1/〈u〉2)〈[u - 〈u〉][ω - 〈ω〉]〉
(3)
where ω ) ω(rb) is a function satisfying the equation
(2)
In eq 2, the dependency on the surface potential, ζ, is embedded in the empirical interpolation function g(ζ), whose analytical form is given by the authors.15 The dimensionless factor K, is defined by the cross-section shape. For circular and plane-parallel cross sections, Griffiths and Nilson15 obtained K ) 1/2 and K ) 1/3, respectively. For the circular cross section, they reported that, within an accuracy of ∼5%, the parameter d given by (2) coincides with d ) β/21/2, which corresponds to the theory of Gas et al.13 The above-discussed studies6-15 have only addressed the dispersion in the circular and plane parallel cross sections. Results that can describe microchannels with high interfacial potentials have only been obtained for the symmetric electrolytes11-13,15 (although, the initial problem was formulated by Gas et al.13 and Andreev and Lizin11,12 for general-type electrolytes). The abovediscussed studies,7-10,14 which are valid for an arbitrary electrolyte, can only be used for low interfacial potentials. It should be noted that integrated devices employ microchannels with cross sections that substantially differ from circular and plane-parallel types.1-5,19,20 Note also, that microchannels with high 902
interfacial potentials (more than 50 mV) are often employed in the analytical techniques.13 It is important that, at high interfacial potentials, the electric double-layer properties are sensitive to changes in the electrolyte content of the equilibrium solution. Thus, a theory is required that is capable of describing the electroosmotic dispersion in the channels: (i) having different cross sections, (ii) bearing arbitrary electric potential at the wall/ solution interface, and (iii) containing a solution of the general electrolyte type. In the general case, developing such a theory is associated with a complicated mathematical problem that, for each given geometry, should be solved numerically. However, the analysis can be drastically simplified while considering the double-layer thickness to be much smaller than the cross-sectional dimensions. Such a consideration is quite reasonable for flow in microchannels. For example, for distilled water, κ-1 = 0.1 µm, whereas the shortest characteristic dimension of the cross section is usually more than 5 µm.1-5,19,20 Because of its small thickness, the double layer can be locally considered as quasi-flat. The quasi-flat model of the double layer bears substantial simplifications in the procedure of obtaining the dispersion coefficient. In the present paper, using this model, we will describe electroosmotic dispersion in microchannels with different geometries. First, we will obtain some general relations that are valid for arbitrary cross-section geometry, interfacial potential, and electrolyte type. Then, we will specify the results considering some particular examples.
∇ 2ω ) 〈u〉 - u
(4)
At zero solute adsorption, eq 4 is subject to the boundary conditions
∇ B ω‚n b ) 0 at the cross-section perimeter line
(5)
The radius vector, br ) bıxx + bıyy, belongs to the cross-section plane; the gradient operator is given by ∇ B ) bıx∂/∂x + bıy∂/∂y, where x and y are the 2D Cartesian coordinates in the cross-section plane; b n is the unity outward normal vector as shown in Figure 1. The electroosmotic velocity distribution, u( b), r within the cross section of a straight long channel was studied by many authors14,15,21-30 who considered circular, plane-parallel, and rect(19) Bianchi, F.; Wagner, F.; Hoffmann, P.; Girault, H. H. Anal. Chem. 2001, 73, 829-836. (20) Dutta, D.; Leighton T. D., Jr. Anal. Chem. 2001, 73, 504-513. (21) Burgreen, D.; Nakache, F. R. J. Phys. Chem. 1964, 68, 1084-1091. (22) Rice, D.; Whitehead, R. J. Phys. Chem. 1965, 69, 4017-4024. (23) Sorensen, T. S.; Koefoed, J. J. Chem. Soc., Faraday Trans. 2 1974, 70, 665-675.
Introducing the function
Ω ) Ψ - (η/E)u
(11)
and combining (8), (9), and (11), we obtain
∇2Ω ) 0
(12)
Using (7), (10), and (11), we derive the boundary condition for eq 12 Figure 1. Microchannel in electric field: (a) side view; (b) cross section.
angular cross sections. For each of these particular geometries, the authors employed a relation first proposed by Burgreen and Nakache,21 who studied the plane-parallel cross section:
u(b) r ) -(E/η)[ζ - Ψ(b)] r
(6)
where η and are the solution viscosity and dielectric permittivity, respectively, and E is the magnitude of the external electric field strength, which is a uniform vector parallel to the longitudinal axis Z (Figure 1a). The function Ψ(b) r describes the potential distribution within the cross section in the absence of the external electric field. At the wall, the potential is assumed to take a constant value, which is referred to as the electrokinetic potential, ζ, i.e.,
Ψ ) ζ at the cross-section perimeter line
(7)
The significance of eq 6 is that, when electric potential Ψ is known at a given location, the velocity becomes known at the same location. Equation 6 was shown to be valid for the abovementioned particular geometries.14,15,21-30 Let us demonstrate that eq 6 is valid for any cross-section geometry. For a long straight channel, Ψ satisfies the two-dimensional version of the Poisson equation
∇ 2 Ψ ) -(F/)
(8)
The internal problem (11) and (12) has the unique finite solution Ω(b) r ) ζ implying that Ω(b) r is uniform over the cross section. Combining this solution with (11), we immediately obtain (6). Thus, we showed that eq 6 describes any cross section regardless of its geometry. Returning to eqs 3-5 together with (6) and taking into account that 〈ω - 〈ω〉〉 ≡ 0, we obtain
d2 )
〈(ω′ - 〈ω〉′)Ψ〉
(14)
(ζ - 〈Ψ〉)2
∇ 2ω′ ) 〈Ψ〉 - Ψ
(15a)
b ) 0 at the perimeter line ∇ B ω′‚n
(15b)
where ω′ ) ωη/E. Hence, for arbitrary cross-section geometry, to describe the electroosmotic dispersion as defined by eq 1, the parameter d of eq 1 is evaluated using (14), where the function ω′(b) r can be obtained as a solution of the problem given by (15). The electric potential distribution, Ψ(b), r in (14) and (15) is determined from the Poisson-Boltzmann equation. The latter is derived from (8) using the Boltzmann distribution for ions that form the space charge
∇ 2Ψ ) -
F
(
∑C Z exp i i
i
ΨFZi RT
)
(16)
(9)
(10)
APPROXIMATION OF THIN DOUBLE LAYER Now, we will consider how the procedure given by (7) and (14)-(16) is simplified due to the thin double-layer approximation. Let us introduce parameter h as
subject to the boundary condition
u ) 0 at the cross-section perimeter line
(13)
Here, F is the Faraday constant and Zi and Ci are, respectively, the ith ion charge (in Faraday units) and molar concentration, which is attributed to an equilibrium solution. Accordingly, the potential Ψ is defined with reference to the equilibrium solution. Equation 16 is subject to boundary condition 7.
where F is the volumetric charge density. In the absence of a pressure gradient, the velocity satisfies the relevant version of the Stokes equation
η∇ 2u ) -FE
Ω ) ζ at the cross-section perimeter line
(24) Yang, C.; Li., D. J. Colloid Interface Sci. 1997, 194, 95-107. (25) Andreev, V. P.; Dubrovsky, S. G.; Stepanov, Y. V. J. Microcolumn Sep. 1997, 9, 443-450. (26) Arulanandam, S.; Li, D. Colloids Surf., A 2000, 161, 89-106. (27) Arulanandam, S.; Li, D. J. Colloid Interface Sci. 2000, 225, 421-428. (28) Dutta, P.; Beskok, A. Anal. Chem. 2001, 73, 1979-1986. (29) Santiago, J. G. Anal. Chem. 2001, 73, 2353-2365. (30) Hsu, J.-P.; Kao, C.-Y. J. Phys. Chem. 2001, 105, 8135-8142.
h ) S/l
(17)
where l is the cross-section perimeter. For a wide class of crosssection geometries, under the necessary condition
κh . 1 Analytical Chemistry, Vol. 75, No. 4, February 15, 2003
(18) 903
the solution Ψ of the Poisson-Boltzmann problem given by (16) and (7) can be approximated using the function Φ(ξ), where
Ψ ) (RT/F)Φ(ξ)
(19)
˜ | exp(-ξ). (RT/Fκ)(l/S)∫κ∆ 0 Φ dξ. It can be shown that |Φ| e |ζ Therefore, at κ∆ . 1, the above given integral very slightly depends on the upper integral limit, which can be replaced by infinity. Using the parameter h defined by (17), the cross-sectional mean potential is finally represented as
which is a solution of the Poisson-Boltzmann problem for the semi-infinite region
∑ν Z
d2Φ/dξ2 ) -
exp(-ΦZk)
k k
Ψ(0) ) ζ˜
(21)
Φ(∞) ) 0
(22)
Here, ζ˜ ) Fζ/RT; νk ) Ck/∑kCkZk2; the normalized coordinate ξ is defined with reference to a point at the cross-section perimeter line (Figure 1b).
ξ ) -κ(b r l - b)‚n r b
(23)
where brl is the radius vector of the ξ axis origin, which belongs to the perimeter line (Figure 1b), the Debye parameter κ is given by
∑C Z
2
i i
/RT
x
(26)
F 1 1 〈Ψ(ω′ - 〈ω′〉)〉 ≈ I ω′ dl - 〈ω′〉 RT κh l l
(
)∫ Φ(ξ) dξ (27) ∞
0
where the contour integral in the brackets is taken over the crosssection perimeter line. While deriving (27), we assumed that, in the close vicinity of the perimeter line, along the normal ξ-axis, ω′ is almost constant. This assumption will be verified after obtaining the expression for ω′. Now, we will obtain an approximate solution of the problem given by (15). It is convenient to represent the unknown function, ω′, as
ω′ ) h2〈Ψ〉W -(RT/κ2F)w
(28)
where the function w satisfies the equation
κ2∇ 2w ) (F/RT)Ψ
∑ν [exp(-ΦZ ) - 1] k
k
(25)
k
For a symmetric electrolyte (Z1 ) -Z2 and ν1 ) ν2 ) 1/2), the analytical result of the second integration, which is known as the Guy-Chapman expression, can be found elsewhere, for example, in ref 31. The error produced due to the above approximation is associated with the nonplanar local geometry of the wall and with the influence of the opposite part of the channel wall (overlap of the double layers). The nonzero local curvature (1/a) results in a relative error of order 1/κa, as stated by Van Aken et al.32 Very rough estimations show that the double-layer overlap results in an error that does not exceed exp(-κh). Let us briefly describe the procedure for obtaining an approximate expression for the cross-sectional mean potential, 〈Ψ〉 ) (1/S)∫SΨ ds, of (14) and (15a). One can single out a thin stripe adjacent to the perimeter line (Figure 1b). The stripe thickness, ∆, should satisfy the inequality κ-1 , ∆ , h. One can also assume that the potential Ψ takes a zero value out of the stripe and is described by (19)-(22) within the stripe. Accordingly, the integral over the cross-section area is reduced to the integral over the stripe. Taking this into account, and using (19) and (23), one can approximate the mean potential as 〈Ψ〉 ) (31) Masiyah, J. H. Electrokinetic transport phenomena; AOSTRA: Edmonton, 1994. (32) Van Aken, G. A.; Lekkeberger, H. N. W.; Overbeek, J. Th. G.; de Bruin, P. L. J. Phys. Chem. 1990, 94, 8468-8472.
904
∞
0
A similar procedure enables the following approximation of the mean value in the numerator of the right-hand side of eq 14
Problem 20-22 can be integrated once to give
2
∫ Φ(ξ) dξ
(24)
i
dΦ/dξ ) - sign(ζ)
1 F ≈ RT κh
(20)
k
κ2 ) F 2
〈Ψ〉
Analytical Chemistry, Vol. 75, No. 4, February 15, 2003
(29)
Using the assumptions employed while deriving (26) and taking into account (19), (23), and (29), the function w can be approximated as an arbitrary finite solution of the equation d2w/dξ2 ) Ψ(ξ). Therefore, with the same accuracy as that provided by approximate equalities 26 and 27, one can write
w ) w(ξ) )
∫ ∫ Φ(ξ′′) dξ′′ dξ′ ξ
∞
ξ′
∞
(30)
The integral in (30) converges because |Φ| e |ζ˜ | exp(-ξ). As for W, combining (15a), (26), (28), and (29) yields the following equation
∇ ˜ 2W ) 1
(31)
Also, combining (15b), (23), (26), (28), and (30), we obtain the boundary condition for eq 31, as
B b ) 1 at the perimeter line ∇ ˜ W‚n
(32)
In (31) and (32), we introduced the following normalization
1 b B B; ∇ ˜r ) b r; ∇ ˜ ) h∇ ˜ 2 ) h2∇ 2 h
(33)
Substituting (26) and (30) into (28) gives
ω′ )
[
h W(b ˜r ) κ
∫ Φ(ξ) dξ - κh1 ∫ ∫ Φ(ξ′′) dξ′′ dξ′] ∞
0
ξ
∞
ξ′
∞
(34)
The second term in brackets in (34) yields a small contribution, O(1/κh). Hence, the assumption made while deriving (27) was correct. Since ω′ is defined by the first term in brackets in (34), the length scale for changes in ω′ is h. Due to the inequality 18, in the close vicinity of the perimeter line, one can set ω′ to be independent of the coordinate ξ (Figure 1b). The second term in (34) was important for obtaining boundary condition 32 from (15b), but it can be omitted while substituting (34) into (27). The final expression for the parameter d, which is present in eq 1, is derived by combing (14), (26), (27), and (34). Ignoring the terms of order O(1/κh), the parameter d, is represented in the form suggested by Griffiths and Nilson15 for the circular and plane-parallel cross sections and given by eq 2 of the present communication. Using eq 1, the result can be written in terms of the dispersion coefficient, D/, as
D/ ) D(1 + δDeo)
(35a)
The electroosmotic contribution to the dispersion coefficient, δDeo, is given by
δDeo ) K(ζE/ηDκ)2g(ζ˜ )
(35b)
where the function g(ζ˜ ) is independent of the cross-section geometry and is given by
Figure 2. Function g(ζ˜ ) for different electrolytes: CI signifies coion; CNT signifies counterion. (1) Binary electrolyte: |ZCI| ) 1; |ZCNT| ) 2. (2) Ternary electrolyte |ZCI| ) |ZCNT1| ) 1; |ZCNT2| ) 2; CCNT1 ) CCNT2. (3) Binary symmetric electrolyte |ZCI| ) |ZCNT| ) 1. (4) Ternary electrolyte |ZCI1| ) |ZCNT| ) 1; |ZCi2| ) 2; CCI1 ) CCI2. (5) Binary electrolyte: |ZCNT| ) 1; |ZCI| ) 2.
to the empirical interpolation function, gG-N(ζ˜ ), suggested by Griffiths and Nilson.15 For (Z:Z) electrolyte
gG-N(ζ) ) g(ζ˜ ) )
[∫ 1
ζ˜
∞
0
Φ(ξ) dξ
[
1
∫ ζ˜
]
2
) Φ
ζ˜
0
x∑ 2
νk[exp(-ΦZk) - 1]
k
dΦ
]
2
(36)
We used (25) to obtain the latter expression in (36). The crosssection geometrical factor K takes the form
K ) 〈W〉l - 〈W〉
(37)
where K ) 〈...〉l ) (1/l)Il ... dl is an average over the cross-section perimeter. The function W(b ˜r ), which should be substituted in (37), is obtained as a solution of the problem 31 and 32 where the normalization is given by (17) and (33). For given geometry, the problem 31 and 32 always has a family of solutions that differ by an additive constant. However, using (37) yields a unique result whichever solution of (31) and (32) is employed. DEPENDENCY ON ELECTROKINETIC POTENTIAL AND ON ELECTROLYTE TYPE According to eqs 35a and b, the electroosmotic correction of the apparent diffusion coefficient is proportional to ζ 2 and to the function g(ζ˜ ) given by (36). For different ionic systems, the behavior of g(ζ˜ ) is illustrated in Figure 2. For all electrolytes, at ζ˜ f 0, the function g(ζ˜ ) approaches unity. Curve 3 plotted for symmetric electrolyte (Z1 ) - Z2 ) 1, ν1 ) ν2 ) 1/2) practically coincides (maximal error is 1.3%) with the curve plotted according
2592 + 24(ζ˜ Z)2 2592 + 96(ζ˜ Z)2 + (ζ˜ Z)4
(38)
According to Griffiths and Nilson, eq 38 yields the common form of the function gG-N(ζ˜ ) for both circular and plane-parallel crosssection geometries. Agreement between (36) and (38) holds for an arbitrary symmetric electrolyte (Z1 ) - Z2 ) Z, ν1 ) ν2 ) 1/(2Z2)). For asymmetric electrolytes, the behavior of g(ζ˜ ) is quite different for different electrolyte systems. For the case of multicharged co-ions, the function g(ζ˜ ) has a maximum that is slightly larger than unity (curves 4 and 5). Up to the rather high electrokinetic potentials, the function takes a value close to unity (for example, at ζ˜ = 4, g(ζ˜ ) > 0.9). Another trend is observed for the case of multicharged counterions. For this case, g(ζ˜ ) monotonically decreases with an increase in ζ˜ and, at rather low electrokinetic potentials, takes values noticeably lower than unity. Thus, for a given geometry, the dispersion coefficient is not only defined by the electrokinetic potential, ζ˜ , and by the Debye parameter, κ, but also depends on the electrolyte content. For the cases displayed in Figure 2, at ζ˜ < 0.5, such variations of g(ζ˜ ) and, hence, of δDeo, do not exceed 10%. However, at higher electrokinetic potentials, the dependency of δDeo on the electrolyte type strongly manifests itself. For example, for the case of ζ˜ = 2, values of g(ζ˜ ) and, hence, δDeo evaluated for 1:2 and 2:1 electrolytes differ nearly by a factor 2 (points G1 and G2 in Figure 2). GEOMETRICAL FACTOR K FOR DIFFERENT CROSS SECTIONS Now, we will consider how the difference in the shape of the channel cross section can affect the geometrical factor K in (35b). Analytical Chemistry, Vol. 75, No. 4, February 15, 2003
905
Figure 3. Annulus: sketch of cross section (a) and variation of geometrical factor K with ratio of external to internal radii, γ ) a1/a2 (b).
For any cross-section geometry, the procedure given by (31), (32), and (37) can be conducted numerically. For some geometries, it is easy to obtain the factor K analytically. Below, we consider some examples. Annulus. This type of cross-section geometry corresponds to a channel whose walls are the coaxial circular cylinders with radii a1 and a2 (a1 > a2) (Figure 3a). For this case
h ) (a1 - a2)/2
(39)
A solution for (31) and (32) takes the form
W(b ˜r ) )
|
|
|
|
1b b 2 2γ ˜r - ˜r c ln b ˜r - b ˜r c 4 (γ - 1)2
(40)
Figure 4. Circumscribed polygon: illustration to analysis of eqs 43 and 44.
For all the circumscribed polygons, a solution for (31) and (32) takes a surprisingly simple form
W(b ˜r ) )
|
1 b b ˜r - ˜r c 4
|
2
(43)
Here, brc is the radius vector of the inscribed circle center. It is trivial to show that (43) satisfies eq 31. Using the illustration in Figure 4, one can show that (43) satisfies boundary condition 32 as well. It can be concluded that (rb - brc)‚n b ) |rb - brc| cos θ ) ac. B ˜ W ) (b ˜r - b ˜r c)/2 ) (rb - brc)/2h. Combining At the same time, ∇ (41) with both the latter equalities, one obtains (32). Applying (37) to (43) yields the geometrical factor K for an arbitrary circumscribed N-gon. N
where b ˜r c ) b r c/h, b r c is the radius vector of the ring center; γ ) a1/a2. Making use of (37), after some straightforward transformations, one obtains an expression for K:
[
1 4γ2 2 K) γ 4γ + 1 + ln γ 2(γ - 1)2 γ2 - 1
]
where ac is the inscribed circle radius. 906
Analytical Chemistry, Vol. 75, No. 4, February 15, 2003
3
K)
1 2
i
+
i)1
∑cot(R /2) i
i)1
(41)
(42)
(44)
N
6
The curve given in Figure 3b was plotted according to (41) to display the variation of the geometrical factor K for the annulus cross section as a function of γ ) a1/a2. When γ increases from 1 to ∞, the coefficient K monotonically increases from 1/3 to 1/2. The limiting cases γ f 1 and γ f ∞ correspond to the planeparallel and circular cross sections, respectively. Thus, the result agrees with the limiting cases analyzed by Griffiths and Nilson,15 who obtained K ) 1/2 for the circular cross section and K ) 1/3 for the plane-parallel cross section. Circumscribed Polygon. Now, we will consider cross sections whose perimeter line is a closed polygon, which can be circumscribed around a circle. This class comprises all triangles, the tetragons with equal sums of the opposite side lengths, all regular N-gons, etc. For all such cross sections
h ) ac/2
∑cot (R /2)
where Ri signifies the ith internal angle of the N-gon (Figure 4). While deriving (44), we made use of the identity N (ac/l)∑i)1 cot(Ri/2) ) 1/2. Let us now consider some particular geometric cases described by eq 44. Regular N-gon. Making use of the substitution Ri ) R ) π(N - 2)/N and after some manipulations, one obtains from (44)
K)
1 1 π + tan2 2 6 N
()
(45)
When N f ∞ one obtains a circle, and eq 44 yields K ) 1/2, which coincides with the above given result and with result of Griffiths and Nilson.15 Some particular geometric cases inherent in eq 45 are given Table 1. Rhombus. For such a circumscribed tetragon N ) 4, R1 ) R3 ) R, and R2 ) R4 ) π - R (Figure 5a). Consequently, eq 44 is reduced to
K ) 2/(3 sin2(R))
(46)
Figure 5. Rhombus: sketch of cross section (a) and variation of geometrical factor, K, with internal angle R (b). Figure 6. Triangle: sketch of cross section (a) and variation of geometrical factor, K, with one of the internal angles, R1, for different values of another angle, R2 (b).
Table 1. Values of K for Some Cross-Section Geometries type of geometry plane-parallel ((41), at γ f 1) circle ((41), at γ f ∞) regular N-gon, (45) circle (N f ∞) hexagon (N ) 6) pentagon (N ) 5) square (N ) 5) triangle (N ) 3) triangle, (47) R1 ) π/2; R2 ) R3 ) π/4 R1 ) π/9 R2 ) R3 ) 4π/9 R1 ) π/36 R2 ) R3 ) 35π/72 rhomb, (46) R ) π/18 R ) π/36 rectangle, (51)
K 1/
3
1/
2
1/ 2 5/9 0.588 2/ 3 1 4/3 4.344 80.67 22.11 87.76 2/ 3
Variation of K with the rhombus internal angle, R, is displayed by the curve in Figure 5b. It is remarkable that, always, K g 2/3. The equality is satisfied when R ) π/2 (point A1), i.e., for a square. At sufficiently small (but absolutely realistic) angles R, the geometrical factor K substantially exceeds unity. The latter can be observed from the examples given in Table 1. In Figure 5b, at R ) 0.08π, we observe that K = 10. Triangle. In this case, N ) 3 and R3 ) π - (R1 + R2). Consequently, (43) is simplified to
K)
3 3 3 1 1 cot (R1/2) + cot (R2/2) + tan [(R1 + R2)/2] + (47) 2 6 cot(R1/2) + cot(R2/2) + tan[(R1 + R2)/2]
For a triangular cross section, eq 47 represents the geometrical factor K as a function of two angles R1 and R2. Applying the standard procedure, it can be shown that, at R1 ) R2 ) π/3, the function K (R1, R2) has an absolute minimum of K ) 1. Hence, for a triangular cross section, K g 1, where the equality corresponds to the equilateral triangle. Also, for a given R2, the dependency of K on R1 has a minimum, which corresponds to R1 ) (π - R2)/2. The latter signifies that, for a given R2, the minimum of K is achieved when the two other angles are equal, i.e., for isosceles triangle. The curves given in Figure 6 display the variation of K with R1 for three fixed values of R2. When R1 f 0, all the curves approach a common asymptotic curve given by K = 2/3R12. Each of the curves has a minimum and is symmetrical with
reference to a vertical line passing through the minimum (R1 ) (π - R2)/2). All the curves are located above the line given by K ) 1 to which the curve corresponding to R2 ) π/3 becomes tangent at the point R1 ) π/3. Similar to the case of a rhombus, for the triangular cross section, the geometric factor K can be much higher than unity (Table 1). Generally, with the help of (44), one can conclude that K exceeds unity for any circumscribed polygon that contains a sufficiently small internal angle. For a polygonal cross sections that have a small internal angle and, hence, are characterized by a geometrical factor K . 1 (for example see Table 1), one can expect a substantial increase in the dispersion coefficient. Recall that, according to (35b), the electroosmotic contribution into the dispersion coefficient, δDeo, is proportional to K. It should be noted that there is a restriction in using the present approach for small angles. In the close vicinity of the polygon corners, the approximation of thin double layer fails. Under the condition κac . 1, the role of the regions around the corners can be ignored. When one of the angles (Ri) is small, one can evaluate the inscribed circle radius, as ac = l tan(Ri/2)/2 = lRi/4. Combining this estimation with the above given inequality we obtain
Ri . 4/κl
(48)
While dealing with the polygonal cross sections, the inequality of (48) yields the condition for applicability of our approach. Combining (46) and (48), we conclude that large geometrical factors K, which are obtained using the approach of the present paper, have to satisfy a strong inequality K , (κl)2/24. For example, for the 0.001 N (1:1) aqueous electrolyte solution (κ ≈ 108 m-1) and for l ) 10-5 m, one obtains that K , 4 × 104. The present theory is applicable for observation times, t, which are sufficiently long to form a quasi-stationary concentration distribution within a channel cross section. For the cross sections containing sharp angles, the observation times need to satisfy the inequality t . l 2/4D, where a half of the cross-section perimeter (l/2) was chosen as the largest length scale parameter describing the concentration distribution. Accordingly, predictions of high geometrical factors, K can be linked to rather long times. For Analytical Chemistry, Vol. 75, No. 4, February 15, 2003
907
example, when l ) 10-5 m and D ) 10-10 m2/s, the observation time has to satisfy the inequality t . 0.25 s. Rectangle. When the cross section is a rectangle with the sides a and b, the parameter h is given by
h ) ab/2(a + b)
(49)
We will use the 2D Cartesian coordinates where the X- and Y-axes are parallel to the sides having lengths a and b, respectively. In such a coordinate system, a solution of the problem 31 and 32 can be represented as
W)
b a 1 + (y˜ - ˜y c)2 (x˜ - ˜x c)2 2 a+b a+b
[
(
)]
(50)
where ˜xc and ˜yc are the normalized Cartesian coordinates of the diagonal intersection. Obtaining the geometric factor, K, with the help of (37), we arrive to the surprising result
K ) 2/3
(51)
Thus, for the rectangular cross section, the coefficient K is constant and, hence, does not depend on a and b. This result matches with that obtained for a square (Table 1). However, in the limit a/b f 0, eq 51 does not yield the result of K ) 1/3. The latter result of Griffiths and Nilson15 for the plane-parallel cross section was confirmed in the present paper by making use of the limiting transition γ ) a1/a2 f 1 in eq 41, which was derived for the case of annulus. For dispersion in the Poiseuille flow in a rectangular channel, a similar situation was discussed by Chatwin and Sullivan.33 Assuming that the ratio of the rectangular sides approaches zero, they obtained an asymptotic expression that differs from the earlier result describing the dispersion coefficient in the plane-parallel cross section. It is interesting to note that the result of K ) 1/3 can be obtained using the limiting transition a/b f 0 in (49) and (50). The latter yields h ) a/2 and W ) (x˜ - ˜xc)2/2 + O(a/b). The procedure given by (37), being applied to the zero-order term, W ) (x˜ - ˜xc)2/2, yields 1/3. However, while obtaining the mean values in (37), the contribution of the first-order term, O(a/b), becomes of the same order as that associated with the zero-order term, and hence, this term should not be omitted. Therefore, for the rectangular cross section, a properly conducted Aris-Taylor procedure does not yield the result of K ) 1/3 at a/b f 0. CONCLUSION In the present paper, we suggested an approximate description of the electroosmotic dispersion in the microchannels having a thin electric double layer at the wall/electrolyte solution interface. The result given by (35) can be represented in terms of the plate height, H ) 2D//〈u〉, which is widely used in the literature6-15 dealing with electroosmotic dispersion
H)
2Dη 2ζE g(ζ˜ )K + ζE ηκ2D
(33) Chatwin, P. C.; Sullivan, P. J. J. Fluid Mech. 1982, 120, 347-358.
908
Analytical Chemistry, Vol. 75, No. 4, February 15, 2003
(52)
Figure 7. Variation of plate height with electric field strengths for two cross-section geometries and different electrolyte systems.
The function g(ζ˜ ) is given by eq 36. Accordingly, it does not depend on the cross-section geometry but depends on the electrolyte type. In the case of a symmetric electrolyte, g(ζ˜ ) is well approximated by empirical function suggested by Griffiths and Nilson.15 Using the function g(ζ˜ ), one can describe channels with arbitrary electrokinetic potential, ζ. For each given cross-section geometry, the dimensionless geometrical factor K is determined using any solution of the problem given by (31) and (32) and the procedure given by (37). The examples considered in the present paper show that, for different cross-section geometries, K takes noticeably different values. For the case of a symmetric electrolyte and circular or planeparallel cross sections, the present theory predictions coincide with the relevant limiting cases analyzed by Griffiths and Nilson.15 For other electrolyte type and for different cross-section geometries, the present report yields the plate height, H, which varies in a wide range. To illustrate the possibility of large variations of H, in Figure 7, we present two curve families plotted for the rectangular and rhomboidal cross sections. The curves display variation of plate height with the electric field strength. Each family contains three curves plotted for 1:2, 1:1, and 2:1 electrolytes, respectively. The curves given in Figure 7 are qualitatively similar: each of them contains a minimum. Such a behavior was discussed in the literature.8,11,12 The minimum exists since, at low electric fields, the dispersion due to molecular diffusion dominates, whereas at high electric fields, the electroosmotic dispersion dominates. Despite the qualitative similarity, we observe a substantial quantitative difference in the plate height, H, evaluated for the same set of parameters but for different electrolyte and crosssection geometry. For example, at E ) 1 kV/cm, the plate height for a 2:1 electrolyte solution in a rhomboidal channel (point B2) is approximately in 30 times more than that for a 1:2 electrolyte solution in a rectangular channel (point B1).
Thus, the longitudinal spreading, which occurs in microchannels due to electroosmotic flow, strongly depends on the electrolyte type and channel cross-section geometry. For any given system, this dependency can be modeled using the results of the present paper. ACKNOWLEDGMENT The authors gratefully acknowledge the financial support from the Natural Sciences and Engineering Research Council of Canada. GLOSSARY
b, r brc
radius vector
S
cross-section area
T
absolute temperature
u
local velocity of liquid
x, y
Cartesian coordinates
Zk
electric charge of the kth ion expressed in the Faraday units
Greek Letters Ri
internal angle of polygon
β
double-layer effective length introduced by Gas et al.13
dielectric permittivity
Φ
normalized solution of the Poisson-Boltzmann problem for semi-infinite region ratio of external to internal ring radius (γ ) a1/a2)
a
radius
a1,2
external and internal radiuses of a ring
ac
radius of the inscribed circle
γ
Ck
concentration of kth ion in equilibrium electroneutral solution
η
liquid viscosity
κ, κ-1
Debye parameter and Debye length
D
molecular diffusion coefficient
F
volumetric density of electric charge
D/
apparent diffusion coefficient due to the longitudinal dispersion
ω, W, w
functions determined within the framework of the Taylor-Aris procedure
d
length-scale parameter determined using Taylor-Aris procedure
ξ
normalized coordinate
Ψ
equilibrium electric potential
E
electric field strength magnitude
ζ
electrokinetic potential
F
Faraday constant
g(ζ˜ )
function of the electrokinetic potential
H
plate height
h
length-scale parameter (h ) S/l)
bıx,y
unity vectors
K
geometrical factor
L
channel length
l
cross-section perimeter
Miscellaneous ∇ B ∇2
2D gradient operator )∇ B ‚∇ B 2D Laplace operator
〈‚‚‚〉
average over cross-section area
〈‚‚‚〉l
average over cross-section perimeter line
N
number of polygon sides
b n
the unity vector normal to the channel wall
Received for review May 30, 2002. Accepted November 11, 2002.
R
gas constant
AC0203591
Analytical Chemistry, Vol. 75, No. 4, February 15, 2003
909
