Table Of Content

Taylor dispersion with absorbing boundaries: A Stochastic Approach Rudro R. Biswas∗ Department of Physics, Harvard University, Cambridge, MA 02138, USA† 7 Pabitra N. Sen 0 Schlumberger-Doll Research, Ridgefield, CT 06877-4108, USA 0 2 We describe how to solve the problem of Taylor dispersion in the presence of absorbing boundaries using an exact stochastic formulation. In addition to providing a clear stochastic picture of n Taylordispersion,ourmethodleadstoclosed-formexpressionsforallthemomentsoftheconvective a displacement of the dispersing particles in terms of the transverse diffusion eigenmodes. We also J find that the cumulants grow asymptotically linearly with time, ensuring a Gaussian distribution 6 in the long-time limit. As a demonstration of the technique, the first two longitudinal cumulants 1 (yielding respectively the effective velocity and the Taylor diffusion constant) as well as the skewness (a measure of the deviation from normality) are calculated for fluid flow in the parallel plate ] h geometry. We find that the effective velocity and the skewness (which is negative in this case) are c enhanced while Taylor dispersion is suppressed dueto absorption at theboundary. e m PACSnumbers: 47.27.eb,05.40.-a,47.55.dr,47.57.eb - t a In a laminarly flowing fluid with flow velocity varying even in their scheme, obtaining higher than second mo- t s perpendicular to the flow direction, suspended particles mentsandthusestimatingdeviation from Gaussianityof . t will jump randomly acrossfast andslow streamlines due the particle concentration profile along the direction of a m to transversediffusion. This fluctuating convectionalve- flow, becomes prohibitively tedious and hence was not locityarisingfromarandomsamplingofvariousstream- addressed. - d lines disperses these particles along the direction of flow Ourapproachgoesastepfurtherthanthiscornerstone n in a well-studied process known as Taylor dispersion[1– work. We provide the correct stochastic picture of the o 4]. The particles drift along the flow with an average underlying process. When decay is present, the usual c velocity v and disperse linearly with time t along the probabilistic techniques[12, 13] fail. We formulate, for e [ flow, i. e. the r.m.s displacement in that direction is the first time, the correct technique — ultimately pro- 1 (D +D)twheretheconvection-inducedTaylordis- viding a visually appealing and computationally simple taylor v persion coefficient D is found to scale as ℓ2v2/D, D stochastic method to find all moments/cumulants of the 2 pbeing the molecular tdaiyffloursion coefficient and ℓ beeing the longitudinal(i.e.alongtheconvectiveflowdirection)dis- 9 transverse dimension. Taylor dispersion has turned out placementforuniform(linear)laminarflow. Wefindthat 3 to be important both fundamentally in hydrodynamics thecumulantsgrowasymptoticallylinearlywithtimeen- 1 and in its applications to diverse fields ranging from bi- suringaGaussiandistributioninthelong-timelimit,asis 0 7 ological perfusion, chemical reactors and lab-on-chips to characteristicofdriftingrandomwalkersonalattice[18]. 0 soil remediation and oil recovery[4–11]. Consider diffusing particles being carried along by a / fluid flowing unidirectionally in the x-direction through t Inspiteofhavingbeenaroundformanydecades,Tay- a a uniform cross-section A. The particle position in the m lor dispersion in the presence of absorption,evenin sim- transversesectionAwillbedenotedby~y,whilethecom- - ple geometries (such as Poiseuille flow in tubes), has plete position will be specified by x (x,~y). The con- d proven to be a mathematically challenging problem[12– centration of the particles, N(x,t), ev≡olves according to n 16]. As Taylor dispersion under Poiseuille flow in tubes the convective diffusion equation inside the fluid volume o is now routinely used to measure molecular diffusion c V: coefficients[17],itisimportanttocomputetheerrorsdue Xiv: taoppalbicsaotripotniosnoauttltihneedwainlls.thTeheprseavmioeuhsopldasratgruraepfhor. otThheer ∂N∂(xt,t) =D∇2xN(x,t)−v(x).∇N(x,t) x∈V (1) methods used originally[1–3] workedonly in the absence r where D is the molecular diffusion constant of the par- a ofabsorption. Thiswaspointedoutinthemost‘modern’ ticles in the still fluid and v(x) v(~y)xˆ is the x- treatmentavailable—byGill,LunguandMoffatt[15,16] ≡ independent convective velocity field of the fluid. The dating back to over 40 years ago. They showed that the presence of absorbing walls is taken into account by the effective velocity is enhanced while Taylor dispersion is boundary condition suppresseddue to absorptionatthe boundary. However, Deˆ.∇xN(x,t)+ρN(x,t)=0, x ∂V (2) ∈ where ρ is a measure of the surface absorption rate that ∗Electronicaddress: [email protected] varies from zero (reflecting walls) to infinity (perfect ab- †AlsoatSchlumberger-DollResearch,Ridgefield,CT06877, USA sorption), while eˆ(x) is the normal to the surface ∂V at 2 x. Since the cross-section is uniform, eˆ(x).xˆ = 0. This formula[19] enables us to replace eˆ.∇x by eˆ.∇y in (2). By integrating over x in (1) and (2), a diffusion equa- n(~y,t)= dy′G(~y,t t′ ~y′)n(~y′,t′). (5) tion for n(~y,t) = dxN(x,t) – the concentration of − | ZA particles in the transverse section A – is obtained with absorbing boundaryRconditions (when N,∂ N 0 as for t t′. In other words, it gives us the fraction of par- x ): x → ticles≥starting out from an arbitrary initial state (~y′,t′) →±∞ and surviving till a specific final state (~y,t). ∂n(~y,t) 2 To illustrate the technique of obtaining the moments =D n(~y,t) ~y A (3a) ∂t ∇y ∈ of x(t), we start with the first one – the mean. It is ob- D eˆ.∇ n(~y,t)= ρn(~y,t) ~y ∂A (3b) tained by averaging (4) over the particles that survive y − ∈ the absorbing walls till time t. Following (5), the parti- Thisyieldsthe~y-statisticsofthedispersingparticlesthat cle density at (~y′,t′) is given by G(~y,t′ ~y )n(~y ,0)dy . i i i we will use in the formalism below involving the equa- However, in the remaining time (t t′) o|nly a fraction tions (4) and (5). G(~y ,t t′ ~y)ofthesewillsurviveRtil−lafinalstate(~y ,t). f f − | Wenowproceedtofindthemomentsofx(t). We shall Thus, the number density of particles that pass through not show explicitly the noise that is responsible for the (~y′,t′) and also survive till t is usual molecular diffusion along the x direction and con- centrateonthedominantconvectivenoise[21]encodedin ν(~y,t′ t)= dy dy G(~y ,t t′ ~y)G(~y,t′ ~y )n(~y ,0) f i f i i the stochastic equation that gives rise to the convective | − | | ZA ZA part in (1). This equation is (6) Multiplying this density by the corresponding displace- dx(t) t mentv(~y)dt′,averagingoverallpossibleintermediatepo- =v(~y(t)) x(t)= dt′v(~y(t′)) (4) dt ⇒ 0 sitions and finally adding up all these averageddisplace- Z ments, we find that Intrinsic diffusive noise along the x-direction may be in- corporated easily — in particular, the mean will be un- 1 t x(t) = dt′ dy v(~y)ν(~y,t′ t) (7) affected while Dtaylor will be additively augmentedbyD. h i N(t)Z0 ZA | We are considering the case x(0) = 0 here. To proceed wherethe(normalizing)denominatorN(t)isequaltothe further, we introduce the Green’s function/propagator total number of particles that survive till t. Similarly, G(~y,t~y′) – the t > 0 solution to (3) when n(~y,0) = squaring (4) and averaging over the surviving particles δ(~y |~y′). It’s most important use is encoded in the one can show that the second moment is given by − 2 2 (x(t)) = dt1dt2 dy1dy2v(~y1)v(~y2)ν(~y1,t1 ~y2,t2 t) (8a) h i N(t)Z0≤t1≤t2≤t ZA⊗A | | ν(~y1,t1 ~y2,t2 t)= dyidyfG(~yf,t t2 ~y2)G(~y2,t2 t1 ~y1)G(~y1,t1 ~yi)n(~yi,0) (8b) | | − | − | | ZA⊗A The density of walkers who survive till t after passing eigenmodes[19]: through (~y1,t1) and (~y2,t2) is given by ν(~y1,t1 ~y2,t2 t), | | analogous to (6). All higher moments may be obtained ∞ by repeating the argumentsleading to the above results. n(~y,t)= n ψ (~y)e−λkt (9) k k k=0 X Thepreceedingargumentsthusprovidethecorrectex- position of the stochastic technique to be used when ab- where n0 =0 (necessarily) and 6 sorptionis alsopresent. This will be the mainstay ofthe remaining results to be derived in this paper. We note 2 (D +λ )ψ (~y)=0 ~y A (10a) that the formulation is exact in the framework we con- ∇y k k ∈ sider — viz. the statistics obeyed by ~y is independent of (Deˆ.∇y+ρ)ψk(~y)=0 ~y ∂A (10b) ∈ x. (λ0 <λ1 <λ2...) dy ψ (~y)ψ (~y)=δ (10c) j k jk We can solve (3) by expanding in terms of diffusion ZA 3 The propagatoris given by [19] Usingthesamearguments,highercumulantsmaybesim- plifiedtoformsanalogoustoequations(14)or(15)inthe ∞ long time limit. Because of the straightforward algebra G(~y,t~y′)= ψ (~y)ψ (~y′)e−λnt (11) | n n involved,computeralgebrasystems(weusedMathemat- nX=0 ica) may be programmed to do this correctly and elimi- The lowest eigenvalue λ0 for this problem becomes non- nate corrections O(e−t/τ). zero for ρ=0 — this is the main difference between the Wehaveanalyticallycomputeduptilthefifthcumulant absorbing6andnon-absorbingcases. Atlongtimesallthe thiswayandfoundthattheyareallproportionaltotime higher modes decay with inverse mean lifetimes of in the long-time limit. We assert that this holds to all orders. This implies that the p.d.f of the longitudinal λ n2/τ (12) displacement tends to a Gaussian with deviations that n ∼ die as 1/t — similar to the case of reflecting boundary where τ is the characteristic relaxation time-scale of the conditions[3]. A similar behaviour is also exhibited by transversediffusionprocessandscalesasℓ2/D whereℓis drifting random walkers on a lattice[18]. the characteristic length-scale of the cross-section A. In thelong-timelimitt τ wecanthusworkwithonlythe lowest mode in (9) an≫d the errors will go as (e−t/τ). y=+ℓ(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0) We consider this long time limit in what fOollows. We - find that the moments can be conveniently evaluated in terms of simple ‘matrix elements’ y=0 fluid -vx(y)=v0`1−y2/ℓ2´ - v = dy ψ (~y)v(~y)ψ (~y) (13) y=−ℓ jk j k @@@@@@@@@@@@@ ZA Using(11)and(10)in(7)andkeepingonlythedominant termsasoutlinedinthepreviousparagraphweobtainthe FIG. 1: Taylor dispersion between parallel plates followingexpressionforthefirstmomentinthelongtime limit: Finallyletuscompute aspecificexample–Taylordis- persion of particles suspended in fluid flowing unidirec- ve x(t) t≫2τ v00 t+ (e−t/τ) (14) tionally between two infinite parallel plates (see fig. 1). h i−−−→ O Absorption at the walls is characterised by the dimen- This resultis verycounterz-i}n|t{uitivesince one na¨ıvelyex- sionless parameter α = ρℓ/D (ρ is defined in (2)). For pects, at late enough times, the effective velocity v to convenience, we have used the dimensionless coordinate e be the same as the asymptotic instantaneous mean ve- η =y/ℓbelow. Theeigenfunctions(10)forthiscaseturn locity that is proportional to vn vψ0 — since out to be: A ∝ A n(~y,t → ∞) ∼ ψ0(~y). HoweveRr, absorptRion at the walls α2+φ2 nπ modifies the probability that a particle starting out at ψ (η)= n cos φ η+ (18a) ~y survives after a long time (i. e. τ) from being ~y- n s(α+α2+φ2n) (cid:18) n 2 (cid:19) ≫ independent to being proportional to ψ0(~y) — this pro- α nπ vides the additional power of ψ0 in (14) (recall that φ =tan φn+ 2 , n=0,1,... (18b) v00 = vψ02). n (cid:18) (cid:19) Analogously, this time using (8) instead of (7) above We obtained, for α 1 R ≪ we calculate the variance in x to be α 11α2 φ =√α 1 + + (α3) , n=0 (δx(t))2 t≫3τ t 2 ∞ vk20 + (e−t/τ) (15) n (cid:18) − 6 360 O (cid:19) h i−−−→ × k=1λk−λ0 O (19a) X 2 3 nπ 2α 2 2α 2α Dtaylor = + + , n=1... 2 nπ − nπ nπ O nπ | {z } (cid:18) (cid:19) (cid:18) (cid:19) while the third cumulant simplifies to: (19b) ∞ ∞ For α , 3 t≫4τ v0j(vjk δjkv00)vk0 →∞ (δx(t)) 6t − (16) h i−−−→ j=1k=1 (λj −λ0)(λk−λ0) φ = (n+1)π, n=0,1,... (20) XX n 2 showing that the skewness, which is a measure of the Finally, the previously introduced eigenvalues λ (see deviation from Gaussianity, dies out with time: n (10)) are obtained from the relation: (δx(t))3 τ φ2 γ = h(hδx(t))2i3i/2 ∼rt (17) λn =Dℓ2n (21) 4 Using our formulæ (14) and (15) we have calculated showing that it is negative and enhanced due to absorp- the effective velocity and the Taylor dispersion constant tion at the walls. The zero-absorptioncase of this result to be: (23) has been independently verified using the method presented in [3]. 2v0 2 2 v = 1+ α+ (α ) ,α 1 e 3 15 O ≪ (cid:18) (cid:19) Inconclusion,wehaveapicturesqueandcomputation- 2v0 3 ally simple solutionto a problemof generalinterest. For = 1+ , α (22a) 3 π2 →∞ example, we may now estimate the error in measuring (cid:18) (cid:19) 8v2ℓ2 4 D using Taylor dispersion[17]. The value of λ0 and con- 0 2 D = 1 α+ (α ) ,α 1 sequently α may be obtained from the observed expo- taylor 945D − 15 O ≪ (cid:18) (cid:19) nential decay of the particle numbers and subsequently = 14% of value at α=0, when α (22b) used in the corrected formula (15) to estimate D from →∞ theexperimentally-observedvariance. Alsonotethatthe When there is no absorption (α = 0), v is equal to the e momentsinthecaseofanarbitrarycross-sectionmaybe asymptotic mean instantaneous velocity of the particles computed numerically to the necessary precision by in- since ψ0 is constant. However, since the slower particles volving a finite number of eigenmodes. This is possible near the walls are also the ones that are preferentially because terms involving higher eigenmodes get progres- absorbed, v increases with α as only surviving parti- e sivelylessimportantasmaybeeasilyseenfromtheforms cles contribute to the moments. Similarly, since more of (14) and (15). particles at the walls get removed for increasing α, the dispersionin the convectionfield asseenby the particles goesdownandsodoesD whicharisesduetoit. Our taylor results (22) agree with those derived in [16] by a more tedious method. One may reproduce (15) from [16] with some effort. However, as already noted before, it is very difficulttoobtainsimpleformulæforthehighermoments using the method presented there. Finally, the skewness, which is a measure of the devi- Acknowledgments ation from gaussianity, is found to be ℓ2 2 γ = 0.186 1+1.7α+ (α ) ,α 1 We would like to thank Prof. B. I. Halperin for his in- − Dt O ≪ r sightfulcomments andSrividyaI.Biswasforstimulating ℓ2 (cid:0) (cid:1) discussions. RRBwouldalsoliketothankSchlumberger- = 0.635 , α (23) − Dt →∞ Doll Research for supporting this work. r [1] G. I. Taylor. Proc. R. Soc. London, 219:186, (1953). [14] W. N. Gill. Proc. R. Soc. Lond. A, 298:335, 1967. [2] G. I. Taylor. Proc. R. Soc. London, 225:473, (1954). [15] W. N. Gill. Proc. R. Soc. Lond. A, 333:115, 1973. [3] R.Aris. Proc. R. Soc. London, 235:67, (1956). [16] E.M.LunguandH.K.Moffatt. J.Engg. Math.,16:121, [4] P-G. deGennes. J. Fluid Mech, 136:189, (1983). 1982. [5] P.G. Saffman. J. Fluid Mech., 6:321, (1959). [17] E.L.Cussler.Diffusion: MassTransferinFluidSystems. [6] A. Ajdari, N. Bontoux, and H. A. Stone. Anal. Chem, Cambridge University Press, second edition, 1997. 78:387, 2006. [18] J. B. Haldane. Biometrika, 31:392, (1940). [7] N. Bontoux, A. Pépin, Y. Che, A. Ajdari, and H. A. [19] P. M. Morse and H. Feshbach. Methods of Theoretical Stone. Lab Chip, 6:930, 2006. Physics I, Chap. 7. McGraw-Hill Book Company, 1953. [8] P.G. Saffman. J. Fluid Mech., 6:194, (1960). [20] M.TodaR.KuboandN.Hashitsume.StatisticalPhysics, [9] U. M. Scheven and P. N. Sen. Physica Review Letters, II. Springer-Verlag N.Y., 1991. 89:254501, (2002). [21] Thistypeofprocessiscommoninmanyareasofphysics [10] C. P. Lowe and D. Frenkel. Phys. Rev. Lett., 77:4552, like spectral diffusion[20]. In line-shape problems of (1996). NMR/ESR, x ≡ iφ = lnM would be the ‘phase’, ~y the [11] H.D.Pfannkuch. Rev.Inst.Fr.Petrol.,318:215,(1963). positionofthediffusingparticleswhilev(~y)wouldbethe [12] C. Van den Broeck. J. Appl. Math. Phys., 34:489, 1983. local Larmor frequency ω(~y). [13] C. Van den Broeck. Physica A, 168:677, (1990).