18 Groundwater Hydrology V (Advection, Dispersion, Diffusion and Sorption)

Prabhas K Yadav

Objectives:

This module will present an introductory topics relevant to transport problems in which mass(es) are reacting during the transport in groundwater. Visualization of few conservative problem is to be learned. The target groups are higher level undergraduate students and the first year PG students. The specific objectives of the module are:

1.      Systematic introduction to transport problems in groundwater

2.      Recognizing factors and processes affecting reactive transport problems

3.      Quantifying and comparing different transport processes with inclusion of sorption

4.      Developing a systematic approach to solve reactive transport problems.

18 Introduction:

In the study of hydrogeology or groundwater, mass transport or simply transport (also used in this module) is a standard term that describes flow of mass along with groundwater. With the transport we deal with quality aspect of groundwater, as compared to flow in which quantity is dealt. The transport question becomes very important when groundwater is used as a source of potable water, which is required to be maintained at a standard quality but more importantly which is to be protected against the quality deterioration or contamination. Thus with transport we deal with not only the physical aspect of hydrogeology, both subsurface physical properties (e.g., homogeneity versus heterogeneity) and flow dynamics as it is with a flow problem, but we also interlink these with chemical and biological aspects of hydrogeology. A fairly representative transport problem will often require that we set-up at least a 2D system.

The transport problem is inherently a multi-disciplinary one, and thus it is important that we first identify different processes that are important in describing the problem. In this module we will focus on non-reacting systems, also called a conservative system, which is purely a physical transport problem.

18.1 Transport Processes

Aquifers as we know already are extensive, at least in the horizontal and lateral directions. However, it’s properties such as conductivity, transmissivity, often changes in a very small scale, in many cases also at a pore-scale. As such to understand aquifer properties and underlying processes, we define a Representative Control Volume (RCV, see Fig 1). Ideally a RCV is very much smaller than the aquifer that is investigated, and at the same it is very much larger than individual grain or pore. The essential condition is that Darcy’s law has to be applicable in the RCV. A big advantage of RCV is that in the aquifer processes analysis it can be oriented as per the coordinate system, e.g., rectangular for Cartesian coordinate.

Before we step forward to learn the transport processes, let us perform a couple of simple laboratory column experiments. This should help us understand more clearly the transport processes. Let us assume that we have a column (say diameter = 3 cm, and length = 10 cm) packed homogeneously with aquifer material, e.g., coarse sand. Let us continuously introduce a tracer, say NaCl (concentration C0), from the entry end of the column instantly at time, t0. For tracking the movement of the tracer, we will periodically collect effluent from the column and analyse for NaCl concentration. The results of our hypothetical experiments are presented in Fig. 2. We observe that mass particles in the column was transported in two different ways. In the first case (Fig. 2b) the effluent and inlet concentrations are always equal for all times except the initial lag time, when there are no NaCl leaving the column.

Clearly, the mass particles are uniformly pressed out of the column in this case. The particles obtained in effluent are called advected mass, and the transport process is termed advection. The other case (Fig. 2c), suggests that mass particles are travelling with different speeds. In this case the mass particles get spread and therefore they exit column at different times. We refer to this type of transport as dispersive transport. It is to be noted that reaction during transport, e.g. adsorption, significantly impacts transport. We will now quantify these transport processes including the reaction with transport.

18.1.1 Advection

As we learned from our above experiment, advection (or also called convection) is the transport of mass with the movement of a moving medium. The moving medium can be termed carrier, and in the groundwater case it is usually fluid. To quantify advection process, we will consider our laboratory column (Fig. 2) which has a constant cross-sectional area, A [L²], the steady-state discharge, Q [L³T^-1]. Further we will assume ne [ – ] the effective porosity of the column packing and v as average linear velocity of water flow. Note that v= q/ne , where q [LT^-1] is the Darcy velocity or also called specific discharge. With these information, we quantify the advective mass flow rate Jadv [MT^-1] as

We have so far considered a uniform cross-section but because continuity equation has to be satisfied, it is not required that flow area remain uniform. This is better understood from Fig 3. Assuming that ne is spatially constant, the mass (shown in grey patch) is the same in both flow sections. During the transition from the smaller to the larger cross-section, the area covered by the solute mass expands laterally and its extension along the flow direction is reduced accordingly, and therefore the concentration is not spatially changed.

18.1.2 Mechanical Dispersion

There is a general agreement in describing the mechanisms behind mechanical dispersion. The three most common mechanisms (also see Fig. 4) that can be found in standard literatures (e.g., Domenico and Schwartz, 1998) are:

I.       Non-uniform flow in an individual pore: From fluid mechanics it is know that a parabolic velocity profile results for a flow passing between two parallel solids (Fig. 4a). Thus mass passing through the centreline of the profile will be fastest and the one at the edge will be the slowest.

II.    Effect of Pore size distribution: Non-uniformity in pores size leads to different flow velocities for transport of mass in wider pores (higher speed) compared to narrower pores (Fig 4b).

III. Effect of flow paths: Due to irregular distribution of pores in the porous media, mass particles are likely to take flow paths varying in lengths thus leading to meandering of particles (Fig 4c).

Next, we quantify dispersive mass flow rate, Jdis [MT-1] using our experiment (Fig 2). Quite clearly the last two mechanisms (point II and III, above) suggest that dispersion can be explained as a random phenomena, which are explained in standard texts, such as Domenico and Schwartz, (1998). Based on observations, the dispersive mass flow due to mechanical dispersion in porous media is found to be

I.   proportional to the cross-sectional area, A, of the flow

II. proportional to the linear velocity, v

III. proportional to the difference in concentration at different location, ΔC = C0 – C

IV. inversely proportional to the transport distance, L

The four proportionality relations mentioned above lead to

with α a constant of proportionality called dispersivity [L] . The negative sign is usedin equality relation to identify that transport is from a higher concentration towards the lower one. Dispersivity is a property of porous media alone and is not dependent on fluid type or the flow characteristics. In a more formal manner the ratio ΔC/L [ML^-4] is termed as concentration gradient. This let us quantify dispersivity as a quantity whose value equals the dispersive mass flow rate through a unit cross-section area for a unit concentration gradient and a unit linear velocity.

As mass transport is influenced by both porous medium and flow properties, a product of dispersivity and linear velocity called the mechanical dispersion coefficient, Dmech [L²T^-1], is more often used in analysing dispersive transport. This changes eq. (5) to

 

Fig. 4: The threeWatmechanismsrResourcesleadingandManagementtodispersion (Adapted from Kinzelbach, 1992)

Since α characterizes spread length, its value is dependent on the flow coordinate direction under consideration. Furthermore, an excellent compilation and analysis by Gelhar et al., (1992) shows a scale (travel distance) dependence of α. The compilation, which is highly referred, suggest the longitudinal α (along the direction of groundwater flow) varies between 10^-3 m and 10³ m with bulk of values between 10^-1 m and 10² m.

Advection and the mechanical dispersion are two processes that are generally used in analysing groundwater transport problems. Nevertheless, it is important for us to introduce a very important transport mechanism called diffusion that is not flow or mechanically driven but rather a concentration gradient driven.

18.1.3 Diffusion

Diffusion can be considered a physiochemical process. The quantification of diffusive mass transport, [ MT^-1] is based on the work in Fick, (1855). In that work the mass flow due to diffusion (1D) is described to be

I.        proportional to the cross-sectional area of the flow, A

II.     proportional to the concentration difference between the transport distance, ΔC = C0 – C

III.   inversely proportional to the transport distance, L

The above three points are combined to constitute the Fick’s first law of diffusion, which is

18.1.4 Sorption

In the first section we distinguished between conservative and the reactive transport systems. Essentially, a reactive system is one in which the mass entering the domain undergoes changes resulting to either increase or decrease of the original mass. Both cases are possible due to some form of (bio)chemical reactions taking place within the domain. Thus we are now exploring a very broad topic on aquatic chemistry and biology. For this introductory module, we will restrict to sorption type reaction only.

Equilibrium Sorption

In the study of groundwater, at least at the introductory level, we use a common term sorption for adsorption and absorption. By sorption, we will actually be learning on adsorption reaction. To begin let us start with common terminologies. Adsorption is the process of accumulation of dissolved chemicals on the surface of a solid, e.g. accumulation of a chemical dissolved in groundwater on the surface of the aquifer materials. In most cases, adsorption is a reversible process, i.e. adsorbed chemicals can become dissolved again. This process is termed desorption. Two components of adsorption-desorption are adsorbent and adsorbate. The first one is the solid offering adsorption sites and the latter one is the chemical that get accumulated on the adsorbent. Interactions between adsorbent and adsorbate depend on properties of both substances. Fig. 5 presents a schematic representation of a sorption process along with its various components.

An equilibrium between the adsorbent and adsorbate on the surface and in the solution is approached with time. The equilibrium condition can be used to quantify the sorption process. The so-called isotherms, which provide a relationship between the solute concentration C [ML-3] in the dissolved phase and the adsorbate mass toadsorbent mass ratio Ca [ ] . The concept of isotherms is derived from the Chemical Engineering studies, which provide several different types of isotherms. In groundwater transport problems Linear, Freundlich and Langmuir isotherms are the common ones. Next we will learn about these isotherms.

The Linear isotherm

The Linear isotherm also called the Henry isotherm or the Kd model is based on the linear relationship between C and Ca. Mathematically it is represented by

in which Kd [L³M^-1] is called the distribution or portioning coefficient. Although simple but this isotherm has several limitations. Among the most important ones are that this isotherm assumes that sorption is unlimited, i.e., an adsorbent has infinite number of sites for sorption. More often in groundwater transport problems the available number of sites are not a limiting factor, however that linearity is always followed may not be the case.

Largely because of its mathematical simplicity (linearity), the linear isotherm is extensively used in analysing groundwater transport problems (see EPA, (1999) for further learning). Let us perform a simple experiment to obtain a in the laboratory. Let us equilibrate six different concentrations (Ci, mg/L) of a solute in six different flasks each with 10 g of an identical adsorbent. The results of the experiment is presented in Table 1. We then plot the Ca (mg/g) against C (mg/L) and linearly fit the data. The slope of the fit, which is 1.404 is the Kd (see Fig. 6).

The Langmuir isotherm

The main limitation of the linear isotherm that the absorbent has infinite capacity to adsorb is avoided in the Langmuir isotherm. This is achieved based on the assumptions that ions are adsorbed as a monolayer on the surface, and the maximum adsorption occurs when the adsorbent surface is completely occupied by the adsorbate. Mathematically the Langmuir isotherm is

in which KL [L³ M^-1] is the Langmuir coefficient and Ca,m is the maximum number of ions that can be adsorbed on the surface. As we did previously, we will perform a hypothetical experiment to estimate Ca,m and obtain KL. In this experiment we will add 250 mg/L of solute to six flasks containing different amounts of absorbent. The experiment results and calculations are tabulated in Table 2. Fig. 7 presents the plot of results in which the Ca,max is found approximately at 0.16 mg/mg. We are still required to find KL. This will require fitting the data that we have in Fig. 7. Although it is possible to fit the results plot with eq. (13), a much simpler approach is to linearise the equation. This is easily achievable by first inverting both side of eq. (13) and rearranging the terms. The linearised form of eq. (13) is

By comparing the linear fit equation (see right plot Fig. 7) of experimental results with eq. (13), both C_a,m and K_L can be easily obtained. Inverting the intercept of the fit line provides C_a,m = 0.18 mg and then KL =0.11 mg can be obtained by equating the slope of the fit with 1/ (K_L C_a,m).

Freundlich isotherm

The Freundlich isotherm is a modified form of KD model and it observes the non-linearity involved with sorption process. It is based on the idea of a power law relating the C and Ca. The isotherm is given as

in which K_F is the Freundlich constant and n is the measure of the non-linearity involved. As is with the K_D model, the Freundlich isotherm is not able to provide any information on maximum possible sorption. The Freundlich isotherm becomes equivalent to K_D model for n =1. For illustration and obtaining K_F, we will work on a similar experiment that we performed for demonstrating the Langmuir

isotherm (see Table 3, and Fig. 8). As we did earlier we will linearise eq. (15), and then use a linear fit to obtain the KF and n from fit parameters. Taking a logarithm (base 10) on both side of eq. (15) and using the properties of logarithm, we can linearise it. The linearised form of eq. (15) is

Fig. 8 and Table 3 illustrates the processes for obtaining KFand n using eq. (16). From the analysis we see that K_F = 0.023 and n = 0.62

Retardation factor

The most significant effect of the sorption is that it leads to retardation of the mass being transported. This basically means that mass will be transported slower than the average groundwater velocity. The coefficient of retardation R is used to estimate retarded contaminant velocity. We first derive R and then find the retarded mass velocity. For this purpose we will use our above discussion on linear isotherm (or K_D model). Further we will apply mass balance equation.

Let us assume aquifer of volume V with an effective porosity ne and material density ρ. Thus we have:

total water volume: n_e V

mass of dissolved chemical: n_e V C

volume of solid: (1 – n_e)V

mass of solid: ρ(1 – n_e)V

mass of adsorbate:    ρ(1 – ne)V Ca = (1 – ne)V ρ KD C

Total mass:  ne V C + (1- ne)  ρ V KD C = ne R V C

with R = 1+ ρ KD(1-ne)/ ne. Using bulk density (mass of solid/total volume) ρb = (1 – ne)V ρ, R becomes

18.2 Quantifying and visualization of Transport processes

18.2.1 Joint action of advection, dispersion and diffusion

The transport of conservative solute in an aquifer can be understood as a superposition of advection, mechanical dispersion and diffusion. To quantify the transported mass over the time, we sum all the mass flow rate, J [MT^-1] that we have considered. i.e., eqs. (1, 6 and 9). Thus we get a combined equation for the mass flow

The spreading of mass due to advection and dispersion can be quantified by combining a mass budget and the corresponding laws of motion. This results in a transport equation called Advection – Dispersion (AD) or convection-dispersion equation. We will attempt to quantify the relative importance of different mass flow rate using eq. (20) with the following 1D column experiment:

Cross-section area of flow, ne A = 1 m²

Linear velocity, v = 1 m/d

Concentration,  C = 1 mg/L = 1 g/m³

Transport distance, L = 1 m

Concentration gradient, ΔC/L = 1 g/m⁴

Dispersivity, ɑ = 0.1 m (ɑ= L/10 is a rough estimate)

Diffusion coefficient, Ddif = 10^-10 m²/s

The contribution of advection can be obtained from eq. (1), we get

in which d[L]is mean grain diameter, and other quantities was defined earlier. As can also be observed in Fig. 9, the following three regimes are defined:

in which d [L] is mean grain diameter, and other quantities was defined earlier. As can also be observed in Fig. 9, the following three regimes are defined:

I.   For Pe < 0.02, diffusion dominates

II.  For 0.02 < Pe < 6, diffusion and dispersion dominates

III.   For Pe > 6, advection and dispersion dominates

18.2.2 Visualizing and quantifying transport processes

Since advection is simply a transport along with groundwater flow, the advective mass at any point is
the initial mass, and the position at any time instant is the product of linear flow velocity, v , and t .

This approach is used if advective mass flow is to obtained in transverse directions. However, for most of the practical problems groundwater flow is considered unidirectional and along the horizontal direction. The quantification of dispersive mass flow, more often the dispersion and diffusion coefficients are rather complicated. We have already summed the two processes in eq. (20), furthermore we have hinted that these can be better described using statistical approach. In fact hydrodynamic dispersion coefficient (we will call it dispersion) has been related to a Guassian distribution statistics as (see, Domenico and Schwartz, 1998):

where σ refers to square root of standard deviation or variance. If position is intended to be obtained, then t can be replaced with v/x. It is to be noted that eq. (22) represents the spatial (or space) distribution of the mass, but when mass evolution is to be obtained from temporal distribution, the following space-time interchange (from Robbins, 1989) has to be applied

Fig. 10 presents a schematic of mass evolution from a 1D column. As can be observed in the Fig 10c, breakthrough curves are representations of the solute concentration as a function of time at specified observation locations. Concentration profiles (Fig. 10b) are representations of the solute concentration as a function of a space coordinate at fixed time levels. Examples of using these plots for calculating dispersion coefficients can be found in Appelo and Postma (1999). Readers are encouraged to use an spreadsheet that can simulate and plot breakthrough curves and concentration profile developed by Prof. J. Craig (Uni. Waterloo, Canada. http://www.civil.uwaterloo.ca/jrcraig/pdf/OgataBanks.xlsm).

While advective transport can essentially be 1D along groundwater flow, tedious efforts are required to understand and quantify dispersion coefficients in 2D and 3D. We can extend our above discussions to quantify transverse dispersivities. As a rule of thumb, the horizontal transverse dispersivities (ɑTh) can be assumed to be one-tenth of longitudinal dispersivity (ɑL) and, likewise, the vertical transverse dispersivity (ɑTv) can be approximated as one-tenth of ɑTh.

18.3 Transport Models and Solutions

So far our discussions have been focussing on processes involved in mass transport in aquifers. In the natural transport problems several processes affect the system simultaneously. Thus our next step will be to combine several processes and attempt a solution of such system. We begin with derivation of the 2D transport model.

18.3.1 Derivation of transport model

Let us consider the same approach as we used in Advection and Dispersion sub-module. The difference here will be in the reactive part. For the reactive part we will consider sorption represented by the KDmodel (eq. 12). The mass balance equation in the RCV is then (see Fig. 11)

Eq. (32), a 2D transport equation also called Advection-Dispersion-Reaction (ADR), can be converted to 1D transport equation by removing the y-components and likewise it can be converted to 3D transport equation by adding the z-components. These are explained in previous sub-module. Next we learn about other essential requirements for transport problems.

18.3.2 Initial and boundary conditions

A complete set of transport problem includes the transport equation, e.g., eq. (32), the initial conditions (time dependent) and the boundary conditions (space dependent). The initial condition describesthe distribution of mass or specifies concentration in the entire area of investigation at some starting time, usually t = 0 is considered. The initial conditions are required for the transport problems that are time dependent. These are so called transient problems and they comprise of general transport problems. A mathematical statement of initial condition is

i.e., att=0, C₀ is the concentration in the entire investigated area.

surrounding just adjacent to the area under investigation. As the name suggests, the boundary conditions have to be specified for the entire boundary of the investigation area, even for the cases when the investigated area in unlimited or extend are infinite. The boundary conditions may or may not be time dependent. For transport problems three types of boundary conditions have been defined.

A first type of boundary condition, also called Dirichlet type, specifies the value of mass or concentration at the boundary. The mathematical statement of the first type boundary condition is

i.e., at any t, the concentration at specified space (x,y,z) is C0. Tracer injection, effluent concentration from polluted area can be few examples of first type boundary conditions. In a second type boundary condition, also called Neumann type condition, the component of the concentration gradient perpendicular to the boundary is specified. The mathematical statement of the second type boundary condition along x-axis is

Since the gradient of concentration is proportional to diffusive flux, the second type boundary condition is also called a flux specified boundary. Second type boundary conditions are also used for specifying the no-flow condition, i.e., no concentration gradient exist across the boundary. However it has to be noted that absolute no-flow can only exist when there is no concentration gradient and the velocity vector is zero. Rock formation at the bottom of aquifer will result to a no flow type of boundary condition.

The third type, Cauchy or Robin type, conditions combines the first two conditions, i.e. both a specific value and gradient is specified. This condition generalizes the other two boundary conditions. Mathematical statement of the third type boundary condition (along x-axis) is

18.4 Few analytical solutions of transport problems

In this section we will attempt to solve the AD model. The techniques for solving the transport problems varies with the complexity of the problems, and therefore we will restrict to simpler problems. We will very briefly get introduced to complex problems (higher dimensions).

Direct integrating of the ADR equation (e.g. eq. 32) is only possible in very simplified cases (e.g., only advection case). However, transformation methods (Laplace, Fourier etc.) can be applied to directly solve ADR equation under some conditions. Thus obtained solutions are called analytical models, which are exact and often a closed-form. The conditions very generally include that aquifer is of regular geometry, is homogeneous and isotropic and the chemical reactions are linear models based. Numerical methods, e.g. Finite-Difference Method (FDM), Finitie-Element Method (FEM), that provide approximate solution of the AD equation are very commonly used to solve transport problems. In this module we will restrict to few analytical models.The solution provided in the succeeding paragraphs generally assumes the following:

Fluid of constant density and viscosity

Flow in -direction only, and velocity is constant.

The longitudinal dispersion coefficient Dx is constant.

In addition we assume aquifer is homogeneous and isotropic. Additional assumptions is provided when required.

18.4.1 The Adevetion-Dispersion (AD) equation

The simplest AD equation can be a transient 1D equation is

and boundary conditions

The conditions implies that our investigation area is semi-infinite and contain C0 concentration initially and Cin concentration is continuously injected to the investigation area from the location x=0. The following analytical solution of this problem is provided in Ogata and Banks, (1961):

Readers should check the wikipedia site (https://en.wikipedia.org/wiki/Error_function) to find methods to approximate the value of erfc.

In order to visualize eq. (38) we perform a column (3 m long) experiment with C0=0 mg/L, v =1 m/h and D =1 m²/h. Let ourCin = 1, 3, 5 mg/L. We will check the concentration in the column at different times. Fig. 12 provides the results of our experiment. The result suggest that it will require 5 hours for our homogeneous distribution of mass in the column.For further visualisation of the Ogata and Banks solution, check GNU-Octave (code (ogata1D.m) provided at https://github.com/prabhasyadav/UGC-Transport for additional simulations.

18.4.2 The AD equation with first-order decay and linear sorption

Next step will be to introduce sorption to transport equation. In this case our transport equation becomes

Fig. 13 Visualizing transport model with decay and sorption using Kinzelbach (1992)

The 1D analysis can can be extended to the 2D and 3D cases. In these cases we bring in use the properties of normal distribution. For a 2D case, we get

Summary:

The goal of the module is to introduces the groundwater transport problem, which is key to analysing, quantifying and solving the groundwater quality problems. To list the module provided:

1. Introduced the conservative and reactive transport problem in the groundwater.
2. Introduced several processes, advection, diffusion/dispersion and reaction (sorption) that influences the transport problem.
3. The derivation of governing equation for the conservative transport problem including the initial and boundary condition
4. Few analytical solution of conservative transport problem with their visualization.

References:

1. Anderson, M. P., and W. W. Woessner (1992), Applied groundwater modeling: simulation of flow and advective transport, Academic Press, San Diego.
2. Appelo, C. A. J., and D. Postma (1999), Geochemistry, groundwater and pollution, Balkema, Rotterdam.
3. Domenico, P. A., and F. W. Schwartz (1998), Physical and chemical hydrogeology, Wiley, New York.
4. EPA (1999), Understanding variation in partition coefficient, Kd, values, U.S. Environmental Protection Agency.
5. Fick, A. (1855), Ueber Diffusion, Ann. Phys. Chem., 170(1), 59–86, doi:10.1002/andp.18551700105.
6. Gelhar, L. W., C. Welty, and K. R. Rehfeldt (1992), A critical review of data on field-scale dispersion in aquifers, Water Resour. Res., 28(7), 1955–1974, doi:10.1029/92WR00607.
7. Genuchten, M. T. van, and W. J. Alves (1982), Analytical solutions of the one-dimensional convective-dispersive solute transport equation, U.S. Department of Agriculture.
8. Kinzelbach, W. (1992), Numerische Methoden zur Modellierung des Transports von Schadstoffen im Grundwasser, Oldenbourg, München.
9. Ogata, A., and R. B. Banks (1961), A solution of the differential equation of longitudinal dispersion in porous media, U.S. Geological Survey.
10. Robbins, G. A. (1989), Methods for determining transverse dispersion coefficients of porous media in laboratory column experiments, Water Resour. Res., 25(6), 1249–1258, doi:10.1029/WR025i006p01249.
11. Wexler, E. (1992), Analytical solutions for one-, two-, and three-dimensional solute transport in groundwater systems with uniform flow, in Techniques of Water-Resources Investigations of the United States Geological Survey, p. 190.