Deposition of Colloids in Porous Media - American Chemical Society

Chapter 3

Deposition of Colloids in Porous Media Theory and Numerical

Solution

Menachem Elimelech and Lianfa Song

Downloaded by UNIV LAVAL on July 14, 2016 | http://pubs.acs.org Publication Date: May 13, 1992 | doi: 10.1021/bk-1992-0491.ch003

Civil and Environmental Engineering Department, School of Engineering and Applied Science, University of California, Los Angeles, CA 90024

The convective diffusion equation describing the transport of sus pended colloidal particles in the vicinity of stationary surfaces in porous media was rigorously formulated by incorporating fundamental theories of mass transport, hydrodynamics of particles in porous media, and colloidal interactions. A finite difference numerical scheme with a variable step size was developed to solve the equation for the concentration distribution of colloidal particles over solid surfaces in porous media, from which the colloid deposition (capture) rates can be evaluated. With the mesh refinement technique employed in the numerical scheme, no oscillations were observed in the numerical solution.

Colloidal particles are ubiquitous in groundwater aquifers. The fate of many pollutants in groundwater is determined in large part by the fate of the colloids with which they are associated (7). Field studies indicate that colloids are mobile in subsurface environments. Mobile colloids can act as a third phase that can enhance transport of adsorbed pollutants. The enhancement in the transport of pollutants sorbed to colloids is usually referred to as facilitated transport or colloid-mediated transport. This process has not yet been considered in pre dictive models of contaminant transport. Colloids interact with stationary surfaces in porous media. The extent of migration of colloids and associated pollutants in groundwater is determined in large part by the rate of their capture by aquifer stationary media. The capture of colloidal particles from flowing suspensions by stationary surfaces is usually referred to as particle deposition. In addition to the role of colloids in pollutant transport, captured colloids can reduce the hydraulic conductivity of porous media due to pore clogging (2,3). Artificial recharge of groundwater, for instance, can be seriously hampered due to pore clogging caused by captured colloidal materials (4). Control of pore clogging is also of considerable interest in the field of petroleum extraction (5). A n excellent review of the various colloidal processes in groundwater and their significance was presented by McDowell-Boyer et al. (7).

0097-6156/92/0491-0026$06.00Α) © 1992 American Chemical Society

Sabatini and Knox; Transport and Remediation of Subsurface Contaminants ACS Symposium Series; American Chemical Society: Washington, DC, 1992.


3. ELIMELECH & SONG

Deposition of Colloids in Porous Media

27

The literature already contains quantitative formulations of colloid depo sition on various surfaces, including rotating disks, stagnation-point flow cells, and parallel-plate channels (5-7). These geometries are convenient for modeling, since the fluid flow field over these surfaces is simple and well defined. A s a result, numerical solutions of the governing differential equations of these systems are straightforward. Furthermore, these simple surface geometries can easily be developed in the laboratory to test the governing equations. The capture of colloids (submicrometer in size) in porous media, however, is more complex and has not been treated rigorously. This study presents a formulation and a numerical solution of the particle transport equation for colloids in porous media. The equation was developed by incorporating fun damental theories of mass transport, colloidal hydrodynamics, and colloidal interactions. In this model, it was assumed that the colloidal particles are spherical and that the porous medium is comprised of uniform spheres (referred to as collectors). A finite-difference numerical method using a mesh refinement technique is developed to solve the transport equation and to calculate the rate of colloid deposition in the porous medium.

Theory The Particle Transport Equation. The transfer of colloidal particles from flowing suspensions toward stationary surfaces is governed by the convective diffusion equation. The equation, in its general form, is given by

^

+

V-J =

fi

(1)

where C is the concentration of particles, t is the time, J is the particle flux vector, and Q is a source/sink term. The particle flux vector is given by (8,9) J = - D · V C + uC + (1/kT) D · F C

(2)

Here D is the particle diffusion tensor, u is the particle velocity vector induced by the flow of the suspending medium, k is the Boltzmann constant, Τ is the temperature, and F is the external force vector. The first, second, and third terms on the right hand side of equation 2 represent the transport of particles induced by diffusion, convection, and external forces, respectively. The external force vector can be derived from the interaction potential energy function, φ, as fol lows: F = -V(|>

(3)

When equations 2 and 3 are substituted into equation 1, the convective diffusion equation under steady conditions and in the absence of a source term, is reduced to


28

TRANSPORT AND REMEDIATION OF SUBSURFACE CONTAMINANTS


V-(uC) = V . ( D - VC) + V

(4)

kT

Surface Interaction Potentials. The potential energy between suspended colloidal particles and collector surfaces resulting from surface interaction forces can be described by the Derjaguin-Landau-Verwey-Overbeek ( D L V O ) theory (8,10,11). In this approach, the total interaction potential is the sum of van der Waals and electric double layer interactions. Quantitative expressions for these interaction potentials, which were used in the theoretical calculations appearing in this work, are described below. The universal van der Waals attraction forces have long been known to play a part in the capture of colloidal particles by surfaces. These forces arise from the fluctuating electromagnetic field between molecules of the colloids and collector surfaces (11). There are several analytical expressions which describe these surface forces between particles and stationary collectors. A n excellent summary of these expressions was given by Gregory (12) and Jia and Williams (9). In this work, the approximation of Gregory (12) for van der Waals attraction will be used. This expression is in good agreement with exact calculations at short distances (up to 20% of the particle radius) and is given by Aa

B

§VDW

~

(5)

6y(l + Uy/X)

where A is the Hamaker constant of the interacting media, a is the particle radius, y is the separation distance (surface to surface) between particles and collectors, and λ is the characteristic wavelength of the interaction, often assumed to be 100 nm. For larger separations, the exact expression for retarded van der Waals interaction derived by Czarnecki (13) w i l l be used. This expression is given by p

§VDW —

A.

2.45λ 60π

y-a

y+3a

p

(y+2a f p

0.59λ (y-3a 3

p

5040π

3

2 . 1 7 λ ί y_ •2a 2

p

n

2

' 72071 [

(y+2a f p

y+5a

p

4

(y+2a )

(6)

p

Czarnecki's expression is inaccurate for distances shorter than about 10 nm; at such short range, however, equation 5 is adequate. Most surfaces and colloids in aqueous media are charged and have electric double layers associated with them (10,11). The surface charge is balanced by an equivalent number of counterions, some of which are located very close to the surface, in the so-called Stern layer, while the remainder are distributed in the diffuse layer. When a charged particle approaches a similarly charged surface, their diffuse layers overlap, and, as a result, a repulsive force develops.


3. ELIMELECH & SONG


29

Analytical expressions for this force as a function of the separation distance are available in the literature (10,14-16). The widely used analytical expression of Hogg et al. (14) for electrical double layer interaction w i l l be used in this study. It gives the potential energy as 1 +exp(-K>Q

Φ£^ = π ε ε α ^ 2 ψ ψ 1 η 0

Γ

1

2

+ (ψΐ + Mfy In [1 - exp (-2Ky )]

1 -exp(-Ky)

(7)

Here ε, and ZQ are the relative dielectric permittivity of water and the permittivity under vacuum, respectively; ψ ! and ψ are the surface potentials of particles and collectors, respectively; and κ is the reciprocal Debye length. This expression is valid for: (1) interaction at constant potential (ii) 1:1 electrolytes; and (iii) surface potentials smaller (in absolute value) than 60 m V . Downloaded by UNIV LAVAL on July 14, 2016 | http://pubs.acs.org Publication Date: May 13, 1992 | doi: 10.1021/bk-1992-0491.ch003

2

Hydrodynamics of Particles i n Porous M e d i a . The role of hydrodynamics in the transport of colloids to stationary surfaces depends upon the geometry of the collectors, the flow field, and the physical properties of the particles and the liquid phase. This study is concerned with the deposition of submicron spherical particles on stationary surfaces in porous media at low Reynolds number flow. This condition (i.e., low Reynolds number flow) is typical of groundwater flow. Several models describing the flow field in porous media at low Reynolds numbers are available (17-19). Among these, Happel's sphere-in-cell flow model (17) is the most commonly used. Happel's model has also been used successfully in particle filtration models in which the transfer of particles to stationary collectors is involved. In Happel's model, the porous medium is treated as an assemblage of identical spherical collectors, each of which is enveloped in a shell of fluid. The thickness of the shell is determined so that the overall porosity of the porous medium is maintained for the single collector (17,20). The undisturbed fluid flow field in Happel's model is obtained from a solution of the Navier-Stokes equation. The stream function for axisymmetric, steady creeping flow around a spherical collector in this model is given by (19)

Ψ =

1-υαϊύη θ 2

2

+ K

+ K

A

0

(8)

r

V J

where a is the radius of the collector; U is the approach velocity of the fluid; r is a radial coordinate originating from the center of the spherical collector; and Κ χ, K , K and K are constants that depend on the porosity of the porous medium. The radial and tangential components of the fluid velocity around the spherical collector (v and ν , respectively) can be derived from equation 8 as follows: c

2

3

4

r

θ

1 θ

3Ψ

r sin θ dr 1

3Ψ

2

r s i n 6 3Θ


(9)

(10)


30

where r and θ are the radial and tangential coordinates of the spherical coordinate system used in this work. Even in the absence of surface interaction forces, the trajectories of particles at short distances from collector surfaces do not follow the fluid streamlines. This deviation is caused by hydrodynamic (viscous) interactions (27). In a viscous fluid, such as water, the approach of a particle to a surface is hindered by the slow drainage of water from the narrowing gap between the particle and surface. A s the gap narrows to zero, this drainage becomes infinitely slow because of the no-slip conditions at the surfaces and, as a result, the drag on the particles significantly increases. Considering hydrodynamic interaction, the radial and tangential particle velocity components (u and u , respectively) are related to the fluid velocity components by (6,22)


r

«,=/i(ff)/ (ff)v 2

"θ=/ (^)ν 3

Q

(11)

r

(12)

θ

Here f (H\ / (7/), and / (7/) are universal correction factors for hydrodynamic interaction; these depend on the dimensionless distance (7/ = y la ) between the particle and the surface of collectors (27, and references therein). The decrease in the mobility of particles in the vicinity of the collector surface due to hydrodynamic interactions results in a decrease in their diffusion coefficient. In this case, the diffusion coefficients are expressed as (9) x

2

3

D =f (P)D r

x

(13)

m

D =f (fl)D B

4

(14)

m

where/ι(77) is a hydrodynamic correction factor, and is the particle diffusion coefficient at large distances from the collector. The latter is a scalar and can be obtained from the Stokes-Einstein equation (19): kT 6πμα

ρ

Mathematical Formulation Particle Transport Equation i n Spherical Coordinates. For a porous medium comprised of spherical collectors, the dimensionless convective diffusion equation in spherical coordinates can be written as

Vm

. 3Φ

. SAW) 3Φ

se" 3Φ

.

Ϊ


3. ELIMELECH & SONG


31

Here C* = C / C , C being the concentration of particles in the approaching fluid; 77 is a dimensionless separation distance from the collector surface (defined previously); Φ = φ/£Γ is a dimensionless potential energy, φ being the total interaction energy between particles and collector; V = v /U and V = VQ/U, v and ν being the radial and tangential velocity components of the fluid defined by equations 9 and 10; N = (2a U)/D„ is a particle Peclet number; and iV/j = (a + a )la . In the derivation of this equation from equation 4, the fol lowing assumptions were made: (i) colloidal forces act perpendicular to the surfaces; and (ii) particle diffusive fluxes tangent to the collector surface are much smaller than the radial (perpendicular) diffusive fluxes. For convenience in numerical calculations, equation 16 can be rewritten as 0

0

R

R

Q

r

θ

Pe

c

p

p

p

t)C* TfC* 7)C* ^ - = α (Η θ) ^ + a (H Θ) + a,(H 0)C* Downloaded by UNIV LAVAL on July 14, 2016 | http://pubs.acs.org Publication Date: May 13, 1992 | doi: 10.1021/bk-1992-0491.ch003

λ

9

2

9

9

(17)

where the functions a (H Q) a (H Q) and a (H Q) can be determined from equation 16. The boundary conditions in our case are x

9

9

2

9

9

3

9

C*(// = 0,9) = 0

(18a)

C*(//->oo,9) = l

(186)

In the first boundary condition, the concentration of particles in contact with the collector is taken to be zero because these particles are no longer part of the dispersed phase. The deposition at H=0 is assumed to be irreversible (the socalled "perfect sink" model). The second boundary condition states that at large separations from the surface of the collectors, the concentration of particles is equal to that of the approaching fluid. In Happel's porous medium model, this boundary condition can be taken at the outer surface of the fluid envelope; the thickness of the fluid envelope is determined by the porosity and collector diameter (79). The third boundary condition arises from the symmetry around the forward stagnation path of the spherical collector (22). Calculation of Particle Deposition Rates. Once the dimensionless concen tration distribution of particles around the collector, C*(77, Θ), is determined, the perpendicular flux of particles at the collector surface can be calculated. When the local particle flux at the surface (i.e., at H=0) is integrated over the entire surface of the collector, the overall rate of particle deposition is obtained. The local particle flux perpendicular to the surface can be expressed as

J(H,Q) =

UC J\H Q) 0

9


(19)

32



where / * ( / / , Θ) is a dimensionless flux given by

In this equation, the first term describes the diffusive flux of particles, the second term describes the convective flux of particles, and the third term describes the flux of particles due to migration (drift) velocity resulting from the interaction potentials (i.e., van der Waals and electrical double layer). The dimensionless local flux of particles at the collector surface is evaluated from equation 20 with H-+0. The overall rate of particle deposition on the collector, I, can be obtained from integration of the local flux (at H=0) over the entire surface of the collector as follows: π 2

I = 2na jj(H ο c

= 0, Θ) sin QdQ

(21)

In studies concerned with deposition of colloidal particles in porous media, it is convenient to use a dimensionless deposition rate, also known as the single collector efficiency,

(

^T^TJF

2

2

)

This can be viewed as the ratio of the overall deposition rate of particles on the collector to the convective transport of upstream particles towards the projected area of the collector. Numerical Solution General Considerations. Equation 17 is similar in form to the one-dimensional advection-dispersion equation with a first-order chemical reaction:

3C

dC D

— -k C

— — -v dx dt = d dx d

2

Vx

x

r

(23)

Here D is the dispersion coefficient, v is the average fluid velocity in the direction x, and k is the reaction rate constant. The advection-dispersion equation is used to describe the mass transport of solutes in various processes. The features of the equation have been the subject of numerous reviews (23,24). In the numerical solution of equation 17, some concepts and results w i l l be adapted from previous studies of the advection-dispersion equation. When the ratio of advection (the second term in the right hand side of equation 23) to dispersion (the first term in the right hand side of equation 23) is small to moderate (so-called "dispersion-dominated transport"), no numerical d

x

r


3. ELIMELECH & SONG

33

difficulties in the solution of this equation are encountered. Several numerical methods to solve this equation are available for this case (23,24)· On the other hand, when this ratio is large (so-called "advection-dominated transport"), severe numerical difficulties are encountered when these methods are employed. These numerical difficulties are of two types: overshoot and numerical dis persion (23,25). A mesh Peclet number is usually used in numerical calculations to characterize the ratio of advection to dispersion. The mesh Peclet number is defined as v Ax/D , Ax being the step size in the direction x. Most numerical methods give an accurate solution when the mesh Peclet number is smaller than unity. However, when the mesh Peclet number increases, oscillations appear in the numerical solution. In the problem of capture of colloidal particles by stationary surfaces, the interaction potential near the surface is extremely large due to the nature of the universal van der Waals attraction forces. Consequently, the velocity of the articles at the collector surface becomes infinitely large. In addition, the difjsion coefficient near the surface approaches a very small value due to the retardation of particle mobilities caused by hydrodynamic interactions. A s a result, at the collector surface, the coefficient a (H, Θ) in equation 17 approaches zero while a (H,Q) approaches infinity. Thus, close to the collector surface equation 17 behaves as the advection-dispersion equation with high Peclet number; oscillations in the numerical solution of equation 17 are obtained in this region when conventional methods are used. These oscillations affect the accuracy of the numerical solution for the concentration distribution of particles, and consequently for the particle deposition rates in the porous medium. In order to overcome the numerical difficulties in the solution of equation 17, a particular mesh refinement technique near the collector surface is used in the numerical scheme. x



d

Ε

x

2

Numerical Discretization. Equation 17 with the boundary conditions described by equations 18a-18c w i l l be solved in two stages. In the first stage, by applying boundary condition 18c, the equation is reduced to an ordinary differential equation:

(H,0)

ai

+ a (H, 0) 2

+ a,(H, 0) C* = 0

(24)

Equation 24 with the boundary conditions 18a and 18b is a typical stiff two-point boundary value problem. A central difference scheme on a non-uniform mesh, developed by Pearson (26), w i l l be used to solve this problem. The solution yields the concentration distribution of particles at the forward stagnation point of the spherical collector (i.e., at 6=0). In the second stage, the concentration distribution of particles over the entire region of the spherical collector w i l l be calculated (i.e., for 0

Deposition of Colloids in Porous Media - American Chemical Society

Recommend Documents