Modeling of Sample Dynamics in Rectangular Asymmetrical Flow

Anal. Chem. 2000, 72, 4331-4345

Modeling of Sample Dynamics in Rectangular Asymmetrical Flow Field-Flow Fractionation Channels S. A. Suslov* and A. J. Roberts†

Department of Mathematics and Computing, University of Southern Queensland, Toowoomba, Queensland 4350, Australia

We model the evolution of the concentration field of macromolecules in a rectangular asymmetrical flow fieldflow fractionation channel using center manifold techniques. The deviation of the primary flow from a parabolic profile influences the concentration field and this is investigated to complement previously known results. The long-term evolution of the components of the sample is shown to be well described by a one-dimensional advection-diffusion equation. The coefficients of this equation are determined by the rigorous analysis of the complete set of equations governing the two-dimensional fluid flow. This model gives quantitative predictions of the elution time of the samples, the width of the concentration pulse, and the resolution of the apparatus. The influence of initial sample width and effects of the secondary relaxation from focusing to elution conditions are discussed. Reported theoretical predictions are in agreement with experimental results published previously. Flow field-flow fractionation (FFFF) channels of different designs and geometries (e.g., refs 1-4) are used in order to separate different species that are present in a solvent during its elution. The effectiveness of separation depends not only on differences in the diffusion coefficients of eluted components but also on the flow conditions in the channel. Although an approximate theory of fractionation has been developed over the past few decades (e.g., refs 2-7) the detailed analytical description of these flows is rather involved mathematically. Until recently, the refinement of existing models was a challenging task requiring complicated algebra.5 Nowadays this limitation is significantly reduced by the use of symbolic computer software to perform algebraic computations. In addition, the development of the center * Corresponding author: (e-mail) [email protected]. † E-mail: [email protected]. (1) Litze´n, A. Anal. Chem. 1993, 65, 461-470. (2) Hovingh, M. E.; Thompson, G. H.; Giddings, J. C. Anal. Chem. 1970, 42 (2), 195-203. (3) Giddings, J. C.; Yang, F. J.; Myers, M. N. Anal. Chem. 1976, 48 (8), 11261132. (4) Wyatt, P. G. J. Colloid Interface Sci. 1998, 197, 9-20. (5) Giddings, J. C.; Yoon, Y. H.; Caldwell, K. D.; Myers, M. N.; Hovingh, M. E. Sep. Sci. 1975, 10 (8), 447-460. (6) Wahlund, K.-G.; Winegarner, H. S.; Caldwell, K. D.; Giddings, J. C. Anal. Chem. 1986, 58, 573-578. (7) Giddings, J. C. In Treatise on Analytical Chemistry; Kolthoff, I. M., Elving, P. J., Eds; Wiley: New York, 1981; Part I, Vol. 5, Chapter 3, pp 63-164. 10.1021/ac9915022 CCC: $19.00 Published on Web 08/12/2000

© 2000 American Chemical Society

manifold theory8-10 and its applications to the dispersion in channels11,12 provide a systematic way to formalize iterative algorithm for deriving accurate asymptotic models of dynamics such as those found in FFFF channels. In this work, we take advantage of computer algebra to confirm previously known results from the approximate theory of FFFF channels, refine them, and obtain new quantitative information on rectangular asymmetrical FFFF channels. In practice, to accelerate the separation process, trapezoidal asymmetrical FFFF channels are frequently used (e.g., ref 1) Analysis of their operation could be done in a manner very similar to the one discussed in the following sections, but the results would be algebraically more involved due to the necessity of including the third spatial dimension. Thus, in this work we, analyze a simpler rectangular geometry, concentrating on new methodological concepts without being overwhelmed by an excessive algebraic complexity. Analysis of the trapezoidal configuration is left for future studies. Consider the transport of some chemical substance in the fluid flow of an asymmetrical FFFF channel as analyzed by Litzen and others (e.g., refs 1 and 13) and sketched in Figure 1. The horizontal plate forming the upper wall of the channel is impermeable to both the eluted molecules and the carrier fluid. The lower wall is impermeable to the eluted molecules but allows a crossflow of the carrier fluid to be pumped out. This cross-flow distributes the sample preferentially to the lower side of the channel (the accumulation wall) as shown in detail Figure 7. The resulting asymmetric distribution of sample concentration c(x,z,t) creates a differential advection of species with different diffusivities. Thus, different species reach the end of the elution channel at different times. This is used for the separation of different species in the sample.1,13,14 In previous work,15 we used center manifold techniques to analyze the process of separation of different species in a symmetrical FFFF channel. On the basis of a derived advectiondiffusion equation, which accurately models the transport of the (8) Carr, J. Applications Of Centre Manifold Theory; Applied Math. Sci. 35; Springer-Verlag: New York, 1981. (9) Carr, J.; Muncaster, R. G. J. Differential Equations 1983, 50, 260-279. (10) Carr, J.; Muncaster, R. G. J. Differential Equations 1983, 50, 280-288. (11) Mercer, G. N.; Roberts, A. J. SIAM J. Appl. Math. 1990, 50, 1547-1565. (12) Mercer, G. N.; Roberts, A. J. Jpn. J. Ind. Appl. Math. 1994, 11, 499-521. (13) Litze´n, A.; Wahlund, K.-G. Anal. Chem. 1990, 62, 1001-1007. (14) Wahlund, K.-G.; Giddings, J. G. Anal. Chem. 1987, 59, 1332-1339. (15) Suslov, S. A.; Roberts, A. J. J. Math. Chem. 1999, 26, 27-46.

Analytical Chemistry, Vol. 72, No. 18, September 15, 2000 4331

sample, we identified the governing nondimensional parameters of the flow for optimal separation. We deduced that the resolution of a symmetrical FFFF channel is the highest for relatively large cross-flow rates. One of the most difficult problems in manufacturing FFFF channels is ensuring uniform permeability of the walls.14 An asymmetrical FFFF channel has just one porous wall and consequently is easier to make than a symmetrical one. Thus, in this work, we analyze the flow in an asymmetrical FFFF channel. Typical operation of the asymmetrical FFFF channel consists of four stages: injection of the sample; its longitudinal focusing during a period of counterflow and simultaneous cross-channel relaxation to the nearly exponential distribution; quick secondary relaxation to elution conditions after the counterflow is switched off and the sample starts drifting along the channel; and finally the elution itself. The injection procedure and, in particular, the injection flow rate create the resulting initial sample size. It is usually impossible to know the actual distribution of the components in the sample immediately after injection,1 but if the following focusing stage is long enough then it approaches a narrow Gaussian in the direction along the channel. The time required to reach this completely focused sample state is estimated in section 3. We use analytical techniques based upon center manifold theory16 to model the dynamics of the sample. First, in section 1 we solve the steady Navier-Stokes equations to find asymptotic expression for the fluid flow in the channel. Then, in section 2, we analyze the two-dimensional advection-diffusion for a sample to deduce that a model for its evolution in the channel is the onedimensional advection-diffusion equation

[

]

∂〈c〉 ∂〈c〉 ∂ -U〈c〉 + De ) ∂t ∂z ∂z

(1)

where t denotes time, z measures distance downstream along the channel, and 〈c〉 ) ∫w0 c(x,z,t) dx is the local amount of sample substance per unit channel length. For both focusing and elution stages we derive expressions (35) for the effective advection velocity U, as it predominantly determines the time required for focusing and the time of efflux of the sample out of the channel, and the effective diffusivity De, as it determines how wide the sample spreads by the time it reaches the end of the channel. We also study the influence of the primary parabolic flow profile distortion due to the cross-flow on the distribution of sample components across the channel, which is found to be weak for the most practical regimes. Analysis of model (1) in section 2 shows that in general the resolution of the asymmetrical FFFF channel increases with the strength of the cross-flow, but so does the retention time. A desirable balance between these two factors should define the optimal operational regime. Comparison of theoretical results with experiments shows the reliability and accuracy of the developed model (1). It is subsequently used to provide detailed information on the actual distribution of the sample components in an asymmetrical FFFF channel. Finally, in section 3, we provide quantitative information on sample relaxation during the transition between focusing and elution stages. This issue has never previously been discussed. We show that the width of the initial sample always decreases (16) Roberts, A. J. J. Aust. Math. Soc. B 1988, 29, 480-500.

4332

Analytical Chemistry, Vol. 72, No. 18, September 15, 2000

Figure 1. Side view of asymmetrical flow field-flow fractionation (FFFF) channel.

when the cross-flow is stronger during the elution stage than during the focusing. If the cross-flow weakens during the switching from focusing to elution, then the sample will either shrink or spread depending on its initial size and the ratio of flow velocities. 1. FLUID FLOW IN A RECTANGULAR ASYMMETRICAL FFFF CHANNEL Consider a rectangular asymmetrical FFFF channel as discussed in refs 1 and 13 and depicted schematically in Figure 1. The dynamics take place between two flat plates located at x* ) 0 and x* ) w. The fluid flow between the plates is driven predominantly by a pressure gradient p*z parallel to the plates. Being that of a Newtonian fluid with kinematic viscosity ν and density F, the velocity field is close to that of parabolic Poiseuille flow except that there is a cross-flow, of velocity -u0, through the accumulation wall (u0 is positive). The cross-flow varies with height x but is assumed to be independent of longitudinal distance z. (The assumption of cross-flow uniformity is valid if the pressure drop across the accumulation wall is much larger than 12Fνv0L/ w2sthe estimated pressure drop along the channel. Experimental study by Litzen1 shows that this is generally the case.) The accumulation wall is permeable to the fluid in order for this crossflow to occur, but it is impermeable to the sample molecules. Within the fluid, the sample substance of concentration c(x,z,t) is advected by the flow and diffuses with coefficient D. In this section, we nondimensionalize the governing differential equations, deduce the advecting fluid velocity field, and confirm that it is nearly parabolic. For order of magnitude estimates of quantities we use the geometry of Litzen.1 The density of the carrier fluid (water) is F ) 1 g/cm3 and the kinematic viscosity ν ≈ 0.01 cm2/s. When the sample molecules are the cow pea mosaic virus (ref 1, p 464), this configuration gives parameters as listed in Table 1. The equations governing the fluid motion are the NavierStokes and continuity equations

∂q* 1 + q*‚∇q* ) - ∇p* + ν∇2q* ∂t* F

(2)

∇‚q* ) 0

(3)

for the incompressible fluid velocity field q* ) u*i + v*j and for the pressure p* (asterisks denote dimensional field quantities). The boundary conditions on the plates are those of no longitudinal flow,

v* ) 0 at x* ) 0 and x* ) w

(4)

Table 1. Typical Set of Physical Parameters for FFFF and the Consequent Parameters Appearing in the Analysisa parameter

v ) 0 at x ) 0 and x ) 1/λ

u ) -1 and u ) 0 at x ) 0 and x ) 1/λ, respectively (12)

value

channel thickness channel length distance to the focusing point kinematic viscosity mean inlet velocity cross-flow velocity molecular diffusivity boundary layer (BL) thickness cross-BL diffusion time downstream advection distance Schmidt number retention parameter velocity ratio ratio of outlet and inlet flow rates

w L z/st ν v0 u0 D l τ ξ Sc λ K 1-M

0.013 cm 28.5 cm 2.3 cm 0.01 cm2/s 2 - 10 cm/s 0.001 - 0.005 cm/s 2 × 10-7 cm2/s (1 - 4) × 10-5 cm 0.01 - 0.1 s 0.02 - 1 cm (2.5 - 5) × 104 0.003-0.015 (1-5) × 10-4 0.1-0.7

a The data are for the cow pea mosaic virus used in the experiments by Litzen (ref 1, p 464).

where λ ) D/(u0w) is the retention parameter and K ) u0/v0 is the ratio of the characteristic cross-flow and longitudinal velocities. Integrating the continuity equation (10) across the channel and using boundary conditions (12), we obtain

d〈v〉/dz + λ ) 0

u* ) 0 at x* ) w

(5)

and specified flow through the accumulation wall,

u* ) -u0 at x* ) 0

v ) 6〈v〉λx(1 - λx) +

2

(7)

0

where l ) D/u0 is the characteristic thickness of the distribution of sample substance in a boundary layer near the accumulation wall, τ ) l/u0 ) l2/D ) D/u02 is the cross-boundary layer advection time (l/u0) or equivalently the cross-boundary layer diffusion time (l2/D), v0 is the mean fluid speed at the entrance to the channel (at z ) 0), ξ ) v0τ is the mean downstream fluid displacement in a cross-boundary layer diffusion time, and Sc ) ν/D is the Schmidt number. Typical values of all these quantities are recorded in Table 1. Substituting these scalings into equations (2)-(6) for an assumed steady fluid flow we obtain

∂u ∂2u 1 ∂u 1 ∂p ∂2u +v + 2+K 2 2 u )Sc ∂x ∂z K ∂x ∂x ∂z 2 2 ∂p 1 ∂v ∂v ∂v ∂v u +v )+K 2 2 + Sc ∂x ∂z ∂z ∂x2 ∂z ∂u ∂v + )0 ∂x ∂z

(

)

(

with boundary conditions

)

〈v〉x (1 - λx)(28λ4x4 - 56λ3x3 + 70Sc

49λ2x2 - 91λx + 38) + O((λSc)-2) (15) u ) -(1 - λx)2(1 + 2λx) -

λx2 (1 - λx)2(4λ3x3 - 6λ2x2 + 70Sc 5λx - 19) + O((λSc)-2) (16)

p*

z* t* u* x* v* , z) , t) , u) , v) , p) l ξ τ u0 v0 Sc Fv

(14)

where ∫1/λ 0 v1(x) dx ) 1/λ as 1/λ is the nondimensional thickness of the channel. Assuming a large Schmidt number (see Table 1), we find the solution of the Navier-Stokes equations (9) and (10) as a series in powers of 1/Sc:

(6)

The nondimensionalization we adopt is chosen to reflect the fact that, for the regime of most effective separation of species, the sample is concentrated near the accumulation wall due to the cross-flow.15 Introduce the following nondimensional variables:

(13)

Thus, the mean flow 〈υ〉 ) 1 - zλ in the chosen nondimensionalization. The nondimensional length of the channel is L/ξ and consequently in order to maintain nonzero fluid flux at the end of the channel the condition L/ξ < 1/λ or equivalently the fraction of fluid lost in the cross-flow M ) Lu0/(wυ0) < 1 must be satisfied. The solution for the fluid flow is found in the following functional form:

v ) 〈v〉υ1(x), u ) u(x), p ) p1(z) + p2(x)

no flow through the upper wall,

x)

(11)

(8) (9) (10)

dp1 81λ〈v〉 ) -12λ2〈v〉 + + O(Sc-2) dz 35Sc

(17)

dp2 λK ) 6λ2K (1 - 2λx) (504λ5x5 - 1260λ4x4 + dx 35Sc 840λ3x3 - 81λx - 19) + O(Sc-2) (18) Note that the above expressions are valid when 1/(λSc) ) u0w/ν is small, that is, for very thin channels or very viscous fluids. Thus, the leading terms in the above expressions correspond to the Stokes limit when the inertia terms in the momentum equations are neglected. Our leading order approximation of the transverse velocity component, u ≈ -(1 - λx)2(1 + 2λx), agrees with eq 1 reported by Wahlund and Giddings.14 The downstream velocity profile is nearly parabolic (solid line in Figure 2a). Corrections of O((λSc)-1) are due to the coupling with the cross-flow. They distort the Poiseuille flow by accelerating it near the accumulation wall and decelerating away from it (dashed line in Figure 2a), but do not affect its mean velocity. The vertical velocity is slightly decelerated in the middle part of the channel because of the coupling (the correction shown by the dashed line in Figure 2b is always positive whereas the primary vertical velocity shown by the solid line is always negative). Higher-order corrections in Analytical Chemistry, Vol. 72, No. 18, September 15, 2000

4333

O(10-3) (see Table 1), and we expect cross-channel advection to keep the sample close to the bottom plate. 2.1. Nearly Exponential Distribution of Sample Substance. It is useful to discuss first the steady distribution of the sample material suggested by Wahlund and Giddings (ref 14, eqs 3-5) and Litzen (ref 1, eq 1). In our notation, their expressions for the concentration are written as

1 c ) C exp -x + λ2x3 - λ3x4 2

(

Figure 2. Primary flow velocities (solid lines) and their distortions due to coupling (dashed lines). For plotting purpose, the distortion profiles are multiplied by a factor 10λSc.

powers of 1/(λSc) can be easily obtained by executing the symbolic algebra code, but their contribution is negligible for practical conditions. The form of the velocity field corrections shown in Figure 2 suggests that in realistic conditions the crosschannel relaxation of the sample will require slightly longer time while elution will be slightly quicker than predicted by any approximate theory that neglects the coupling between primary and cross-flows. An important conclusion can be made from expression (15) in the limit of very strong cross-flow λ f 0 and large, but finite Schmidt number. Then in the concentration boundary layer (x of the order 1), the leading term for the longitudinal velocity is v ∼ 〈υ〉(19x)/(35Sc) while it is of the order λ for a symmetrical FFFF channel. It will be shown later that this fact leads to the qualitative difference in retention ratio for an asymmetrical flow channel in contrast to a symmetrical one. 2. SAMPLE DYNAMICS We consider the dynamics of the sample in the steady flow field (eqs 15 and 16). Here we show how to model these dynamics by the one-dimensional advection-diffusion equation (1). Inside the channel, the sample evolves according to the nondimensional two-dimensional advection-diffusion equation

∂c ∂c ∂2c ∂c +v +u ) +K ∂t ∂z ∂x ∂x2

2

∂2c ∂z2

(19)

where c ) c(x,z,t) is the sample concentration field. Herein we assume that the molecules are neutrally buoyant. Equation (19) is a differential form of the mass conservation principle for the eluted substance and accounts for the flux given by Fick’s first law of diffusion and flux associated with migration of the molecules with the carrier fluid. The nondimensional boundary conditions for the sample substance are those of no flux through the channel walls

c+

∂c ∂c ) 0 and ) 0 at x ) 0 and x ) 1/λ, respectively ∂x ∂x (20)

The main nondimensional parameter λ, appearing as the inverse nondimensional thickness of the channel, is typically small, 4334


)

(21)

where C is the concentration at the accumulation wall. It is straightforward to check by direct substitution that if C ) const, then this distribution does not satisfy (19) and, consequently, violates the sample mass conservation principle. There is a simple physical reason for this. The substance is carried by a decelerating fluid which is removed through the porous accumulation wall. Since the sample material cannot penetrate through the accumulation wall, it must necessarily build up; i.e., its concentration must increase with time and/or z. A steady z -independent distribution of eluted material in an asymmetrical FFFF channel is physically impossible. Consequently, one must allow for the variation of C along the channel when developing any approximation for asymmetrical channels: if C ) C(z) and in addition K is very small, then by substituting (21) in (19) we obtain that this equation is satisfied to error of order K 2 if

C(z) )

C0 〈v〉

)

C0 1 - λz

(22)

where C0 is the accumulation wall concentration at the entrance to the channel. (Discussion of this distribution was suggested by the reviewer.) Therefore, approximation (21) for the cross-channel sample distribution is relatively accurate only for a very special case when the concentration at the accumulation wall grows according to (22). This is a poor approximation in practice because samples are injected near the entrance in the channel and so their concentration distribution is typically localized in the z direction. Consequently, now we proceed to develop formulas describing realistic distributions of the sample material. 2.2. Derivation of Advection-Diffusion Model. Here we summarize briefly the computer algebra algorithm developed earlier.15 Under the action of the cross-flow balanced by diffusion, the sample distribution near the accumulation wall approaches the exponential distribution c ) C0 exp(-x), which follows from eqs (21) and (22) in the limit of small λ. The shear velocity, different at different x, will smear this sample substance cloud out along the channel, while cross-flow and diffusion continue to act to push the cross-channel distribution toward the exponential distribution. The net effect is that the cloud has a concentration that is slowly varying along the channel and is approximately exponential across it. Thus, after the quick decay of cross-stream transients, we justify the relatively slow long-term evolution of a sample substance cloud for which z derivatives of C, ∂nC/∂zn, are small and C is generally different from (22) so that (21) is not exactly applicable. The application of center manifold theory to dispersion in channels and pipes has been developed by Roberts, Mercer, and

Watt.11,12,16-18 The essence of the technique is that in a FFFF channel we separate rigorously the slow transport along the channel from the relatively fast quasi-equilibrization across the channel. In these conditions, the local concentration field c(x,z,t) can be written as a function of the accumulation wall concentration C(z,t) and its along-channel derivatives; i.e., c ) h(C,x). All dependence on time and the longitudinal coordinate z is “absorbed” by C. Thus, evolution of the local field is parametrized by C so that

∂C/∂t ) g(C)

(23)

Center manifold theory puts these physical intuitive ideas on a rigorous foundation16 in the context of dynamical systems and ensures that such a model exists, retains major features of the original physical problem, and can be constructed by satisfying the governing equations to some order of accuracy (e.g., ref 8). The function h, C exp(-x) to leading approximation, describes the details of the sample field throughout space and time in terms of the concentration C of sample component at the accumulation wall. The function g represents the time evolution of the concentration at the accumulation wall. We approximate the functions h and g by assuming that the concentration field varies slowly along the channel; that is, ∂/∂z is a small operator. Formally we determine h and g in the form of the following asymptotic series ∞

g∼

∑

n)1

gn

∂nC ∂zn

∞

and h ∼

∑

hn(x)

n)0

∂ nC ∂zn

(24)

These series were shown to converge for the problems of dispersion in channels,11,12 and thus, we expect that the first few terms in the series representation of g give a model for the evolution of the sample in the FFFF channel. (For simplicity, the approximate solutions will be derived in the following sections in the form of the power series in retention parameter λ. They provide good accuracy for λ < 0.1.) To find the asymptotic expansions (24) we implement an iterative algorithm19 in computer algebra. (The computer algebra code can be obtained by e-mailing the authors at ssuslov@ usq.edu.au or [email protected].) Assume that some approximate solution of the advection-diffusion equation (19) with boundary conditions (20) is found in the center manifold form (23); for example, the iteration is initiated with the approximation c ) C exp(-x) and g ) 0. Each iteration refines such an approximation by finding a correction h′ to the shape of the center manifold and a correction g′ to the evolution thereon. As established by Roberts,19 the corrections are found by solving

∂2h′ ∂h′ ) R + g′exp(-x) + ∂x ∂x2

h′ +

respectively, and

h′ ) 0 at x ) 0

∫

1/λ

0

R + g′ exp(-x) dx ) 0

(28)

so that the boundary conditions (26) and (27) are satisfied. Then the differential equation (25) is solved to find the correction h′. The iterations continue until the desired terms in powers of small λ and small z derivatives are found in the asymptotic approximation. Computer algebra performs all the algebraic details. Note that the left-hand side of (25) is a simple approximation of the cross-channel mixing operator ∂2/∂x2 - u ∂/∂x entering the complete eq (19). The advantage of this iteration technique is that it does not require using the exact operator on the left-hand side of (25) since the accuracy of the solution is of the same order as the residual R which is driven to zero. The quality of approximation of the operator in the left-hand side of (25) affects only the number of iterations required to reduce the residual to the specified error level, but not the accuracy of the solution once this error level is reached. From the computer algebra results, the concentration field is to low order

{[

]

x2 c ) e-x C 1 + λ2x2(3 + x) - λ3 (12 + 4x + x2) - λ〈v〉 × 2 ∂C 2 x [3 - 2λ(3 + x) + 3λ2(54 + 18x + 6x2 + x3)] + λ2〈v〉2 × ∂z ∂2C x2 [72 + 24x + 9x2 - 12λ(60 + 20x + 6x2 + x3)] + 2 2 ∂z λ2K

2

∂2C x2 λx2 2 [6 + 2x λ(24 + 8x + x × )] + C 840Sc ∂z2 2

[76(3 + x) - 129λ(12 + 4x + x2) + 4λ2(4680 + ∂C x2 [228 ∂z 840Sc

1560x + 447x2 + 135x3 + 19x4)] - 〈v〉

516λ(3 + x) + 12λ2(2472 + 824x + 263x2 + 38x3) λ3(532152 + 177384x + 50358x2 + 10272x3 + 1169x4)] + ∂2C x2 [190(24 + 8x + 3x2) - 1670λ(60 + 20x + 2 700Sc ∂z

6x2 + x3) + λ2(4826880 + 1608960x + 460320x2 + 86268x3 + 11885x4 + 855x5)] + λK

where R is the residual of eq (19), with boundary conditions

(27)

This last boundary condition arises because we seek a solution parametrized by the concentration at the accumulation wall: C(z,t) ) c|x ) 0. The correction to the evolution g′ is determined from the solvability condition

λ〈v〉2 (25)

∂h′ ∂h′ ) 0 and ) 0 at x ) 0 and x ) 1/λ (26) ∂x ∂x

2

∂2C x2 [76(3 + x) ∂z2 420

129λ(24 + 8x + x2) + 4λ2(35244 + 11748x + 2433x2 + (17) Watt, S. D.; Roberts, A. J. J. Aust. Math. Soc. B 1994, 38, 101-125. (18) Watt, S. D.; Roberts, A. J. SIAM J. Appl. Math. 1995, 55 (4), 1016-1038. (19) Roberts, A. J. Comput. Phys. Commun. 1997, 100, 215-230.

420x3 + 38x4)] + O

(

∂3C 1 4 , ,λ ∂z3 Sc2

)}


(29) 4335

where the evolution of the sample concentration along the bottom plate is described to low order by

equivalently described by

λ ∂C ) 6λ2C(1 - 2λ + 30λ2) + C (19 - 129λ + ∂t 35Sc ∂C 1560λ2) - 6λ〈v〉 (1 - 2λ + 54λ2 - 516λ3) ∂z 〈v〉 ∂C (19 - 129λ + 2472λ2 - 44346λ3) + 35Sc ∂z

c ) e-x 〈c〉 1 - λ2(12 - 3x2 - x3) +

2

2λ ∂ C [6〈v〉2(38 - 835λ + 40224λ2) + 35Sc ∂z2

(

]

4x3 - x4) + λ〈v〉

∂2C [72λ2〈v〉2(1 - 10λ + 604λ2) + K 2(1 + 12λ2 ∂z2 48λ3 + 3312λ4)] +

{[

λ3 (72 - 12x2 2

∂〈c〉 [6 - 3x2 - 2λ(12 - 3x2 - x3) + ∂z

∂2〈c〉 3 3λ2(420 - 30x2 - 16x3 - 6x4 - x5)] - λ2〈v〉2 2 [144 2 ∂z 12x2 - 8x3 - 3x4 - 4λ(456 - 42x2 - 18x3 - 6x4 - x5)] λ2K

2

∂2〈c〉 [24 - 6x2 - 2x3 - λ(120 - 24x2 - 8x3 - x4)] + ∂z2

(

)

∂3C 1 λ4 K (19 - 258λ + 11748λ )] + O 3 , 2, , λ5 (30) ∂z Sc Sc 2

2

2.3. Conservative Form of the Model. The distribution of the sample material given by (29) and the evolution equation (30) are used to express the concentration field and evolution equation in terms of the experimentally observed quantity 〈c〉 ) ∫1/λ 0 c dx, which represents the local amount of the sample material per unit channel length. First, integrating (29) over the channel thickness gives

(

λ (38 -387λ + 35Sc

〈c〉 ) C 1 + 12λ2 - 36λ3 + 1296λ4 +

)

9888λ2) - 6λ〈v〉

∂C (1 - 4λ + 240λ2 - 2640λ3) ∂z

〈v〉 ∂C (19 - 258λ + 10380λ2 - 225480λ3) + 35Sc ∂z 2

C [3〈v〉2(7 - 84λ + 6828λ2) + 2K 2(1 - 5λ + 12λ ∂z2 2∂

2λ ∂2C 558λ )] + [6〈v〉2(133 -3507λ + 223224λ2) + 35Sc ∂z2 2

(

)

∂3C λ2 λ4 K (38 -645λ + 46188λ )] + O 3 , 2, , λ5 (31) ∂z Sc Sc 2

2

O

∂3〈c〉 λ 4 , ,λ ∂z3 Sc

)}

(33)

where for brevity we omit terms inversely proportional to the Schmidt number and terms of order λ4. Second, differentiating (31) with respect to time and using (30) and (32) we obtain the evolution equation

[

] ( )

∂〈c〉 ∂3〈c〉 ∂〈c〉 ∂ ) -U〈c〉 + De +O ∂t ∂z ∂z ∂z3

(34)

where

U ) f λ〈v〉 and De ) AK

2

+ Bλ2〈v〉2

(35)

are effective downstream advection speed and effective diffusion coefficients, respectively, where

f ) 6(1 - 2λ + 30λ2 - 276λ3) + 1 1 (19 - 129λ + 1560λ2 - 24306λ3) + O 2 , λ4 35Scλ Sc λ

(

A ) 1 + 360λ4 +

(

)

)

1 456λ3 + O 2, λ5 7Sc Sc

Inverting (31) we obtain

[

λ (38 - 387λ + C ) 〈c〉 1 - 12λ + 36λ - 1152λ 35Sc 2

3

4

]

∂〈c〉 8976λ ) + 6λ〈v〉 (1 - 4λ + 210λ2 - 2424λ3) + ∂z 2

B ) 72(1 - 10λ + 604λ2) + 1 12 (38 - 835λ + 40224λ2) + O 2 , λ3 35Scλ Sc λ

(

)

∂3〈c〉 λ2 λ4 5 , , , λ (32) ∂z3 Sc2 Sc

One great advantage in using 〈c〉 to parametrize the model is that model evolution equation (34) is then necessarily in conservative form and represents the flow of a material with a density 〈c〉 per unit channel length and downstream flux U〈c〉 - De ∂〈c〉/∂z. Since 〈υ〉 decreases downstream so does the effective diffusivity De. As a result, zone-broadening effects (see ref 20 and references therein) are expected to be weaker in the asymmetrical FFFF channel and thus the restriction on the maximum width of the inital sample in asymmetrical FFFF channels is substantially relaxed in comparison with symmetrical flow channels.

Substitute this into (29) to find that the concentration field is

(20) Schure, M. R.; Barman, B. N.; Giddings, J. C. Anal. Chem. 1989, 61, 27352743.

〈υ〉 ∂〈c〉 (19 - 258λ + 9240λ2 - 207444λ3) 35Sc ∂z 24λ2

∂2〈c〉 [3〈v〉2(3 - 38λ + 3025λ2) + K 2(1 - 5λ + 2 ∂z

522λ2)] -

2 2λ ∂ 〈c〉 [114〈v〉2(6 - 167λ + 10485λ2) + 35Sc ∂z2

K 2(38 - 645λ + 43452λ2)] + O

4336

(

)


2.4. Retention Time and Sample Width. Using the advection-diffusion equation (34), we now derive approximate expressions for the retention time tr and the retention ratio and investigate zone-broadening effects in a rectangular asymmetrical FFFF channel. These are the major characteristics of the apparatus determining its performance. First integrate (34) with respect to z to see that the total mass ∞ is conserved, dm0/dt ) 0, where m0 ) ∫-∞ 〈c〉 dz is the total amount of the sample substance. The non-dimensional coordinate of the center of sample pulse is

zc )

∫

1 m0

∞

-∞

z〈c〉 dz

|

tsr v ) v s 0 0

ln(1 - M ) M

≈ -(1 - 6λ2 + 24λ3)

(37)

tr

( )

λ(f - 2Bλ )

ln

(41)

|

(

)

≈ (1 - 2λ + 8λ2 - 80λ3)

( ) 〈v*〉st

〈v*〉out

(39)

where 〈v*〉out and 〈v*〉st denote the mean longitudinal fluid velocities at the exit from the channel and at the initial location of the center of the sample (〈v*〉st ) v0 if zst ) 0). Recalling that M ∝ u0, from (39) we deduce that the sample retention time increases with the strength of the cross-flow for two reasons. First, stronger crossflow sweeps the sample closer to the accumulation wall where the effective advection is slower (the retention ratio decreases in inverse proportion to the value of the cross-flow velocity). Second, intensifying the cross-flow while keeping the inlet flux fixed decreases the average flow speed along the channel. This leads to even smaller effective sample advection speed and, consequently, larger retention time. Note also that the values of coefficients f, A, and B increase when the coupling of the primary and cross-flow is taken into account (when the Schmidt number is finite). Thus, flow coupling decreases slightly the retention time (39), as we already indicated in section 1, and stronger zonebroadening effects associated with the increase of the effective diffusivity (35). Using the expression for the retention time tsr in a symmetrical FFFF channel (eq 38 in ref 15) and assuming that

∆D D

(42)

confirming the experimental fact that the resolution of a channel increases with increasing cross-flow. A similar expression is obtained for a symmetrical FFFF channel,15 namely

∆t sr

)

2

)

t0s

i.e., the retention time in an asymmetrical FFFF channel is slightly shorter then in a symmetrical one. If two species with close diffusivities D1 and D2 such that 2|D1 - D2|/(D1 + D2) ) ∆D/D , 1 are separated in an asymmetrical FFFF channel, then the relative difference in retention time associated with different diffusivities of species is

t sr w/u0

≈ 1 - 6λ2 + 24λ3

(38)

where z/st is the initial location of the center of the sample. Note that the limiting value of zc ) 1/λ corresponds to the point of zero mean flow velocity in a sufficiently long channel. Solving for t at which zc ) L/ξ < 1/λ, we obtain an approximate expression for the retention time tr:

)

t0

∆tr 1 ∂tr ( f - 2Bλ2)′λ ∆D ≈ ∆D ) 1 + λ tr tr ∂D D f - 2Bλ2

/ 1 zst -λ2(f - 2Bλ2)t 1 zc ) - e λ λ ξ

λz/st/ξ

|

|

with a solution

(40)

For example, if M ) 0.7, the retention time in an asymmetrical flow channel is ∼1.7 times larger than in a symmetrical flow channel. If λ and the void time are the same then

tsr

dzc d D (z ) ) U(zc) + dt dzc e c

(

tr

(36)

Then multiplying (34) by z and integrating it over the domain after some rearrangement for the particular U(zc) and De(zc) given by (35) we obtain

1τ tr ) ln 2 2 1-M (f - 2Bλ )λ

z/st ) 0 and that λ, the average inlet flow velocity and geometry of symmetrical and asymmetrical flow channels, are the same we obtain

∆D D

≈ (1 - 2λ - 4λ2 - 8λ3)

(43)

If the operational parameters are chosen in such a way that the flow rate through the accumulation wall and the retention time are the same for symmetrical and asymmetrical configurations, i.e., λ ) λs and tr ) tsr , then for the same sample materials we obtain

∆tr ∆tsr

≈ 1 + 12λ2 - 48λ3

i.e., the distance between two peaks will be slightly larger if the asymmetrical configuration is used. In the above expressions, we have assumed an infinite Schmidt number to reduce the amount of algebraic detail. For the set of characteristic parameters listed in Table 1, this approximation does not lead to any noticeable error. The dimensional void time t0 (the time required for a hypothetic particle of carrier fluid drifting with the average flow velocity to pass the complete channel length) is obtained from a simple integral

(

)

( )

〈v*〉st 1-M w τ ln ) ln / λ u0 〈v*〉out 1 - λzst/ξ

(44)


4337

t0 ) τ

∫

L/ξ

dz

z /st/ξ 〈v〉

)-

σ2 ) ξ2

∫

∞

-∞

〈c〉(z - zc)2 dz

∫

)

m0

∞

-∞

〈c〉z2 dz m0

- zc2 )

m2 - zc2 m0 (48)

Multiplying (34) by z2, integrating over the domain, and rearranging, we obtain an evolution equation for m2/m0

() [

m2 d m2 ) 2 -λ2(f - 3Bλ2) + dt m0 m0 2

AK

Figure 3. Retention ratio as a function of retention parameter λ for asymmetrical (solid and dashed lines) and symmetrical (dash-dotted line) FFFF channels. For plotting purpose, the distortion of the retention ratio function for asymmetrical channel (dashed line) is multiplied by a factor of Sc. The retention ratio for an asymmetrical FFFF channel is always greater than that for symmetrical one and tends to a finite limit for small λ.

]

+ Bλ2 + λ(f - 4Bλ2)zc (49)

Substituting the solution of the above equation in (48) we obtain /

2

(1 - λzst/ξ) -2λ2(f - 2Bλ2)t σ2 AK 2 ) 2 e + 2 2 ξ λ (f - 3Bλ ) λ2

[

]

/

σ2st

2

(1 - λzst/ξ) -2λ2(f - 3Bλ2)t AK 2 e (50) + 2 2 2 ξ λ (f - 3Bλ ) λ2

and then the retention ratio R is given by

R≡

0

t ) λ(f - 2Bλ2) ) 6λ(1 - 2λ + 6λ2 - 36λ3) + tr

(

)

1 1 (19 - 129λ + 648λ2 - 4266λ3) + O 2, λ5 (45) 35Sc Sc

where σ2st is the variance of the initial sample centered at z/st. For practical applications, it is more convenient to rewrite the expression for the variance in terms of the distance z* from the channel inlet: 2

σ2 ) The retention ratio for a symmetrical FFFF channel is easily derived from our previous results15

Rs ) 6λ(1 - 2λ)

1 - 2λ 2λ

( ( )

(w - K zst) Aw2 + 2 f - 3Bλ K 2

w - K z* w - K z/st

)

2(1 - Bλ2/(f - 2Bλ2))

(51)

)

The minimum of this expression is

σ2min ) (47)

(for example, see refs 1, 3, 6, 14, and 20) and are plotted in Figure 3. Note the increase in retention ratio for the asymmetrical FFFF channel in comparison with a symmetrical one. Another important feature of the asymmetrical flow channel is that the distortion of the primary parabolic flow leads to further increase in retention ratio. It is nonzero even in the limit of very strong cross-flows λ f 0. In contrast to the symmetric configuration, increased crossflow not only sweeps the sample closer to the accumulation wall but also increases the longitudinal component of the primary flow in the boundary layer. As a result and in contrast to a symmetrical FFFF channel, the sample will eventually elute from an asymmetrical FFFF channel no matter how strong the cross-flow is. The variance σ2 of the sample is proportional to the square of its width and characterizes the effectiveness of the separation process. The variance is 4338

](

/ 2

σ2st -

(46)

Both of these expressions for the retention ratio coincide in the small λ, large Sc limit with the well-known approximate expression

R ) 6λ coth

[

(w - K z*) Aw2 + 2 f - 3Bλ K 2


Aw2 f - 3Bλ2

(52)

and it is reached at z* ) w/K , i.e., at the location where all carrier fluid is lost in cross-flow and the flow stops completely. The variance of the sample at the exit of the channel is 2

σL )

[

Aw2 L2 Aw2 2 2 (1 M ) + σ + st f - 3Bλ2 M 2 f - 3Bλ2

(

) ](

z/st L2 1 M L M 2

2

1-M 1 - M z/st/L

)

2(1 - Bλ2/(f - 2Bλ2))

(53)

which has a minimum given by (52) at M ) 1, i.e., when all carrier fluid is lost in the cross-flow. Thus, the more fluid is lost in an asymmetrical FFFF channel in a cross-flow the weaker the integral zone-broadening effects are. For further simplification, we assume z/st ) 0 and expand the above in λ to obtain an approximate expression for the plate height

[

]

σ2L σ2st B ≈ (1 - M )2 1 - 2λ3 ln(1 - M ) + L L R ln(1 - M ) M B 2D M 1+ λ3 1 + (1 - M )2 A Rυ0 2 R 2 M

H h ≡

[

(

)]

R′λ ∆D ∆tr ≈λ tr R D

where

trv0 ∆trvo ∆D ) ≈ λR′λ (1 - M ) 4σL 4σL D R′λ L 1 - M ∆D ln(1 - M ) (55) R 4σL M D

λ

)

Bλ2 B χ ≈ 2λ3 1 ln(1 - M ) ≈ f - 2Bλ2 f - 2Bλ2 24λ (1 - 8λ + 582λ - 12λ2 ln(1 - M ))

2

and s

The first line in (54) represents contribution H h i due to the initial width of the sample, the second one contribution H h long due to the longitudinal diffusion, and the third one, H h neq, accounts for the nonequilibrium effects.5,13 (Strictly speaking the longitudinal diffusion plate height is that which would be observed if the only acting mechanism is molecular diffusion. It is

H h long )

2Dtr 2D ln(1 - M ) )L Rv0 M

The second line in (54) has a nonequilibrium contribution, and only in the limit M f 0 becomes truly the longitudinal diffusion plate height.) As M f 0 (no fluid is lost in the cross-flow), the first two lines reduce to the well-known expressions for the plate heights in a symmetrical flow channel, H h si ) σ2st/L and H h slong ) 2D/Rυ0 (e.g., ref 2), respectively. The limiting expression for the nonequilibrium plate height formally becomes the same as the one for a symmetrical flow channel except that function χ is different from s

3

tsr

R′sλ ∆D Rs D

≈λ

respectively. Then the resolution functions for asymmetrical and symmetrical channels are

R ≈

3

∆tsr

-

w2v0 ln(1 - M ) χ (1 - M )2 (54) D M

(

and

2

χ ≈ 24λ (1 - 8λ + 12λ )

obtained for a symmetrical flow channel (which is derived from the results reported in ref 15 and coincides with eq 32 in ref 5). This is due to the difference in retention characteristics of the channels. As the portion of fluid lost in cross-flow increases (M increases) H h i and H h neq quickly decrease because of the focusing due to the fluid deceleration while H h long slightly increases due to the longer retention time so that the overall effect is to decrease the plate height. Finally, it is useful to quantify the resolution of the apparatus. The performance of the FFFF channel is assessed by considering the elution of two substances with nearly identical diffusivities and is characterized by the ratio R of the time interval ∆tr between the two concentration peaks to the width in time of an individual peak 4σL/vo estimated at the exit from the channel. Using (37) and (35), we obtain the dimensional zone drift velocity at the channel outlet. It is vo ) v0U(L/ξ) ) λf υ0(1 - M ) ) (R + 2Bλ3)v0(1 - M ). The outlet zone velocity in a symmetrical channel is vso ) vs0Rs. Equations (42) and (43) can be written as

R ≈

∆tsr vso

≈

4σsL

tsvs s r 0 λR′λ s 4σL

R′sλ L ∆D ∆D )λ s s D R 4σ D

(56)

L

respectively, where we neglect terms of order λ3. In typical operational conditions, the influences of the initial sample width and molecular diffusion are negligible in comparison with nonequilibrium effects. Then σL ≈ (LH h neq)1/2 and using (54) we obtain

R ≈

x

λ R′λ 4 R

s

λ R′λ Dt0 ∆D and R s ≈ 2 D 4 Rs χw

x

Dt0s ∆D χsw2 D

(57)

Finally, to compare the resolutions of symmetrical and asymmetrical channels, we consider their ratio for the cases of equal void or retention times assuming that geometry of channels, λ, D, and ∆D are the same

R Rs

|

R Rs

|

≈

tr ) tsr

t0 ) t0s

≈

R′λ R′sλ

x

χsRs ≈ 1 - 258λ2 + 6λ2 ln(1 - M ) (58) χR

x

R′λRs R′sλR

χs ≈ 1 - 261λ2 + 6λ2 ln(1 - M ) (59) χ

Thus, we conclude that the resolution of an asymmetrical FFFF channel is slightly lower than that of a symmetrical one under similar flow conditions. 2.5. Numerical Experiments. To further support the advection-diffusion model (34), we solve it numerically and compare the result with experimental data obtained by Wahlund and Giddings (ref 14, Figure 10). They considered the separation of cytochrome c (DC ≈ 1.2 × 10-6 cm2/s), albumin (DA ≈ 6 × 10-7 cm2/s), and thyroglobulin (DTh ≈ 2.3 × 10-7 cm2/s) in a channel of length L ) 43 cm, width b ) 2.07 cm, and thickness w ) 0.052 cm. Using provided information on flow rates, we deduce that u0 ≈ 6.07 × 10-4 cm/s and v0 ≈ 0.763 cm/s. Initially the sample consisting of the mixture of the above components was focused at z/st ) 4.0 cm. No other information on the initial spreading of the sample or actual concentrations of the components was provided. Thus, in our numerical simulation, we choose the initial distributions to be of the form aj exp[-(z* - z/st)2/(2σ2st)], j ) 1, 2, 3, where values σ2st ) 1/2 cm2 and a1:a2:a3 ) 1.00:0.52:0.16 Analytical Chemistry, Vol. 72, No. 18, September 15, 2000

4339

experimental pulse for thyroglobulin is wider than predicted by the model. Such a broadening of the concentration pulse would occur if, for example, the thyroglobulin component consisted of different species with diffusivities occupying a finite range centered at DTh ) 2.3 × 10-7. Since the experimental thyroglobulin pulse is wider by about ∆trTh ≈ 2 min than the numerical solution for monodisperse substance, we estimate the diffusivity range using (42) as

∆DTh ≈

Figure 4. Comparison of numerical solution of model (34) (thin line) with the experimental results by Wahlund and Giddings14 (thick line). Flow parameters are discussed in the text.

were chosen for cytochrome c, albumin, and thyroglobulin, respectively, to obtain numerical peak concentrations similar to the ones observed experimentally. The comparison is shown in Figure 4. A small hump between the albumin and thyroglobulin peaks on an experimental curve suggests that a fourth component was present in the mixture, but it is not discussed in ref 14, and consequently, we do not model it. The theoretical void time (44) is shown by the thin vertical line in Figure 4. It is slightly longer then the experimental value shown by a vertical arrow. Apart from the errors associated with scanning and manual digitizing of Figure 10 of ref 14, there is a physical reason for this discrepancy. The theoretical value of the void time is based on the average flow speed while the experimental void time is that of the first detector response of the presence of the eluted species in the flow. It occurs when a small fraction of underrelaxed sample material advected close to the horizontal midplane of the channel with the speed larger than the average longitudinal velocity reaches the location of a detector. The incomplete relaxation of the initial sample is also responsible for the “fronting”14 of the experimental curves: the front slopes of experimental concentration profiles are less steep than the tail ones, meaning that a fraction of the sample material advects ahead of the main sample. The second possible reason for the fronting might be associated with an asymmetrical distribution of sample material in the initial sample. Fronting is weaker for highly retained species such as thyroglobulin. In fact, the front concentration slope for thyroglobulin is steeper than the tail one. This situation qualitatively agrees with the predictions of model (34): concentration pulses obtained from numerical simulation are steeper in time at the front than at the back. When the mixture of different species is separated in asymmetrical FFFF channel, this front steepening is stronger for highly retained components since the influence of fronting is less profound than for species with higher diffusivities, which pass the channel faster. This is confirmed by the experimental results: cytochrome c peak has substantial fronting, albumin peak is almost symmetric and thyroglobulin peak has more tailing than fronting. Since model (34) does not take into account an underrelaxed sample fraction and associated fronting effects, all concentration profiles obtained numerically are qualitatively similar to the experimental distribution of thyroglobulin. The model solution gives an excellent prediction of the retention times for all species (locations of experimental and numerical concentration maxima coincide). The model solution also provides very good estimation of the width in time of the pulses for cytochrome c and albumin. The 4340


∆trTh D ≈ 1.4 × 10-8 cm2/s trTh Th

(60)

or DTh ) (2.30 ( 0.07) × 10-7 cm2/s. The deviation of the diffusion coefficient, which causes widening of the experimental pulse, is only ∼3%. The reason for widening the range of the diffusion coefficient for thyroglobulin might be related to the fact that the highly retained sample is swept too close to the accumulation wall so that its local concentration becomes too high and so-called sample-overloading effects21 occur. For example, the sample overloading may result in nonlinear deviations from Fick’s first law of diffusion. They are not accounted for by the basic eq (19) and, subsequently, by model (1). Thus, the application of the above model is limited to regimes and substances for which the hypothesis of linearity of diffusion is valid. Having established the reliability limits of model (34), by comparing its predictions with the experimental results, we use this model to provide detailed information on the dispersion of the sample within the asymmetrical FFFF channel of this experiment. In Figure 5, we present predicted characteristics of the sample components as they advect along the channel. Since dispersion and focusing due to the deceleration of the mean flow oppose each other, the maximum concentration of the samples presented in Figure 5a is a nonmonotonic function of location along the channel. It is inversely proportional to σ whose variation (51) along the channel is shown in Figure 5b. The sign of dσ2(z*st)/dz depends on the initial width of the sample. It follows from (51) that if

σ2st
K. If the injection of sample is made quickly in the vicinity of the stagnation zone, then the sample substance essentially displaces the carrier fluid from the injection region without being swept back by the counterflow. This suggests a rough estimation of the initial width of the sample. For example, if the amount m0 of injected sample is a fraction of a milliliter in a typical FFFF channel, then the width of the zone and the concentration of the contaminant per unit channel length, respectively, are

wi0 ≈

m0 wb

and

〈c〉i0 ≈

m0 ) wb wi0

(63)

where b is the width of the channel. Litzen’s1 experiment had b ) 1.25 cm and m0 ) 0.05 mL, which would produce the initial sample with a width of order wi0 ≈ 3 cm. After injection of the sample is completed, two processes take place. First, cross-channel relaxation leads to redistribution of Analytical Chemistry, Vol. 72, No. 18, September 15, 2000

4341

sample so that its concentration becomes essentially exponential across the channel thickness as given by (33). The time necessary for this relaxation was estimated using cross-channel advection arguments by Wahlund and Giddings14 to be

trel )

( (

w 2 1 + 2x*/w x*/w ln + u0 9 1 - x*/w 3(1 - x*/w)

)

)

(64)

Using this formula and experimental data from ref 1, we find first that during the 30-s relaxation period after injection was completed ∼90% (x*/w ≈ 0.89) of the sample substance was swept toward the accumulation wall, and consequently, the exponential distribution of the sample material was mostly established. Second, the sample is advected toward the stagnation zone by a longitudinal counterflow. If the time passed after injection is long enough, then the advection-diffusion model (34) predicts that the cross-channel concentration evolves to the equilibrium distribution near the stagnation zone approximately obtained by solving

[

∂〈c〉i ∂ ) -Ui〈c〉i + K ∂t ∂zi

2 i

]

∂〈c〉i ∂zi

(65)

the other terms in (34) being small near the focusing point. Thus, expect that after the focusing stage is completed the steady concentration is given by (33) where

〈c〉i ≈ m0 w

x

m0 w

x

( ) fiz/2 i 2

fi exp 2π 2w

(

)

(

x

)

/ 2

(z* - zst) 3 exp -3 π w2

(66) (In practice, injection of the sample is done through a thin pipe, and consequently, its initial distribution varies significantly in the spanwise direction. Nevertheless, spanwise variations may disappear quickly because of the action of the counterflow so that by the end of the focusing stage the distribution of sample material is essentially two-dimensional as shown by Wahlund and Giddings14 in their Figure 9. Note that the steady solution of equation (34) is actually

(

)

2

-fi/(2Biλi ) z/2 i Biλi2 2

w

(67)

due to the spatial variation in the shear dispersion contribution (Biλ2i 〈v〉2i ). For large cross-flow rates solutions, (66) and (67) are very close. In particular, the minimum width of the focused sample as estimated from both expressions essentially does not depend on the flow conditions and is wi ) 2w/fi1/2 ≈ 2w/61/2, which confirms the estimation given in ref 14 for highly retained species. Thus to avoid excessive algebraic complexity we use the simplified eq 65. Expression 66 corresponds to a very sharply peaked concentration with 〈c〉i max/〈c〉i0 ) m0(fi/(2π))1/2/(w2b) ≈ 250 for the experimental parameters cited in ref 1. This concentration 4342

(

∞

〈ci〉(z, t) )

∑C

n

exp λnt -

n)0

) (x )

fi z/2 i 2Ai w2

Hn

fi z/i

Ai w

(68)

where Hn denotes the Hermite polynomials and constants Cn are determined by the distribution of sample material immediately upon injection.22 Eigenvalues λn ) -nfiλ2i , n ) 1, 2, ..., characterize the speed of the focusing process. The slowest decaying component of the solution corresponds to n ) 1. The characteristic time for decay to the focused state is then

tf ∝

τi λ2i fi

≈

w2 6D

(69)

)

fi(z* - z/i )2 m0 fi ≈ exp 2 2π w 2w

〈c〉i ∝ 1 +

might be too high and lead to overloading of the channel.21 Indeed this was noted in ref 1 when low injection flow rates were used and the initial spread of the sample was sufficiently small. In practice, the focusing up to the limiting distribution (66) is never performed (ref 1, p 469) to avoid nonlinear effects associated with the interaction of the sample molecules when their concentration becomes too high.21 Moreover, analysis of (53) shows that there is no need to focus the initial sample too sharply during the counterflow. All one needs to do is to focus initial sample to a size significantly smaller than the zone broadening that occurs during the elution. To estimate the time required for such focusing, we note that the advection-diffusion equation (65) with an arbitrary localized initial condition 〈ci〉(z, 0) ) 〈ci0〉(z) has the solution


Note that when the cross-flow rate is sufficiently large (compare with ref 14), the focusing time and the limiting width of the focused sample practically do not depend on the flow conditions and are determined only by the thickness of the channel and the diffusion coefficient of the sample substance (assuming it is not affected by the sample overloading effects). In the experiments by Litzen,1 ferritin (D ) 4.1 × 10-7cm2/s) and the cow pea mosaic virus (CPMV, D ) 2 × 10-7cm2/s) were used as experimental sample components. The estimated focusing times for these samples are tf ≈ 70 s and tf ≈ 140 s, respectively. The experimental focusing time was 30 s. Apparently, this was not quite long enough to establish a sharply peaked distribution of sample material along the channel. Indeed, the observed sample width was ∼1 cm rather than a fraction of a millimeter as predicted by (66) and (67). Previous estimates of the focusing time are based purely on advection arguments (see, for example, ref 14, eq 34, and ref 1, eq 24). It is rewritten in our notation as

tf ≈

( )

z/u - z/st w2 ln / 6D zf - z/st

(70)

where z/u and z/f denote the location of the boundaries of unfocused and focused samples, respectively. Substituting z/st ) 0 cm, z/u ) 1.5 cm, and z/f ) 1 cm, we obtain tf ≈ 28 s for D ) 4.1 × 10-7 cm2/s, which confirms the experimental observations by (22) Suslov, S. A.; Roberts, A. J. J. Aust. Math. Soc. B, Part E (electronic) 1998, 40, E1-E26; http://jamsb.austms.org.au/V40/E007

Litzen.1 Thus, we conclude that in practice the cross-channel distribution of sample material is typically close to its relaxed state given by (66). On the other hand, its longitudinal spread is far from the equilibrium and depends strongly on the injection procedure so that its width has to be determined experimentally. 3.2. Secondary Relaxation and Initial Conditions for the Model. In general, the cross-flow rates at the focusing/relaxation stage and elution stage are different; that is, λ * λi. If, for example, the cross-flow velocity at the focusing stage is smaller than that in the advection regime (λi > λ),1 then the sample will tend to be swept more toward the accumulation wall as it advects along the channel during the elution before it relaxes to new flow conditions. We call this process secondary relaxation. Although this is a relatively quick process,14,15 it has a long-term effect on the distribution of the sample substance; in particular, it will affect the width of the sample at the outlet. In other words, using the concentration distribution (66) obtained during the focusing directly as an initial condition for the sample evolution, eq (34) would lead to the deviation of the long-term predictions of the model from the experimental results. The correct initial condition for 〈c〉 has to be found by projecting the above distribution on the model appropriate to the elution flow. We do this projection by following the procedure outlined by Roberts.23,24 The proper model initial condition is the projection 〈c〉0 from the physical field corresponding to 〈c〉i to the center manifold (long-time solution) along the isochronsin the state space an isochron is a surface of all the initial states that have the same long-term dynamics on the center manifold (up to an exponentially small error). Consequently, under the inner product defined as

〈a, b〉 ≡

∫ ∫ L/ξ

0

1/λ

0

ab dx dz

(71)

The dual operator D appearing in (73) is obtained from equations adjoint to (19) with respect to the inner product defined as

〈〈a, b〉〉 ≡

DZ ) -

∂Z ∂Z ∂Z ∂2Z +v +u + 2 +K ∂t ∂z ∂x ∂x

2

∂2Z ∂z2

(76)

with boundary conditions

∂Z ) 0 at x ) 0 and x ) 1/λ ∂x

(77)

and normalization condition

〈Z, E 〉 ) 1

(78)

The solution of (73) is quite involved, but since the sample is swept close to the accumulation wall during the focusing in the counterflow regime we only find a leading approximation of this solution in the vicinity of the accumulation wall. The iterative algorithm used to find this solution is largely the same as described in section 2.2. To the lowest order Z ) δ(z - ζ), where δ(z - ζ) is the Dirac delta function and 0 < ζ < L/ξ0 is an arbitrary inner point in the channel. Using this as an initial guess we find the two-term approximation for the projection function

Z ≈ δ(z - ζ) - δ′(z - ζ)λ〈v〉[6 - 3x2 - 12ex - 1/λ -

∂Z ∂Z ≈ 6λ(1 - 2λ)〈v〉 ∂t ∂z

(80)

(72)

where Z is the projection function satisfying the dual equation

DZ ) 〈DZ, E 〉Z

(73)

Here we neglected terms of order exp(-1/λ) and λ3 and higher. Now that the expression for the projection function is found we solve (72) for 〈c〉0. Recalling that

∫

L/ξ

0

g(z)δ(z - ζ) dz ) g(ζ) and

∫

and E is the local tangent vector to the center manifold

(

2

)]

∂ ∂ λ + O 2, , λ3 ∂z ∂z Sc

g(z)δ′(z - ζ) dz ) -g′(ζ) (81)

we obtain

(74)

Note that E is an operator rather than just a function as occurs in the finite dimensional cases discussed in ref 24. Being introduced into the inner product (71), it acts on the other factor involved before the integration is performed. (23) Roberts, A. J. J. Aust. Math. Soc. B 1989, 31, 48-75. (24) Roberts, A. J. Comput. Phys. Commun., in press.

L/ξ

0

∂c ) e-x 1 - λ2(12 - 3x2 - x3) + λ〈v〉(6 - 3x2 ∂〈c〉 2λ(12 - 3x2 - x3))

(75)

2λ(12 - 3x2 - x3)] (79)

〈Z, 〈c〉i - 〈c〉0〉 ) 0

[

t

0

so that here

the model initial conditions are determined to satisfy

E≡

∫ 〈a, b〉 dτ

〈c〉0 ≈ 〈c〉i + 6

d(〈c〉i〈v〉) × dz λi 1 - (λ + λi - 4λ2 - 4λλi - 2λ2i ) (82) λ

( )

Assuming that

〈c〉i )

m0 e-(z* - z/st)2/(2σ2st) σst x2π


(83) 4343

Figure 7. Velocity field and an instantaneous concentration field near the accumulation wall: (a) during the focusing and (b) soon after elution started. Plot c represents cross-channel average concentration during the focusing (solid line) and after elution has started (dashed line). λi ) 0.02, λ ) 0.01, σ ) 1; other parameters correspond to the experiment by Litzen.1

we obtain

[

〈c〉0 ≈ 〈c〉i 1 + 6(λi - λ)(λ + λi - 4λ2 - 4λλi - 2λ2i ) ×

(

1+

z* - z/st w - z* K σ2st

(

)]

)

(84)

We interpret this result from a physical point of view. Assume, for example, that the cross-flow at the focusing stage is weaker than the one during the elution. Indeed, during the focusing in the experiment by Litzen1 the volumetric flow rate through the accumulation wall of area of 34.63 cm2 was 3 mL/min. Using this information, we deduce that for CPMV as an experimental substance v0i ≈ 1.45 × 10-3cm/s and λi ≈ 0.02. On other hand, the elution stage was characterized by λ as small as 0.003 (see Table 1). Thus, as soon as the switching between the focusing and elution regimes occurs the sample is swept much closer to 4344


the accumulation wall and relaxes to the new flow conditions. This is a relatively quick process and thus it is not described by model (33)-(34). Nevertheless, the obtained projection is able to provide quantitative characteristics of this transient process without considering it in detail. In Figure 7a, the concentration field (83) is shown just before the elution is started. The distribution of the sample material is approximately symmetric about the stagnation point, exponential across the channel, and Gaussian in the longitudinal direction (solid line in Figure 7c). Note that due to the cross-flow the sample is concentrated in a thin layer near the accumulation wall with the thickness of ∼0.05w. Figure 7b shows the distribution (84) of the sample after relaxation to the elution conditions. First, observe that the sample material is indeed concentrated closer to the accumulation wall because of the stronger cross-flow. Second, while the sample relaxes to new conditions, it is advected toward the outlet. Thus, by the time cross-flow relaxation is finished, the sample is shifted slightly in the direction of primary flow. This is clearly seen from Figure 7b and c (dashed line).

It is instructive to find the variance σ20 of the concentration distribution 〈c〉0 soon after the start of elution

2

σ0 )

∫

∞

-∞

〈c〉0(z* - z/st)2 dz*

∫

∞

the primary flow, thus providing more favorable conditions for further focusing. Since the sample width decreases during the secondary relaxation, the peak concentration of sample slightly increases as seen in Figure 7c.

)

〈c〉0 dz* -∞

σ2st(1 - 12(λi - λ)(λ + λi - 4λ2 - 4λλi - 2λi2)) (85)

We conclude that switching to higher cross-flow rate (λ < λi) leads to a slight longitudinal focusing of the sample. This is because the increase of the cross-flow rate causes faster deceleration of

ACKNOWLEDGMENT This research is supported by a grant from the Australian Research Council. The authors thank the reviewers for valuable comments. Received for review December 29, 1999. Accepted June 16, 2000. AC9915022


4345

Modeling of Sample Dynamics in Rectangular Asymmetrical Flow

Recommend Documents