- Email: [email protected]
A non-linear theory of high-concentration-gradient dispersion in porous media
0309-1708(95)00012-7
ELSEVIER
Advances in Water Resources, Vol. 18, No. 4, pp. 203-215, 1995 Copyright 0 1995 Elsevier Science Limited printed in Great Britain. All rights reserved 0309~1708/95$09.50+0.00
A non-linear theory of high-concentration-gradient dispersion in porous media S. Majid Hassanizadeh & Anton Leijnse National Institute of Public Health and Environmental Protection (RIVM), PO Box 1, 3720 BA Bilthoven, The Netherlands (Received 23 September 1994; accepted 15 March 1995)
Tie application
of Fick’s law to describe hydrodynamic dispersion in porous media is based on the assumption of a linear dependance of a solute dispersive m:ass flux on its concentration gradient. Both theoretical and experimental studies ha.ve shown that the Fickian description of dispersion is not valid when large concentration variations in the porous medium are encountered. However, an appropriate alternative is still lacking. In this work, based on a theoretical derivation of the Fickian dispersion equation, a non-linear theory of dispersion is suggested. In the non-linear theory, in addition to the longitudinal and transversal dispersivities, a new parameter is introduced. Miscible displacement experiments are carried out in order to investigate the effects of large variations in salt mass fraction and to assess the validity of the new theory. Low-concentration liquid is displaced upwards in a vertical column by a high-concentration liquid. Thus, only hydrodynamically stable flow regimes are considered. The experiments are simulated by means of both classical Fick’s law and the new non-linear theory. It is found that low-concentration-gradient experiments can be simulated satisfactorily using the Fickian-type dispersion equation. However, calculated breakthrough curves for high-concentration-gradient experiments deviate substantially from the measured curves. It appears that a satisfactory fit to high-concentration-gradient data can be obtained only if the value of longitudinal dispersivity is reduced by a factor of three. Using the nonlinear theory, however, it is possible to simulate both low- and highconcentration-gradient experiments with a unique set of parameter values. Key wor& dispersion, solute transport, Fick’s law, Darcy’s law, brine transport, salt transport.
and I is the unit tensor. Although aL and or are usually assumed to be independent of the fluid properties, it is known that they are time and scale dependent.” The Fickian-type eqn (l), according to which the dispersive mass flux is linear in concentration gradient and independent of the concentration itself, has been successfully used in describing the dispersion of lowconcentration solutes in porous media. Indeed in many cases of ground water pollution, solute transport in unsaturated media, or tracer experiments, the concentration of dissolved matter is low so that density effects may be neglected. However, there exist certain situations of practical interest where high salt concentrations in groundwater are encountered. An important example is the presence of concentrated brines in deep aquifers overlying salt formations. The problem is of utmost interest in the safety assessment studies of high-level radioactive waste disposal in salt formations. In Fig. 1
INTRODUCTION
Hydrodynamic dispersion of solutes in porous media is normally described by a Fickian type equation with a dispersion tensor replacing the diffusion coefficient’ J=-D.VC
(1)
where J is the dispersive mass flux vector, D is the dispersion tensor, and C is the mass concentration. The dispersion tensor is assumed to be independent of the solute mass fraction and is given by17 D = Wo + oz.q)I + (QL - oT)qqlq
(2)
where 12 is porosity, DO is the effective molecular diffusion coefficient, (Ye and or are the longitudinal and transversal dispersivities, q is the apparent mean velocity of fluid with respect to the rock (also called Darcian velocity or specific discharge), q is its magnitude 203
204
S. Majid Hmsanizadeh,
50
1
0-I ~.
.I’.
" 2 cd
-50. ., -i.
-250
-
1000
1200
1100 density (kg/m3)
Fig. 1. Distribution of groundwater mass density in strata overlying the Gorleben salt dome in Germany (courtesy of Bundesanstalt fur Geowissenschaften und Rohstoffe, Hannover, Germany). the mass density of groundwater in the strata overlying the Gorleben salt dome in Germany is given.18 It is evident that the mass density increases drastically (due to an increase in salt concentration), from about 1040 kgmP3 to 1180 kg mm3 over a distance of 30 m. The question then arises whether equations such as Darcy’s law and Fick’s law are still valid under such large concentration gradients. In fact, many experiments have shown that the linear relationship (1) does not hold anymore when there exist large differences in density and/or viscosity of the
A. L.eijnse
displacing and displaced fluids2-5,‘3~16920 (for a review see Kempers13). In general, it appears that dispersivity coefficients decrease as the density (or concentration) differences increase. A similar observation has been made in a recent experimental study by Hassanizadeh et al. I1 where low-concentration salt water is displaced by high-concentration brine in an upward direction. Traditionally, in an attempt to allow for these effects, most researchers still employ eqns (1) and (2) but allow for a dependence of oL and or on density differences, viscosity ratios and the mean flow velocity. We believe that this approach is heuristic and a more fundamental description of dispersion under large concentration gradients has to be sought. In this work, the experiments of Hassanizadeh et al.” are briefly described, a non-linear dispersion theory is presented and theoretical and experimental evidences supporting the theory are discussed. Note that this study does not concern unstable displacements. Thus, situations where a heavy fluid overlies a lighter fluid, such that fingering occurs, are not considered here.
EXPERIMENTAL METHODS The experiments were performed in a two-dimensional vertical column made of two plexiglass plates separated by l-cm thick spacers. The inner dimensions of the column were 60 cm wide, 125 cm high, and 1 cm thick (see Fig. 2). The column was packed with glass beads with diameters ranging from O-40 to 0.52mm. The packing of the column was carried out under wet conditions. That is, the column was first half filled with water and then glass beads were placed by means of a
3-wayvalve n fresh water outflow
II
infkw
0 manometers oelectrodes
f
neter
3-wayvalve
Fig. 2. Schematic representation of
the
f&iieter
1;
-“].I
experimental set-up.
High-concentration-gradient pipe at the bottom of the column. An electric vibrator was used to vibrate the column during packing. It turned out that because of the size and uniformity of grains and homogeneity of the packing, the dispersivity of the porous medium was very low. The fluid used in the experiments was prepared by dissolving pure NaCl in distilled water. The concentration of resident salt solution was 1000 kg mm3 or higher. Flow was always directed upward with the heavier solution displacing the lighter resident solution. Thus, all one-dimensional experiments may be considered as being stable displacements. Thle fluid was injected through nine equally-spaced holes at the bottom of the column. Obviously, both flow and transport near the inlet holes had a two-dimensional character. However, both measurements and simnlations have shown that when all holes were open, both flow and transport become onedimensional within a distance of 5 cm into the column. More persistent two-dim.ensional situations were created by closing some of the inlet holes. The water exited the column through nine outlet holes at the top of the column. The upper part of the column was open to the atmosphere. It must be noted that the experiments were carried out in a constant-temperature room where the variations in temperature was about f 0.5”C. Salt mass fraction of the fluid flowing through the column was measured by means of 16 pairs of electrodes embedded into the plexiglass plates at 16 different locations (see Fig. 2). At each location, two electrodes were embedded into the plexiglass plates facing each other. The electrodes consisted of thin platinum disks, with a diameter of 1 cm, mounted on a plexiglass plug screwed into the column wall. The platinum disks were in contact with glass beads and water but they did not protrude into the porous medium; thus they are not expected to have perturbed flow and transport. An alternating current of constant amplitude was supplied to each electrode pair a:nd at the same time the voltage difference between the two electrodes was measured. The voltage difference is a function of the porous medium type and porosity, water composition, temperature, and the measuring circuit. Each and every pair of electrodes was calibrated in situ. With all other factors kept constant, the measured voltage difference could be directly related to the salt mass fraction of the solution. Individual calibration curves were prepared for each electrode pair in the range of 0.001 to 0.24 kg kg-‘. Simultaneous and continuous reading and registration of all the electrodes was made possible by a computerized data acquisition system. Note that because of
dispersion in porous media
205
significant variations in mass density, we worked with mass fraction, instead of mass concentration, in all experimental and modelling activities. The total error in the measurement of mass fraction has been estimated” to be in the range -5.0 to + 3.55%. The mass fraction of the fluid entering and leaving the column was measured by means of a conductivity meter. This meter was calibrated, in the same room where experiments were carried out, for the salt solution used in the experiments. The volumetric rate of flow of the fluid was measured before entering the column by means of an in-line flow meter. The flow meter was also calibrated for the range of mass fractions indicated above. Two sets of manometers along the two vertical sides of the column were used to measure the water head. From the values of water head and flowrate, the permeability of the medium for a given salt concentration could be calculated. In all calculations the effect of salt concentration on mass density and viscosity was taken into account. The total error associated with the estimation of average permeability of the column was about f2%. A peculiar result was that the permeability appeared to be a linear function of mass fraction of the solution; decreasing with an increase in mass fraction. The permeability showed a reduction of about 15% when mass fraction of the solution increased from 0.001 to 0.24 kg kg-‘. Note that in all flow experiments for determination of permeability, about 30 in total, the mass density of the solution was constant for the duration of the experiment. An explanation for this decrease will be offered shortly. It is known that clays show an opposite effect: their permeabilities increase when fresh water is replaced by salt water. This is believed to be due to the swelling of clay as a result of ionic forces. Prior to a displacement experiment, a solution with a known density was circulated through the column for a couple of days until a uniform salt concentration distribution in the column was reached. Then the flow was stopped by closing the valves at the bottom of the column. The reservoir feeding the nine inlet holes was emptied and filled with the displacing solution (which had a mass fraction higher than the resident mass fraction). After flushing the feeding reservoir a few times and ensuring that the salt concentration of fluid in the feeding lines was equal to the concentration of fluid in the high-concentration reservoir, the valves were opened and the displacement started. In this way, the mass fraction of the solution entering the column was known and kept constant during a displacement experiment.
Table 1. Specifications of displacement experiments
Experiment code Resident mass fraction ( 1U3 kg kg-‘) Displacing mass fraction ( 10e3kg kg-‘) Flowrate range (cm3s-l)
LlDOl
LlD02
LlD03
L2DO1
HlDOl
HlD02
H2DOl
1.37 2.88 0.40
1.75 2.78 0.36-0.45
2.97 3.70 0.20-0.23
3.84 4.50 O.Ol-0.05
2.78 81.03 0.04-0.40
3.17 235.60 0.07-0.60
4.65 176.66 0.17-0.23
206
A. Leijnse
S. Majid Hassanizadeh,
5
0.8
‘E 0.7 ,e g 0.6 2 ‘0 .B
0.5 0.4
f
0.3
=
0.2
1000
0
2000
3000
4000
5000
6000
7000
6(
time (s) Fig. 3.
Simulation of the LCG-experiment LlDOl, with Q~ = O-09cm and p = 0.
Two types of displacement experiments were carried out: low concentration gradient (LCG) and high concentration gradient (HCG). The differences in mass fractions of the resident and displacing fluids were in the order of 0.001 kg kg-’ for LCG experiments and up to 0.232 kg kg-’ for HCG experiments. All experiments and the range of mass fractions and flow rates are listed in Table 1. Details of the experiments are found in Hassanizadeh et al.” Typical measured breakthrough curves are shown in Fig. 3 for the central column of electrodes in the LlDOl experiment. The breakthrough curves for the left and right column of electrodes were almost identical in the one-dimensional experiments. However, as one might expect, this was not the case for the two-dimensional experiments. For example, in the
--*-.---
0
electrode 1L ----electrode 2C --electrode 3C ---
1000
2000
3000
case of H2DOl (where only the right half of the inlet holes were open), the front arrives at the right electrode earlier than at the central and left electrodes (see Fig. 4). However, it is apparent that even here, after the third row a one-dimensional situation prevails.
NON-LINEAR
THEORY
It is well known that when large concentration variations occur in groundwater flow, the fluid density will also change accordingly. In the absence of density effects, groundwater flow is independent of the solute transport. However, when large concentration variations occur, the two problems are coupled and have to
electrode 1C _-_ electrode 2R -_-_ electrode 3R
4000
DISPERSION
5000
electrode 1R electrode 3L
6000
7000
8000
9000
time (s)
Fig. 4. Breakthrough curves for the two-dimensional experiment H2DOl.
10
High-concentration-gradient
be solved simultaneously. The basic governing equations for coupled flow and transport are the equations of conservation of mass for fluid and salt. For non-reacting solutes, these equations
aw
v,,+ w
* .Vw+V.J=O
(4)
where w is the mass fraction of salt (w = C/p) and vf is the average velocity of the fluid phase. The Darcy velocity q is related to vf by q = n(v* - v”‘~)
(5)
The solid velocity vsofidis commonly either neglected or is related (through its divergence) to the rate of change of porosity. The former alternative has been selected in our simulations. The dispersive mass flux, J, is defined as the mass flux of salt relative to the mean motion of fluid J = npw(vS -v*)
(6)
where u’ denotes the average velocity of salt ‘particles’. Equations (3) and (4) am coupled through the equation of state. In our simulati~ons, the following relationship, which holds for a wide range of pressure and mass fraction values, has been employed” P = POew@*
(P - POI + 74
(7)
where pf is the fluid compressibility, equal to 4.45 x lOPi kg m-l K2, and y is the mass fraction coefficient, equal to 0.69. Equations (3)-(7) have to be supplemented with appropriate relationships for q and J. Commonly, the well-known Darcy and Fickian dispersion equations are employed, whose origins are in the equations of motion. For typical conditions encountered in porous media flow, the equations of motion for fluid and salt, are respectively given as7 Vp-pg=+f
(8)
VI;s =;:
(9)
These equations are vialid for the case of creeping isothermal flow and transport, with gravity as the only external body force acting. The term g is a force accounting for the resistance of the porous medium to the fluid flow. Similarly, +,”is a force accounting for the resistance of the porous medium to the salt dispersion; this resistance is offered1 by water molecules and solid grains. Both of these forces are expected to be a function of the velocities of water and salt components. The term fis is the relative chemkal potential of salt, which is the driving potential for the spreading of salt from high concentration areas to low concentration areas. The relative chemical potential accounts for the increase in the macroscopic free energy of the solution as a result of
dispersion in porous media
207
increase in the solute mass fraction;7 it is a function of salt mass fraction and fluid temperature and pressure iX = ixs(P,W, T)
(10)
Equation (8) indicates that the fluid tends to move from areas of high dynamic pressure (i.e. the excess pressure after accounting for the effect of gravity) to areas of low dynamic pressure. This force is opposed by the resistance of the medium to flow, 9, which is assumed to depend on q, J, and state variables. For low concentration situations, and when the flow velocity is small, one may assume a linear dependence of ;/ on q. The resulting equation is Darcy’s law. However, for high velocity flow (as a result of a large pressure gradient and/or high medium permeability), the linear dependance will not be valid any more.‘,17 Instead, a second or higher order relationship for 9 has to be employed. As a second order approximation, the Forchheimer equation is obtained9
(1+ a&l = -(k/p)
. (VP - ~g)
(11)
where p is the dynamic viscosity of the fluid, k is the permeability tensor, and a is the Forchheimer coefficient. The viscosity is also known to be a non-linear function of the salt mass fraction.” A typical relationship is p = /.&l + 1.85~ - 4.1~~ + 44.5~~)
(12)
Now, consider eqn (9). This equation suggests that the solute tends to spread (in addition to its movement due to the mean fluid flow) from areas of high chemical potential (and thus high mass fraction) to areas of low chemical potential. This force is opposed by the resistance of the medium to the solute dispersion,<, which is also assumed to depend on q, J, and state variables. In principle, < may be represented by a Taylor series in terms of J. When J is sufficiently small, such a series may be truncated. If terms higher than second order in the series are neglected, one will obtain 3 = -RS(l + /3J)J
(13)
where J denotes the magnitude of J, R” is a second-order tensor accounting for the resistance of the porous medium to the dispersion of the solute, and p is a new coefficient accounting for the non-linear effects. Both R and /3 are considered to be independent of J, but they may still depend on q. When there exist small concentration gradients, i.e. for small V& the dispersion will not be large and the second order term in eqn (13) may be neglected. Substitute the resulting relationship for < in eqn (9), and use eqn (10) to write Ofi: in terms of VW. Then, for isothermal conditions, and for the range of pressure variations encountered here, one obtains the following equation for the dispersion flux J=
-pD.Vw
(14)
S. Majid Hassanizadeh, A. Leijnse
208
where D, the dispersivity tensor, is related to R” by D
=
f-
(RS)-’
(15)
and is a function of flow velocity q. The empirical eqn (2) is commonly employed to describe this functional dependence. Note that the Fickian dispersion eqn (14) is valid even if the fluid mass density is not constant. The same cannot be said about eqn (1). Obviously, only if the fluid density is constant, the two formulations are equivalent. Now, consider situations where a large concentration gradient exists (and thus a large V$). In such cases the tendency for dispersion, and therefore the dispersive mass flux, could be very large. This also means (intuitively and also as suggested by eqn (9)), that the resistance of the medium to dispersion will grow. Under these conditions, the second order term in eqn (13) may not be neglected any more. Substitution of (13) in (9) yields a non-linear dispersion equation (1 + @J)J = -pD . VW
(16)
where ,0 is a new dispersion coefficient and D is still given by eqn (2) with longitudinal and transversal dispersivities considered to be (constant) properties of the porous medium and independent of the fluid properties and transport processes. It is instructive to note that eqn (16) is identical in form to the Forchheimer equation, (1 l), applicable to high velocity flow. In fact the two processes are quite similar; the coefficient ,B is the counterpart of the Forchheimer coefficient a. In principle, similar to D, the coefficient p may still depend on the flow velocity q. In the next section, both eqns (14) and (16) are employed to model the displacement experiments reported here.
SIMULATION OF EXPERIMENTS A numerical model has been constructed to solve the fluid and salt mass balance equations (3) and (4) supplemented with state equations (7) and (12), Darcy’s equation (1 l), and an appropriate equation for J. Equations (16) and (2) can be combined. In a onedimensional situation, with the z-coordinate taken to be positive upwards, the following non-Fickian equation is obtained (1 +
PJ)J = -pa&
(17)
Here, the molecular diffusion is neglected and both q and J are assumed to be positive upwards. Initial and boundary conditions are needed to solve the flow and transport equations. Initially, the solution is at rest and has a uniform salt mass fraction. Thus at t = 0, q = 0 and w = wi for all z
(lga)
where wi and pi are the mass fraction and mass density of the resident solution. For boundary conditions at the lower end of the column, because both the flow rate and the mass fraction are known, flux-type conditions are chosen at z = 0, pq = pa&) for all t>O
and J + pwq = poweq&) Wb)
where q. is the measured inflow rate, wg is the salt mass fraction of the displacing solution, and p. = p(wo) is given by eqn (7). At the outlet, the column is open to the atmosphere (zero gauge pressure). Moreover, as long as the displacement front has not reached the end of the column, there will be no dispersion. Thus, the following boundary conditions for the upper end are chosen atz=L,p=OandJ=Oforallt>O
(lgc)
It is important to note that because the flow is onedimensional and the measured flow rate is specified as a boundary condition, there would be no coupling between Darcy and Fick equations. That is, in onedimensional simulations, the flow velocity is completely determined by the boundary conditions and, therefore, the specific form of Darcy’s law will have no influence on the calculated salt mass fraction distribution. Therefore, any deviation between measured and calculated results is bound to stem from the deficiencies in either the value of model parameters or the dispersion law being used. Among the model parameters, wi, wo, and q. are measured with relatively high accuracy. Thus, for the classical Fickian formulation, only two parameters, namely porosity, n, and longitudinal dispersivity, Q~ need to be determined. These are evaluated on the basis of LCG-experiments. Two sets of simulations have been performed. In the first set, the classical Fickian relation is employed, i.e. the value of /I is set to zero. First, with the inflow rate known and almost constant for the LCG-experiment LlDOl (see Table l), porosity is determined from the arrival time of the displacement front at different rows of electrodes. It appears that porosity is slightly variable in both longitudinal and lateral directions and has an average value of about 0.4. For dispersivity, an average value of 0.09cm has been obtained by fitting the simulation results to experimental data (see Fig. 3). This value is about two times the mean grain size of the medium, as one would expect for a homogeneous porous medium (cf. Scheidegger,17 p. 306). Using these values of porosity and dispersivity, other low- and highconcentration-gradient experiments have been simulated. A satisfactory agreement between calculated and measured breakthrough curves is obtained for other LCG experiments. Results of simulation of the experiment LlD02 are given in Fig. 5. The slight deviation in the arrival times between measured and calculated breakthrough curves is probably due to errors in flow measurements, or a result of inaccuracies in the
High-concentration-gradient
. -
dispersion in porous media
209
measured calculated
1 0.9
0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0
I do0
2000
3000
4000
5000
Sob0
7600
8( 30
time (s)
Fig. 5. Simulation of the LCG-experiment LlD02, with (Ye= 0.09cm and p = 0. evaluation of porosity. The agreement in shape of calculated and measured curves, however, is quite satisfactory. This confirms the applicability of the Fickian dispersion equation at low-concentration-gradient situations. Keeping the value of porosity and dispersivity unchanged, the classical set of equations are also used to simulate the HCG-experiments. As is evident in Figs 6 and 7, results of simulations for HCGexperiments deviate significantly from the measured breakthrough curves. Interestingly enough, it is possible to obtain a good fit for IXG experiments if the value of dispersivity is lowered by a factor of three (oL= 0.03 cm). One ma:y be tempted, in an attempt to salvage the Fickian relationship, to consider the dispersivity to be a function of mass fraction or mass
fraction gradient. This possibility has been investigated and the answer is found to be negative, as discussed shortly. Instead, following the line of theoretical development presented here, the non-linear formulation of dispersion is considered as the correct alternative. The non-linear formulation is believed to be applicable if the magnitude of dispersive mass flux is high. In order to see if there is actually a difference in the magnitude of dispersive mass flux for different experiments, the maximum value of J in the system (which occurs at the front) is calculated and plotted as a function of time (Fig. 8). It is evident from Fig. 8 that the dispersive mass flux for HCG-experiments is about two orders of magnitude larger than that for LCG-experiments. Therefore, the non-linear term in (17) is bound to be important.
)-ziEiq . 0.9 5 .2
0.8 o.7
: 3 E -0 .g
0.6
E
0.3
:
. .. ...
. . .
??
: ::
0.5 0.4
0.1 0.2 0
I
I 0
2000
I 4000
6000
8000
10000
12000
14000
If 100
time (s) Fig. 6.
Simulation of the HCG-experiment HlDOl, with CQ= 0.09cm and p = 0.
S. Majid Hassanizadeh, A. Leijnse
210
1 0.9
5
0.8
‘S 0 E
0.7
E
;';
n .ij 0.4 a E =
Om3 0.2 0.1 0
0
2000
4000
6000
6000
12000
10000
14000
lf 100
time (s)
Fig. 7. Simulation of the HCG-experiment HlD02, with (Ye = 0.09cm and p = 0.
Thus, in the second set of simulations, the non-linear dispersion equation (17) is employed. Using the values of porosity and dispersivity obtained from the LIDOlexperiment, data from H 1DO2 are employed to determine the coefficient of the non-linear term. A very good agreement between measured and calculated breakthrough curves has been obtained with a value of /3 = 4.0 x 104m2 s kg-’ (Fig. 9). Next, keeping the values of n, cxL, and ,f3 unchanged, the other HCG-experiment, HlDOl (see Fig. 10) as well as all LCG-experiments were simulated and once again very good agreement was obtained. As expected, the addition of the non-linear term does not affect the results of LCG simulations. No attempts were made to simulate the two-dimensional experiments. The main reason has been the expense of computer simulations. Because of the very low
dispersivity and relatively large velocities, a highly refined mesh was necessary to minimize numerical dispersion. Given the size of the domain, this would have resulted in a large mesh and a high demand on computer resources. The data sets, however, are available for release to interested modellers.
SIGNIFICANCE OF THE NON-LINEAR TERM As mentioned above, when employing the standard linear theory, it is possible to simulate an HCGexperiment when employing an apparent dispersivity with a value lower than the dispersivity obtained for LCGexperiments. A relation for an apparent dispersivity, ay,rP, can be obtained from the non-linear dispersion theory if
.. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
.
.
. . .
.
iii
Et
lo98 0
"'I'n'n 1000
2000
3000
4000 time (s)
"'I 5000
6000
7000
8( 30
Fig. 8. Plot of the maximum value of the magnitude of dispersive mass flux, J, as a function of time for different experiments.
High-concentration-gradientdispersion in porous media
Jciiiiq F .
.
.
*
211
.
*
0!3 0.8 0.'7 0.16 0.5 0.4 0.3 0.2 0.1 O,"',"',"'~"','.',.",'.',"' 0 2000 4000 6000
8000
10000
12000
14000
16 100
time (s) Fii. 9. Simulation of the HCG-experiment HlD02, with CQ= 0.09cm and fl= 4.0 x 104m2skg-‘. we set, for a one-dimensional
situation
J = paap&
(19)
Substituting this relation in (17) and solving one obtains the following important relation. aaPP
2
-OL = I+ dmaLq(Aw/Az)
foraapp, (20)
In Fig. 11, a plot of the relative apparent dispersivity as a function of the apparent dispersive mass flux pqaL(Aw/Az) and the coefficient of non-linearity ,8, is provided. Indicated in this figure are the ranges of the apparent dispersive mass flux for field situations and for
our HCG-experiments. Because no information is available on the value of p other than for the experiments reported here, plots have been provided for a wide range of values. It is evident that for certain field situations, the non-linear term can become important, depending on the value of the non-linear coefficient. For example, for the aquifers overlying the Gorleben salt dome, assuming a dispersivity of lOm, a velocity of 1 m year-‘, and a mean mass fraction gradient of 0*25/30 m-’ (see Fig. l), the value of mean apparent dispersive mass flux is estimated to be 3 x 10e6 kgmP2 s-l. From Fig. 11, it is clear that even for such a low value of the apparent dispersive mass flux, the non-linear terms could be important. For example, for a value of ,0 = 10’m2 s kg-‘, the apparent dispersivity
time (s)
Fig. 10. Simulation of the HCG-experiment HlDOl, with CQ= 0.09 cm and p = 4.0 x 104m2skg-‘.
S. Majid Hassanizadeh, A. Leijnse
212
1 0.9
0.8 0.7
d
B
0.6 0.5
aw 0.4
1o-5
IO4
apparent dispersive mass flux pqg(Aw/Az) (kg/m2 s)
Fig. 11. Plot of the normalized apparent dispersivity as a function of the apparent dispersive mass flux and the coefficient of nonlinearity fi.
would be approximately 20% smaller than the actual dispersivity. In eqn (20), one can identify the dimensionless group of parameters &, defined by (21)
PN = %w(WAz)
This dimensionless number is a measure of the significance of non-linear effects. It is evident that for large mass fraction gradients, for high-velocity flow, or for media with a large dispersivity, the non-linear effects may be important. For example, for the low concentrationgradient experiment LlDOl, a value of 10m2m-l for Aw/Az can be estimated based on the breakthrough curve at electrode 4C. With p = 4 x lo4 m2 s kg-‘, the value of PN is approximately 0.1, which has a negligible effect on the apparent dispersivity: CY,,,~= 0.98~~~.Note that the effects will be even smaller for the other LCGexperiments LlD02 and LlD03, for which a lower value of PN is calculated. For the HCG-experiment HlD02, based on the slope of the breakthrough curves, the value of Aw/Az is estimated to vary between 4.5 and 3.0m-’ during the experiment. Thus, the associated value of ,/?N varies between 44 and 16 which corresponds to a value of the apparent dispersivity, oar,,,, varying between 0.26~~~ and 0.39~~. These values agree very well with the fact, as reported in the previous section, that the HCG-experiments can be simulated with a dispersivity of 0.03 cm.
DISCUSSION
OF RESULTS
Results obtained here clearly demonstrate that the Fickian dispersion equation is not valid under highThe non-linear situations. concentration-gradient dispersion eqn (16), however, seems to provide a
satisfactory description of dispersion under these conditions. It should be noted that a number of heuristic modifications to the linear Fickian equation were also attempted (see Hassanizadeh’ for details of the models). For example, the longitudinal dispersivity was specified to be a function of mass fraction (linearly or exponentially) or a function of mass fraction gradient (linearly, second order, or exponentially). Also an extended (yet linear) version of Fick’s law proposed by Hassanizadeh7 and Hassanizadeh and Leijnse” was attempted. But, an agreement between simulated and measured breakthrough curves could not be obtained in any of the cases. The fact that according to eqn (20) and the graph in Fig. 11, the ratio o,nP /oL will be equal to or smaller than one, is in accordance with many experimental observations reported in the literature.2’12-14,20 Krupp and Elrick14 and Slobod and Howlett2’ have performed both stable and unstable displacement experiments. They have found that in the case of stable experiments, the mixing zone (or the dispersion coefficient) becomes smaller as the difference in density of the fluid on the two sides of the front increases. What is interesting in these studies is that, in contrast to our experiments, the flow direction is vertically downward. That is to say that the fresh water displaces salt water. Thus, the non-linear effects described here exist regardless of the direction of the displacement (in stable cases). Kempers’2T’3 has also performed stable experiments and has concluded that dispersion is suppressed as a result of a density difference. Kempers12 provides a graph, obtained from numerical simulations, which gives the ratio of apparent dispersivity to the medium dispersivity as a function of the difference in fluid densities, the ratio of fluid viscosities at the two sides of
High-concentration-gradient
the front, medium properties, and flow velocity. Note that this approach does not constitute a dependence of dispersivity on mass fraction or flow velocity. The value of apparent dispersivity lis determined a priori based on fluid density differences, viscosity ratios at the two sides of the front, and flow velocity. Then, the value of dispersivity is kept const,ant for a given simulation even though both mass fraction and velocity vary in time and space. Also note that determining a single value of apparent dispersivity may work in one-dimensional systems where only one front is present. However, if two fronts with different fluid densities and/or viscosity contrasts are present, it will not be possible to define an apparent dispersivity for the whole system uniquely. This is particularly important in two- or three-dimensional situations where the existence of two or more fronts is even more plausible. In such situations the role of transversal dispersivity becomes important too. Spi@ has studied transversal dispersion across a horizontal boundary between two layers of horizontally moving fluids with different densities and with the lighter fluid overlying the heavier one. These studies show that transverse dispersivity also decreases with increase in concentration difference across a front. The non-linear dispersion eqn (16) avoids the difficulties and empiricism associated with estimating an apparent dispersivity. This equation can be used in a threedimensional flow regime and it does not rely on the existence of a single front in the system. Note that the decrease of the width of the dispersion zone (i.e. smaller apparent dispersivity) in high-concentration-gradient experiments does not mean that dispersive mass flux will be smaller compared to situations with a low concentration gradient. On the contrary, due to much higher mass fraction gradients, dispersive mass flux is high and, as a result, the resistance of the system to dispersion will be also high. In other words, the ease with which such a large amount of mass can be dispersed is decreased (lower apparent dispersivity). This explains the need for keeping the second order term in the constitutive eqn (13). It is worth noting that a numerical experiment was also performed in which a high-concentration solution was displaced by a slightly higher concentration solution and the dipsersive mass flux was calculated. It turned out that the dispersive mass flux has the same order of magnitude as that for LCG-experiments. This indicates that the gradient of concentration (or mass fraction), and not the value of concentration, causes deviations from the linear theory. It is known that dispersion is a scale-dependent phenomenon and its magnitude depends on the scale of our models. The non-linear effects described here are expected to be important at the first macroscopic scale of a porous medium; that is, one scale above the pore scale, also called the D,arcy scale. Thus, in the cases where additional mixing as a result of instabilities or
dispersion in porous media
213
heterogeneities occur, the high-concentration-gradient effects may or may not be important depending on the modelling scale. For example, in situations where a lighter fluid underlies a heavier fluid, instabilities may develop and fingering may occur. In those cases, two approaches are possible. One may choose to include the details of fingering and model the system as two dimensional. In this approach, dispersion is modelled at the Darcy scale and the non-linear effects may prove to be important. Alternatively, one may choose to model the system at a higher scale. Thus, one may lump the effects of fingering into an effective dispersion term and model the system as one-dimensional. In such an approach, basically an additional averaging is carried out and the dispersion is modelled at the scale of instabilities. Thus, the effective dispersion will be large and the non-linear effects discussed in this manuscript, will probably be negligible. We are confident that the deviations from Fick’s law observed in our study are not due to the presence of heterogeneities. This can be shown by using the results of Welty and Gelhar22 who have studied the effects of density and viscosity variations on macro-dispersivity in heterogeneous media. They provide an analytical expression for longitudinal macro-dispersivity, as a function of mean displacement distance, mean concentration and concentration gradient, density and viscosity coefficients, mean velocity, gravity, and correlation scale and variance of the log permeability process. They too find that for stable cases, the value of macrodispersivity decreases as the magnitude of concentration gradient increases. The main dispersion mechanism in the system they are considering is due to the spatial variability in flow velocity caused by local heterogeneities in permeability. The system we are considering, however, is relatively homogeneous and the main dispersion mechanism is due to variations in interstitial flow velocity and concentration gradients. Given the very low value of longitudinal dispersivity we have calculated, for the formula derived by Welty and Gelhar” to be applicable, the porous medium in our study should have a correlation length as small as O-1cm and the permeability must vary by a factor of about 6. These values are highly irrelevant to the homogeneous system in our study whose permeability has a very large correlation length and a small variance. Thus, deviations observed in our studies cannot be due to local heterogeneities, but are caused by microscale velocity variations under the influence of density effects. In this study, only the applicability of the Fickian equation has been explored and no comment can be made on the validity of Darcy’s law. For investigating that question, it would have been necessary to have pressure measurements inside the column and compare them with the calculated pressure distribution; these were unfortunately unavailable for HCG experiments because the manometer readings were useless for displacement
214
S. Majid Hassanizadeh, A. L-e&se
experiments, certainly for high-concentration-gradient experiments. The reason is that during displacement experiments, the density of the water column in manometer tubes varies and, thus, the water pressure cannot be linearly related to the manometer level. However, in a study by Schincariol and Schwartz,” it has been observed that a plume of dense salt solution moving horizontally in an ambient fresh water has a velocity much higher than that of a tracer injected with the salt solution. They have also observed brine movement occurring against the prevailing hydraulic gradient. These observations suggest that the pressure gradient may not be the only important term in Darcy’s law, and a term accounting for the large concentration gradients as proposed by Hassanizadeh7 may also be needed. A peculiar result obtained in our experiments has been the dependence of permeability on the salt mass fraction. This can t;e explained by the fact that the resistivity of the porous medium is actually a function of the fluid as well as the matrix properties. Obviously, the resistance of the medium to the flow of brine is larger than the resistance to the flow of fresh water. Traditionally, the resistivity is assumed to be given by p/k with the additional assumption that k is independent of fluid properties. Although viscosity increases with an increase in salt mass fraction, our results show that the increase of resistance to flow is even larger. Thus, the assumption that k is independent of fluid properties is not necessarily valid for a high salt concentration solution. Additional laboratory and field experiments are required to further investigate the issues and questions raised here. In particular, applicability of the non-linear dispersion equation under more general conditions and validity of Darcy’s equation under large concentration gradient conditions have to be studied. The significance of non-linear effects in heterogeneous media and under field conditions needs to be investigated. Additional experiments are also required to determine the dependence of /3 on flow velocity and the possible dependence of permeability on solute concentration.
CONCLUSIONS The Fickian description of dispersion in a porous medium is based on the assumption that the dispersive mass flux is small such that the second order term in the constitutive relation (13) is negligible. Experiments have shown that this assumption is valid only if small concentration gradients are present in the system. For large concentration gradients, large amounts of mass tend to be dispersed and the above-mentioned assumption breaks down. As a result, the breakthrough curves calculated based on the Fickian-type dispersion eqn (14) deviate substantially from the measured curves. It appears that under high-concentration-gradient condi-
tions, the dispersive mass flux will be large so that the non-linear dispersion eqn (17) has to be employed. Both low- and high-concentration-gradient experiments are simulated successfully with this equation. The model parameters aL. and p appear to be independent of the fluid properties. The coefficient p could, however, be a function of the flow velocity.
ACKNOWLEDGMENTS The experiments reported here have been originally designed and set up in the Geotechnique Laboratory of the Delft University of Technology, The Netherlands. The expertise of Ab Mensinga and Joop van Leeuwen of this laboratory has been instrumental in the success of these experiments. Discussions with Jan van Eijkeren of Center for Mathematical Methods of RIVh4 have been very useful in developing the non-linear theory. Wilco de Vries of RIVM has assisted with the calibration of electrodes and the flowmeter. Their contributions are gratefully appreciated. The final version of this paper was written when the senior author was visiting the Department of Civil Engineering and Geological Sciences of the University of Notre Dame. Financial support for this visit was generously made available by RIVM and the University of Notre Dame.
REFERENCES 1. Bear, J., Dynamics of Flui& in Porous Media. Elsevier, New York, USA, 1972. 2. Ben Salah, M. D., Influence des contrastes de viscositk et de densiti sur le dkplacement en milieu poreux de deux fluides miscibles. Revue de l‘lnstitut Frangais du Petrole, (1965) 1237-55. 3. Bues, M. A., Aachib, M. and Zilliox, L., Incidence of heterogeneities on pollutant transport-density and viscosity contrasts of the liquid phase-structure of solid matrix. In Contaminant Transport in Groundwater, 4. H.E. Kobus & W. Kinzelbach, 1990, pp. 251-7. 4. Bues, M. A. & Zilliox, L., D&placement miscible avec contrastes de densitC et de viscositk en milieu poreux. Identification des parametres de d&placementet stabi& en regime de dispersion mkanique. J. Hydrol., 120 (1990) 125-41. 5. Brigham, W. E., Reed, P. W. & Dew, J. N., Experiments on mixing during miscible displacements in porous media. Sot. Petrol. Engng J., March (1961) 1-8. 6. Hassanizadeh, S. M., Derivation of basic equations of mass transport in porous media. Part 1. Macroscopic balance laws. Adv. Water Resour., 9 (1986) 196-206. 7. Hassanizadeh, S. M., Derivation of basic equations of mass transport in porous media. Part 2. Generalized Darcy’s and Fick’s laws. Adv. Water Resour., 9 (1986) 207-22. 8. Hassanizadeh, S. M., Experimental study of coupled flow and mass transport: A model validation exercise. In ModelCare Modelling,
90: Calibration and Reliability
in Groundwater
IAHS Publ. No. 195, 1990, pp. 241-50. 9. Hassanizadeh, S. M. & Gray, W. G., High velocity flow in porous media. Transport in Porous Media, 2 (1987) 52 1-3 1.
High-concentration-gradient
dispersion in porous media
215
10. Hassanizadeh, S. M. 81 Leijnse, A., On the modeling of brine transport in porous media. Water Resow. Res., 24 (1988) 321-30. 11. Hassanizadeh, S. M., Leijnse, A., De Vries, W. J. & Stapper, R. A. M., Experimental Study of Brine Transport in Porous Media. Report 728514005, RIVM, Bilthoven, The Netherlands, 1990. 12. Kempers, L. J. T. M., Effects of Fluid Properties on
17. Scheidegger, A. E., The Physics of Flow through Porous Media. 3rd edn, Univ. of Toronto Press, Toronto, Canada, 1974. 18. Schelkes, K. & Vogel, P., Paleohydrological information as an important tool for groundwater modeling of Gorleben site. In Paleohydrological Methods and Their
Hydrodynamic Dispersion: Comparison of Analytical Models to Numerical Simulations. Proc. 1st SPE/IMA Europ. Conf Mathematics of Oil Recovery. Cambridge,
19. Schincariol, R. A. & Schwartz, F. W., An experimental investigation of variable density flow and mixing in homogeneous and heterogeneous media. Water Resour. Res., 26 (1990) 2317-29. 20. Slobod, R. L. & Howlett, W. E., Effects of gravity segregation in laboratory studies of miscible displacement in vertical unconsolidated porous media. Sot. Petrol. Engng J., March (1964) 1-8. 21. Spitz, K., Dispersion in Porosen Medien: Einth_$ von Inhomogeniteiten und Dichteunterschieden. PhD Thesis, University of Stuttgart, Germany, 1985. 22. Welty, C. & Gelhar, L. W., Stochastic analysis of the effects of fluid density and viscosity variability on macrodispersion in heterogeneous porous media. Water Resour. Res., 27 (1991) 2061-75.
UK, 25-27 July 1989. 13. Kempers, L. J. T. M., Dispersive Mixing in Stable and Unstable Miscible Displacements. PhD Thesis, Delft University of Technology, The Netherlands, 199 1. 14. Krupp, H. K. & Elrick, D. E., Density effects in miscible displacement experiments. Soil Sci., 107 (1969) 372-80. 15. Matheron, G. & De Marsily, G., Is transport in porous media always diffusive? A counter example. Water Resour. Res., 16 (1980) 901-17. 16. Newberg, M. A. & Foh, S. E., Measurement of Longitudinal Dispersion Coefficients of Gas Flowing Through Porous Media. SPE Paper 17731, 1988.
Applications for Radioactive Waste Disposal. Proc. OECD/ NEA Wkshop. OECD, Paris, France, 1992, pp. 237-50.
ELSEVIER
Advances in Water Resources, Vol. 18, No. 4, pp. 203-215, 1995 Copyright 0 1995 Elsevier Science Limited printed in Great Britain. All rights reserved 0309~1708/95$09.50+0.00
A non-linear theory of high-concentration-gradient dispersion in porous media S. Majid Hassanizadeh & Anton Leijnse National Institute of Public Health and Environmental Protection (RIVM), PO Box 1, 3720 BA Bilthoven, The Netherlands (Received 23 September 1994; accepted 15 March 1995)
Tie application
of Fick’s law to describe hydrodynamic dispersion in porous media is based on the assumption of a linear dependance of a solute dispersive m:ass flux on its concentration gradient. Both theoretical and experimental studies ha.ve shown that the Fickian description of dispersion is not valid when large concentration variations in the porous medium are encountered. However, an appropriate alternative is still lacking. In this work, based on a theoretical derivation of the Fickian dispersion equation, a non-linear theory of dispersion is suggested. In the non-linear theory, in addition to the longitudinal and transversal dispersivities, a new parameter is introduced. Miscible displacement experiments are carried out in order to investigate the effects of large variations in salt mass fraction and to assess the validity of the new theory. Low-concentration liquid is displaced upwards in a vertical column by a high-concentration liquid. Thus, only hydrodynamically stable flow regimes are considered. The experiments are simulated by means of both classical Fick’s law and the new non-linear theory. It is found that low-concentration-gradient experiments can be simulated satisfactorily using the Fickian-type dispersion equation. However, calculated breakthrough curves for high-concentration-gradient experiments deviate substantially from the measured curves. It appears that a satisfactory fit to high-concentration-gradient data can be obtained only if the value of longitudinal dispersivity is reduced by a factor of three. Using the nonlinear theory, however, it is possible to simulate both low- and highconcentration-gradient experiments with a unique set of parameter values. Key wor& dispersion, solute transport, Fick’s law, Darcy’s law, brine transport, salt transport.
and I is the unit tensor. Although aL and or are usually assumed to be independent of the fluid properties, it is known that they are time and scale dependent.” The Fickian-type eqn (l), according to which the dispersive mass flux is linear in concentration gradient and independent of the concentration itself, has been successfully used in describing the dispersion of lowconcentration solutes in porous media. Indeed in many cases of ground water pollution, solute transport in unsaturated media, or tracer experiments, the concentration of dissolved matter is low so that density effects may be neglected. However, there exist certain situations of practical interest where high salt concentrations in groundwater are encountered. An important example is the presence of concentrated brines in deep aquifers overlying salt formations. The problem is of utmost interest in the safety assessment studies of high-level radioactive waste disposal in salt formations. In Fig. 1
INTRODUCTION
Hydrodynamic dispersion of solutes in porous media is normally described by a Fickian type equation with a dispersion tensor replacing the diffusion coefficient’ J=-D.VC
(1)
where J is the dispersive mass flux vector, D is the dispersion tensor, and C is the mass concentration. The dispersion tensor is assumed to be independent of the solute mass fraction and is given by17 D = Wo + oz.q)I + (QL - oT)qqlq
(2)
where 12 is porosity, DO is the effective molecular diffusion coefficient, (Ye and or are the longitudinal and transversal dispersivities, q is the apparent mean velocity of fluid with respect to the rock (also called Darcian velocity or specific discharge), q is its magnitude 203
204
S. Majid Hmsanizadeh,
50
1
0-I ~.
.I’.
" 2 cd
-50. ., -i.
-250
-
1000
1200
1100 density (kg/m3)
Fig. 1. Distribution of groundwater mass density in strata overlying the Gorleben salt dome in Germany (courtesy of Bundesanstalt fur Geowissenschaften und Rohstoffe, Hannover, Germany). the mass density of groundwater in the strata overlying the Gorleben salt dome in Germany is given.18 It is evident that the mass density increases drastically (due to an increase in salt concentration), from about 1040 kgmP3 to 1180 kg mm3 over a distance of 30 m. The question then arises whether equations such as Darcy’s law and Fick’s law are still valid under such large concentration gradients. In fact, many experiments have shown that the linear relationship (1) does not hold anymore when there exist large differences in density and/or viscosity of the
A. L.eijnse
displacing and displaced fluids2-5,‘3~16920 (for a review see Kempers13). In general, it appears that dispersivity coefficients decrease as the density (or concentration) differences increase. A similar observation has been made in a recent experimental study by Hassanizadeh et al. I1 where low-concentration salt water is displaced by high-concentration brine in an upward direction. Traditionally, in an attempt to allow for these effects, most researchers still employ eqns (1) and (2) but allow for a dependence of oL and or on density differences, viscosity ratios and the mean flow velocity. We believe that this approach is heuristic and a more fundamental description of dispersion under large concentration gradients has to be sought. In this work, the experiments of Hassanizadeh et al.” are briefly described, a non-linear dispersion theory is presented and theoretical and experimental evidences supporting the theory are discussed. Note that this study does not concern unstable displacements. Thus, situations where a heavy fluid overlies a lighter fluid, such that fingering occurs, are not considered here.
EXPERIMENTAL METHODS The experiments were performed in a two-dimensional vertical column made of two plexiglass plates separated by l-cm thick spacers. The inner dimensions of the column were 60 cm wide, 125 cm high, and 1 cm thick (see Fig. 2). The column was packed with glass beads with diameters ranging from O-40 to 0.52mm. The packing of the column was carried out under wet conditions. That is, the column was first half filled with water and then glass beads were placed by means of a
3-wayvalve n fresh water outflow
II
infkw
0 manometers oelectrodes
f
neter
3-wayvalve
Fig. 2. Schematic representation of
the
f&iieter
1;
-“].I
experimental set-up.
High-concentration-gradient pipe at the bottom of the column. An electric vibrator was used to vibrate the column during packing. It turned out that because of the size and uniformity of grains and homogeneity of the packing, the dispersivity of the porous medium was very low. The fluid used in the experiments was prepared by dissolving pure NaCl in distilled water. The concentration of resident salt solution was 1000 kg mm3 or higher. Flow was always directed upward with the heavier solution displacing the lighter resident solution. Thus, all one-dimensional experiments may be considered as being stable displacements. Thle fluid was injected through nine equally-spaced holes at the bottom of the column. Obviously, both flow and transport near the inlet holes had a two-dimensional character. However, both measurements and simnlations have shown that when all holes were open, both flow and transport become onedimensional within a distance of 5 cm into the column. More persistent two-dim.ensional situations were created by closing some of the inlet holes. The water exited the column through nine outlet holes at the top of the column. The upper part of the column was open to the atmosphere. It must be noted that the experiments were carried out in a constant-temperature room where the variations in temperature was about f 0.5”C. Salt mass fraction of the fluid flowing through the column was measured by means of 16 pairs of electrodes embedded into the plexiglass plates at 16 different locations (see Fig. 2). At each location, two electrodes were embedded into the plexiglass plates facing each other. The electrodes consisted of thin platinum disks, with a diameter of 1 cm, mounted on a plexiglass plug screwed into the column wall. The platinum disks were in contact with glass beads and water but they did not protrude into the porous medium; thus they are not expected to have perturbed flow and transport. An alternating current of constant amplitude was supplied to each electrode pair a:nd at the same time the voltage difference between the two electrodes was measured. The voltage difference is a function of the porous medium type and porosity, water composition, temperature, and the measuring circuit. Each and every pair of electrodes was calibrated in situ. With all other factors kept constant, the measured voltage difference could be directly related to the salt mass fraction of the solution. Individual calibration curves were prepared for each electrode pair in the range of 0.001 to 0.24 kg kg-‘. Simultaneous and continuous reading and registration of all the electrodes was made possible by a computerized data acquisition system. Note that because of
dispersion in porous media
205
significant variations in mass density, we worked with mass fraction, instead of mass concentration, in all experimental and modelling activities. The total error in the measurement of mass fraction has been estimated” to be in the range -5.0 to + 3.55%. The mass fraction of the fluid entering and leaving the column was measured by means of a conductivity meter. This meter was calibrated, in the same room where experiments were carried out, for the salt solution used in the experiments. The volumetric rate of flow of the fluid was measured before entering the column by means of an in-line flow meter. The flow meter was also calibrated for the range of mass fractions indicated above. Two sets of manometers along the two vertical sides of the column were used to measure the water head. From the values of water head and flowrate, the permeability of the medium for a given salt concentration could be calculated. In all calculations the effect of salt concentration on mass density and viscosity was taken into account. The total error associated with the estimation of average permeability of the column was about f2%. A peculiar result was that the permeability appeared to be a linear function of mass fraction of the solution; decreasing with an increase in mass fraction. The permeability showed a reduction of about 15% when mass fraction of the solution increased from 0.001 to 0.24 kg kg-‘. Note that in all flow experiments for determination of permeability, about 30 in total, the mass density of the solution was constant for the duration of the experiment. An explanation for this decrease will be offered shortly. It is known that clays show an opposite effect: their permeabilities increase when fresh water is replaced by salt water. This is believed to be due to the swelling of clay as a result of ionic forces. Prior to a displacement experiment, a solution with a known density was circulated through the column for a couple of days until a uniform salt concentration distribution in the column was reached. Then the flow was stopped by closing the valves at the bottom of the column. The reservoir feeding the nine inlet holes was emptied and filled with the displacing solution (which had a mass fraction higher than the resident mass fraction). After flushing the feeding reservoir a few times and ensuring that the salt concentration of fluid in the feeding lines was equal to the concentration of fluid in the high-concentration reservoir, the valves were opened and the displacement started. In this way, the mass fraction of the solution entering the column was known and kept constant during a displacement experiment.
Table 1. Specifications of displacement experiments
Experiment code Resident mass fraction ( 1U3 kg kg-‘) Displacing mass fraction ( 10e3kg kg-‘) Flowrate range (cm3s-l)
LlDOl
LlD02
LlD03
L2DO1
HlDOl
HlD02
H2DOl
1.37 2.88 0.40
1.75 2.78 0.36-0.45
2.97 3.70 0.20-0.23
3.84 4.50 O.Ol-0.05
2.78 81.03 0.04-0.40
3.17 235.60 0.07-0.60
4.65 176.66 0.17-0.23
206
A. Leijnse
S. Majid Hassanizadeh,
5
0.8
‘E 0.7 ,e g 0.6 2 ‘0 .B
0.5 0.4
f
0.3
=
0.2
1000
0
2000
3000
4000
5000
6000
7000
6(
time (s) Fig. 3.
Simulation of the LCG-experiment LlDOl, with Q~ = O-09cm and p = 0.
Two types of displacement experiments were carried out: low concentration gradient (LCG) and high concentration gradient (HCG). The differences in mass fractions of the resident and displacing fluids were in the order of 0.001 kg kg-’ for LCG experiments and up to 0.232 kg kg-’ for HCG experiments. All experiments and the range of mass fractions and flow rates are listed in Table 1. Details of the experiments are found in Hassanizadeh et al.” Typical measured breakthrough curves are shown in Fig. 3 for the central column of electrodes in the LlDOl experiment. The breakthrough curves for the left and right column of electrodes were almost identical in the one-dimensional experiments. However, as one might expect, this was not the case for the two-dimensional experiments. For example, in the
--*-.---
0
electrode 1L ----electrode 2C --electrode 3C ---
1000
2000
3000
case of H2DOl (where only the right half of the inlet holes were open), the front arrives at the right electrode earlier than at the central and left electrodes (see Fig. 4). However, it is apparent that even here, after the third row a one-dimensional situation prevails.
NON-LINEAR
THEORY
It is well known that when large concentration variations occur in groundwater flow, the fluid density will also change accordingly. In the absence of density effects, groundwater flow is independent of the solute transport. However, when large concentration variations occur, the two problems are coupled and have to
electrode 1C _-_ electrode 2R -_-_ electrode 3R
4000
DISPERSION
5000
electrode 1R electrode 3L
6000
7000
8000
9000
time (s)
Fig. 4. Breakthrough curves for the two-dimensional experiment H2DOl.
10
High-concentration-gradient
be solved simultaneously. The basic governing equations for coupled flow and transport are the equations of conservation of mass for fluid and salt. For non-reacting solutes, these equations
aw
v,,+ w
* .Vw+V.J=O
(4)
where w is the mass fraction of salt (w = C/p) and vf is the average velocity of the fluid phase. The Darcy velocity q is related to vf by q = n(v* - v”‘~)
(5)
The solid velocity vsofidis commonly either neglected or is related (through its divergence) to the rate of change of porosity. The former alternative has been selected in our simulations. The dispersive mass flux, J, is defined as the mass flux of salt relative to the mean motion of fluid J = npw(vS -v*)
(6)
where u’ denotes the average velocity of salt ‘particles’. Equations (3) and (4) am coupled through the equation of state. In our simulati~ons, the following relationship, which holds for a wide range of pressure and mass fraction values, has been employed” P = POew@*
(P - POI + 74
(7)
where pf is the fluid compressibility, equal to 4.45 x lOPi kg m-l K2, and y is the mass fraction coefficient, equal to 0.69. Equations (3)-(7) have to be supplemented with appropriate relationships for q and J. Commonly, the well-known Darcy and Fickian dispersion equations are employed, whose origins are in the equations of motion. For typical conditions encountered in porous media flow, the equations of motion for fluid and salt, are respectively given as7 Vp-pg=+f
(8)
VI;s =;:
(9)
These equations are vialid for the case of creeping isothermal flow and transport, with gravity as the only external body force acting. The term g is a force accounting for the resistance of the porous medium to the fluid flow. Similarly, +,”is a force accounting for the resistance of the porous medium to the salt dispersion; this resistance is offered1 by water molecules and solid grains. Both of these forces are expected to be a function of the velocities of water and salt components. The term fis is the relative chemkal potential of salt, which is the driving potential for the spreading of salt from high concentration areas to low concentration areas. The relative chemical potential accounts for the increase in the macroscopic free energy of the solution as a result of
dispersion in porous media
207
increase in the solute mass fraction;7 it is a function of salt mass fraction and fluid temperature and pressure iX = ixs(P,W, T)
(10)
Equation (8) indicates that the fluid tends to move from areas of high dynamic pressure (i.e. the excess pressure after accounting for the effect of gravity) to areas of low dynamic pressure. This force is opposed by the resistance of the medium to flow, 9, which is assumed to depend on q, J, and state variables. For low concentration situations, and when the flow velocity is small, one may assume a linear dependence of ;/ on q. The resulting equation is Darcy’s law. However, for high velocity flow (as a result of a large pressure gradient and/or high medium permeability), the linear dependance will not be valid any more.‘,17 Instead, a second or higher order relationship for 9 has to be employed. As a second order approximation, the Forchheimer equation is obtained9
(1+ a&l = -(k/p)
. (VP - ~g)
(11)
where p is the dynamic viscosity of the fluid, k is the permeability tensor, and a is the Forchheimer coefficient. The viscosity is also known to be a non-linear function of the salt mass fraction.” A typical relationship is p = /.&l + 1.85~ - 4.1~~ + 44.5~~)
(12)
Now, consider eqn (9). This equation suggests that the solute tends to spread (in addition to its movement due to the mean fluid flow) from areas of high chemical potential (and thus high mass fraction) to areas of low chemical potential. This force is opposed by the resistance of the medium to the solute dispersion,<, which is also assumed to depend on q, J, and state variables. In principle, < may be represented by a Taylor series in terms of J. When J is sufficiently small, such a series may be truncated. If terms higher than second order in the series are neglected, one will obtain 3 = -RS(l + /3J)J
(13)
where J denotes the magnitude of J, R” is a second-order tensor accounting for the resistance of the porous medium to the dispersion of the solute, and p is a new coefficient accounting for the non-linear effects. Both R and /3 are considered to be independent of J, but they may still depend on q. When there exist small concentration gradients, i.e. for small V& the dispersion will not be large and the second order term in eqn (13) may be neglected. Substitute the resulting relationship for < in eqn (9), and use eqn (10) to write Ofi: in terms of VW. Then, for isothermal conditions, and for the range of pressure variations encountered here, one obtains the following equation for the dispersion flux J=
-pD.Vw
(14)
S. Majid Hassanizadeh, A. Leijnse
208
where D, the dispersivity tensor, is related to R” by D
=
f-
(RS)-’
(15)
and is a function of flow velocity q. The empirical eqn (2) is commonly employed to describe this functional dependence. Note that the Fickian dispersion eqn (14) is valid even if the fluid mass density is not constant. The same cannot be said about eqn (1). Obviously, only if the fluid density is constant, the two formulations are equivalent. Now, consider situations where a large concentration gradient exists (and thus a large V$). In such cases the tendency for dispersion, and therefore the dispersive mass flux, could be very large. This also means (intuitively and also as suggested by eqn (9)), that the resistance of the medium to dispersion will grow. Under these conditions, the second order term in eqn (13) may not be neglected any more. Substitution of (13) in (9) yields a non-linear dispersion equation (1 + @J)J = -pD . VW
(16)
where ,0 is a new dispersion coefficient and D is still given by eqn (2) with longitudinal and transversal dispersivities considered to be (constant) properties of the porous medium and independent of the fluid properties and transport processes. It is instructive to note that eqn (16) is identical in form to the Forchheimer equation, (1 l), applicable to high velocity flow. In fact the two processes are quite similar; the coefficient ,B is the counterpart of the Forchheimer coefficient a. In principle, similar to D, the coefficient p may still depend on the flow velocity q. In the next section, both eqns (14) and (16) are employed to model the displacement experiments reported here.
SIMULATION OF EXPERIMENTS A numerical model has been constructed to solve the fluid and salt mass balance equations (3) and (4) supplemented with state equations (7) and (12), Darcy’s equation (1 l), and an appropriate equation for J. Equations (16) and (2) can be combined. In a onedimensional situation, with the z-coordinate taken to be positive upwards, the following non-Fickian equation is obtained (1 +
PJ)J = -pa&
(17)
Here, the molecular diffusion is neglected and both q and J are assumed to be positive upwards. Initial and boundary conditions are needed to solve the flow and transport equations. Initially, the solution is at rest and has a uniform salt mass fraction. Thus at t = 0, q = 0 and w = wi for all z
(lga)
where wi and pi are the mass fraction and mass density of the resident solution. For boundary conditions at the lower end of the column, because both the flow rate and the mass fraction are known, flux-type conditions are chosen at z = 0, pq = pa&) for all t>O
and J + pwq = poweq&) Wb)
where q. is the measured inflow rate, wg is the salt mass fraction of the displacing solution, and p. = p(wo) is given by eqn (7). At the outlet, the column is open to the atmosphere (zero gauge pressure). Moreover, as long as the displacement front has not reached the end of the column, there will be no dispersion. Thus, the following boundary conditions for the upper end are chosen atz=L,p=OandJ=Oforallt>O
(lgc)
It is important to note that because the flow is onedimensional and the measured flow rate is specified as a boundary condition, there would be no coupling between Darcy and Fick equations. That is, in onedimensional simulations, the flow velocity is completely determined by the boundary conditions and, therefore, the specific form of Darcy’s law will have no influence on the calculated salt mass fraction distribution. Therefore, any deviation between measured and calculated results is bound to stem from the deficiencies in either the value of model parameters or the dispersion law being used. Among the model parameters, wi, wo, and q. are measured with relatively high accuracy. Thus, for the classical Fickian formulation, only two parameters, namely porosity, n, and longitudinal dispersivity, Q~ need to be determined. These are evaluated on the basis of LCG-experiments. Two sets of simulations have been performed. In the first set, the classical Fickian relation is employed, i.e. the value of /I is set to zero. First, with the inflow rate known and almost constant for the LCG-experiment LlDOl (see Table l), porosity is determined from the arrival time of the displacement front at different rows of electrodes. It appears that porosity is slightly variable in both longitudinal and lateral directions and has an average value of about 0.4. For dispersivity, an average value of 0.09cm has been obtained by fitting the simulation results to experimental data (see Fig. 3). This value is about two times the mean grain size of the medium, as one would expect for a homogeneous porous medium (cf. Scheidegger,17 p. 306). Using these values of porosity and dispersivity, other low- and highconcentration-gradient experiments have been simulated. A satisfactory agreement between calculated and measured breakthrough curves is obtained for other LCG experiments. Results of simulation of the experiment LlD02 are given in Fig. 5. The slight deviation in the arrival times between measured and calculated breakthrough curves is probably due to errors in flow measurements, or a result of inaccuracies in the
High-concentration-gradient
. -
dispersion in porous media
209
measured calculated
1 0.9
0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0
I do0
2000
3000
4000
5000
Sob0
7600
8( 30
time (s)
Fig. 5. Simulation of the LCG-experiment LlD02, with (Ye= 0.09cm and p = 0. evaluation of porosity. The agreement in shape of calculated and measured curves, however, is quite satisfactory. This confirms the applicability of the Fickian dispersion equation at low-concentration-gradient situations. Keeping the value of porosity and dispersivity unchanged, the classical set of equations are also used to simulate the HCG-experiments. As is evident in Figs 6 and 7, results of simulations for HCGexperiments deviate significantly from the measured breakthrough curves. Interestingly enough, it is possible to obtain a good fit for IXG experiments if the value of dispersivity is lowered by a factor of three (oL= 0.03 cm). One ma:y be tempted, in an attempt to salvage the Fickian relationship, to consider the dispersivity to be a function of mass fraction or mass
fraction gradient. This possibility has been investigated and the answer is found to be negative, as discussed shortly. Instead, following the line of theoretical development presented here, the non-linear formulation of dispersion is considered as the correct alternative. The non-linear formulation is believed to be applicable if the magnitude of dispersive mass flux is high. In order to see if there is actually a difference in the magnitude of dispersive mass flux for different experiments, the maximum value of J in the system (which occurs at the front) is calculated and plotted as a function of time (Fig. 8). It is evident from Fig. 8 that the dispersive mass flux for HCG-experiments is about two orders of magnitude larger than that for LCG-experiments. Therefore, the non-linear term in (17) is bound to be important.
)-ziEiq . 0.9 5 .2
0.8 o.7
: 3 E -0 .g
0.6
E
0.3
:
. .. ...
. . .
??
: ::
0.5 0.4
0.1 0.2 0
I
I 0
2000
I 4000
6000
8000
10000
12000
14000
If 100
time (s) Fig. 6.
Simulation of the HCG-experiment HlDOl, with CQ= 0.09cm and p = 0.
S. Majid Hassanizadeh, A. Leijnse
210
1 0.9
5
0.8
‘S 0 E
0.7
E
;';
n .ij 0.4 a E =
Om3 0.2 0.1 0
0
2000
4000
6000
6000
12000
10000
14000
lf 100
time (s)
Fig. 7. Simulation of the HCG-experiment HlD02, with (Ye = 0.09cm and p = 0.
Thus, in the second set of simulations, the non-linear dispersion equation (17) is employed. Using the values of porosity and dispersivity obtained from the LIDOlexperiment, data from H 1DO2 are employed to determine the coefficient of the non-linear term. A very good agreement between measured and calculated breakthrough curves has been obtained with a value of /3 = 4.0 x 104m2 s kg-’ (Fig. 9). Next, keeping the values of n, cxL, and ,f3 unchanged, the other HCG-experiment, HlDOl (see Fig. 10) as well as all LCG-experiments were simulated and once again very good agreement was obtained. As expected, the addition of the non-linear term does not affect the results of LCG simulations. No attempts were made to simulate the two-dimensional experiments. The main reason has been the expense of computer simulations. Because of the very low
dispersivity and relatively large velocities, a highly refined mesh was necessary to minimize numerical dispersion. Given the size of the domain, this would have resulted in a large mesh and a high demand on computer resources. The data sets, however, are available for release to interested modellers.
SIGNIFICANCE OF THE NON-LINEAR TERM As mentioned above, when employing the standard linear theory, it is possible to simulate an HCGexperiment when employing an apparent dispersivity with a value lower than the dispersivity obtained for LCGexperiments. A relation for an apparent dispersivity, ay,rP, can be obtained from the non-linear dispersion theory if
.. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
.
.
. . .
.
iii
Et
lo98 0
"'I'n'n 1000
2000
3000
4000 time (s)
"'I 5000
6000
7000
8( 30
Fig. 8. Plot of the maximum value of the magnitude of dispersive mass flux, J, as a function of time for different experiments.
High-concentration-gradientdispersion in porous media
Jciiiiq F .
.
.
*
211
.
*
0!3 0.8 0.'7 0.16 0.5 0.4 0.3 0.2 0.1 O,"',"',"'~"','.',.",'.',"' 0 2000 4000 6000
8000
10000
12000
14000
16 100
time (s) Fii. 9. Simulation of the HCG-experiment HlD02, with CQ= 0.09cm and fl= 4.0 x 104m2skg-‘. we set, for a one-dimensional
situation
J = paap&
(19)
Substituting this relation in (17) and solving one obtains the following important relation. aaPP
2
-OL = I+ dmaLq(Aw/Az)
foraapp, (20)
In Fig. 11, a plot of the relative apparent dispersivity as a function of the apparent dispersive mass flux pqaL(Aw/Az) and the coefficient of non-linearity ,8, is provided. Indicated in this figure are the ranges of the apparent dispersive mass flux for field situations and for
our HCG-experiments. Because no information is available on the value of p other than for the experiments reported here, plots have been provided for a wide range of values. It is evident that for certain field situations, the non-linear term can become important, depending on the value of the non-linear coefficient. For example, for the aquifers overlying the Gorleben salt dome, assuming a dispersivity of lOm, a velocity of 1 m year-‘, and a mean mass fraction gradient of 0*25/30 m-’ (see Fig. l), the value of mean apparent dispersive mass flux is estimated to be 3 x 10e6 kgmP2 s-l. From Fig. 11, it is clear that even for such a low value of the apparent dispersive mass flux, the non-linear terms could be important. For example, for a value of ,0 = 10’m2 s kg-‘, the apparent dispersivity
time (s)
Fig. 10. Simulation of the HCG-experiment HlDOl, with CQ= 0.09 cm and p = 4.0 x 104m2skg-‘.
S. Majid Hassanizadeh, A. Leijnse
212
1 0.9
0.8 0.7
d
B
0.6 0.5
aw 0.4
1o-5
IO4
apparent dispersive mass flux pqg(Aw/Az) (kg/m2 s)
Fig. 11. Plot of the normalized apparent dispersivity as a function of the apparent dispersive mass flux and the coefficient of nonlinearity fi.
would be approximately 20% smaller than the actual dispersivity. In eqn (20), one can identify the dimensionless group of parameters &, defined by (21)
PN = %w(WAz)
This dimensionless number is a measure of the significance of non-linear effects. It is evident that for large mass fraction gradients, for high-velocity flow, or for media with a large dispersivity, the non-linear effects may be important. For example, for the low concentrationgradient experiment LlDOl, a value of 10m2m-l for Aw/Az can be estimated based on the breakthrough curve at electrode 4C. With p = 4 x lo4 m2 s kg-‘, the value of PN is approximately 0.1, which has a negligible effect on the apparent dispersivity: CY,,,~= 0.98~~~.Note that the effects will be even smaller for the other LCGexperiments LlD02 and LlD03, for which a lower value of PN is calculated. For the HCG-experiment HlD02, based on the slope of the breakthrough curves, the value of Aw/Az is estimated to vary between 4.5 and 3.0m-’ during the experiment. Thus, the associated value of ,/?N varies between 44 and 16 which corresponds to a value of the apparent dispersivity, oar,,,, varying between 0.26~~~ and 0.39~~. These values agree very well with the fact, as reported in the previous section, that the HCG-experiments can be simulated with a dispersivity of 0.03 cm.
DISCUSSION
OF RESULTS
Results obtained here clearly demonstrate that the Fickian dispersion equation is not valid under highThe non-linear situations. concentration-gradient dispersion eqn (16), however, seems to provide a
satisfactory description of dispersion under these conditions. It should be noted that a number of heuristic modifications to the linear Fickian equation were also attempted (see Hassanizadeh’ for details of the models). For example, the longitudinal dispersivity was specified to be a function of mass fraction (linearly or exponentially) or a function of mass fraction gradient (linearly, second order, or exponentially). Also an extended (yet linear) version of Fick’s law proposed by Hassanizadeh7 and Hassanizadeh and Leijnse” was attempted. But, an agreement between simulated and measured breakthrough curves could not be obtained in any of the cases. The fact that according to eqn (20) and the graph in Fig. 11, the ratio o,nP /oL will be equal to or smaller than one, is in accordance with many experimental observations reported in the literature.2’12-14,20 Krupp and Elrick14 and Slobod and Howlett2’ have performed both stable and unstable displacement experiments. They have found that in the case of stable experiments, the mixing zone (or the dispersion coefficient) becomes smaller as the difference in density of the fluid on the two sides of the front increases. What is interesting in these studies is that, in contrast to our experiments, the flow direction is vertically downward. That is to say that the fresh water displaces salt water. Thus, the non-linear effects described here exist regardless of the direction of the displacement (in stable cases). Kempers’2T’3 has also performed stable experiments and has concluded that dispersion is suppressed as a result of a density difference. Kempers12 provides a graph, obtained from numerical simulations, which gives the ratio of apparent dispersivity to the medium dispersivity as a function of the difference in fluid densities, the ratio of fluid viscosities at the two sides of
High-concentration-gradient
the front, medium properties, and flow velocity. Note that this approach does not constitute a dependence of dispersivity on mass fraction or flow velocity. The value of apparent dispersivity lis determined a priori based on fluid density differences, viscosity ratios at the two sides of the front, and flow velocity. Then, the value of dispersivity is kept const,ant for a given simulation even though both mass fraction and velocity vary in time and space. Also note that determining a single value of apparent dispersivity may work in one-dimensional systems where only one front is present. However, if two fronts with different fluid densities and/or viscosity contrasts are present, it will not be possible to define an apparent dispersivity for the whole system uniquely. This is particularly important in two- or three-dimensional situations where the existence of two or more fronts is even more plausible. In such situations the role of transversal dispersivity becomes important too. Spi@ has studied transversal dispersion across a horizontal boundary between two layers of horizontally moving fluids with different densities and with the lighter fluid overlying the heavier one. These studies show that transverse dispersivity also decreases with increase in concentration difference across a front. The non-linear dispersion eqn (16) avoids the difficulties and empiricism associated with estimating an apparent dispersivity. This equation can be used in a threedimensional flow regime and it does not rely on the existence of a single front in the system. Note that the decrease of the width of the dispersion zone (i.e. smaller apparent dispersivity) in high-concentration-gradient experiments does not mean that dispersive mass flux will be smaller compared to situations with a low concentration gradient. On the contrary, due to much higher mass fraction gradients, dispersive mass flux is high and, as a result, the resistance of the system to dispersion will be also high. In other words, the ease with which such a large amount of mass can be dispersed is decreased (lower apparent dispersivity). This explains the need for keeping the second order term in the constitutive eqn (13). It is worth noting that a numerical experiment was also performed in which a high-concentration solution was displaced by a slightly higher concentration solution and the dipsersive mass flux was calculated. It turned out that the dispersive mass flux has the same order of magnitude as that for LCG-experiments. This indicates that the gradient of concentration (or mass fraction), and not the value of concentration, causes deviations from the linear theory. It is known that dispersion is a scale-dependent phenomenon and its magnitude depends on the scale of our models. The non-linear effects described here are expected to be important at the first macroscopic scale of a porous medium; that is, one scale above the pore scale, also called the D,arcy scale. Thus, in the cases where additional mixing as a result of instabilities or
dispersion in porous media
213
heterogeneities occur, the high-concentration-gradient effects may or may not be important depending on the modelling scale. For example, in situations where a lighter fluid underlies a heavier fluid, instabilities may develop and fingering may occur. In those cases, two approaches are possible. One may choose to include the details of fingering and model the system as two dimensional. In this approach, dispersion is modelled at the Darcy scale and the non-linear effects may prove to be important. Alternatively, one may choose to model the system at a higher scale. Thus, one may lump the effects of fingering into an effective dispersion term and model the system as one-dimensional. In such an approach, basically an additional averaging is carried out and the dispersion is modelled at the scale of instabilities. Thus, the effective dispersion will be large and the non-linear effects discussed in this manuscript, will probably be negligible. We are confident that the deviations from Fick’s law observed in our study are not due to the presence of heterogeneities. This can be shown by using the results of Welty and Gelhar22 who have studied the effects of density and viscosity variations on macro-dispersivity in heterogeneous media. They provide an analytical expression for longitudinal macro-dispersivity, as a function of mean displacement distance, mean concentration and concentration gradient, density and viscosity coefficients, mean velocity, gravity, and correlation scale and variance of the log permeability process. They too find that for stable cases, the value of macrodispersivity decreases as the magnitude of concentration gradient increases. The main dispersion mechanism in the system they are considering is due to the spatial variability in flow velocity caused by local heterogeneities in permeability. The system we are considering, however, is relatively homogeneous and the main dispersion mechanism is due to variations in interstitial flow velocity and concentration gradients. Given the very low value of longitudinal dispersivity we have calculated, for the formula derived by Welty and Gelhar” to be applicable, the porous medium in our study should have a correlation length as small as O-1cm and the permeability must vary by a factor of about 6. These values are highly irrelevant to the homogeneous system in our study whose permeability has a very large correlation length and a small variance. Thus, deviations observed in our studies cannot be due to local heterogeneities, but are caused by microscale velocity variations under the influence of density effects. In this study, only the applicability of the Fickian equation has been explored and no comment can be made on the validity of Darcy’s law. For investigating that question, it would have been necessary to have pressure measurements inside the column and compare them with the calculated pressure distribution; these were unfortunately unavailable for HCG experiments because the manometer readings were useless for displacement
214
S. Majid Hassanizadeh, A. L-e&se
experiments, certainly for high-concentration-gradient experiments. The reason is that during displacement experiments, the density of the water column in manometer tubes varies and, thus, the water pressure cannot be linearly related to the manometer level. However, in a study by Schincariol and Schwartz,” it has been observed that a plume of dense salt solution moving horizontally in an ambient fresh water has a velocity much higher than that of a tracer injected with the salt solution. They have also observed brine movement occurring against the prevailing hydraulic gradient. These observations suggest that the pressure gradient may not be the only important term in Darcy’s law, and a term accounting for the large concentration gradients as proposed by Hassanizadeh7 may also be needed. A peculiar result obtained in our experiments has been the dependence of permeability on the salt mass fraction. This can t;e explained by the fact that the resistivity of the porous medium is actually a function of the fluid as well as the matrix properties. Obviously, the resistance of the medium to the flow of brine is larger than the resistance to the flow of fresh water. Traditionally, the resistivity is assumed to be given by p/k with the additional assumption that k is independent of fluid properties. Although viscosity increases with an increase in salt mass fraction, our results show that the increase of resistance to flow is even larger. Thus, the assumption that k is independent of fluid properties is not necessarily valid for a high salt concentration solution. Additional laboratory and field experiments are required to further investigate the issues and questions raised here. In particular, applicability of the non-linear dispersion equation under more general conditions and validity of Darcy’s equation under large concentration gradient conditions have to be studied. The significance of non-linear effects in heterogeneous media and under field conditions needs to be investigated. Additional experiments are also required to determine the dependence of /3 on flow velocity and the possible dependence of permeability on solute concentration.
CONCLUSIONS The Fickian description of dispersion in a porous medium is based on the assumption that the dispersive mass flux is small such that the second order term in the constitutive relation (13) is negligible. Experiments have shown that this assumption is valid only if small concentration gradients are present in the system. For large concentration gradients, large amounts of mass tend to be dispersed and the above-mentioned assumption breaks down. As a result, the breakthrough curves calculated based on the Fickian-type dispersion eqn (14) deviate substantially from the measured curves. It appears that under high-concentration-gradient condi-
tions, the dispersive mass flux will be large so that the non-linear dispersion eqn (17) has to be employed. Both low- and high-concentration-gradient experiments are simulated successfully with this equation. The model parameters aL. and p appear to be independent of the fluid properties. The coefficient p could, however, be a function of the flow velocity.
ACKNOWLEDGMENTS The experiments reported here have been originally designed and set up in the Geotechnique Laboratory of the Delft University of Technology, The Netherlands. The expertise of Ab Mensinga and Joop van Leeuwen of this laboratory has been instrumental in the success of these experiments. Discussions with Jan van Eijkeren of Center for Mathematical Methods of RIVh4 have been very useful in developing the non-linear theory. Wilco de Vries of RIVM has assisted with the calibration of electrodes and the flowmeter. Their contributions are gratefully appreciated. The final version of this paper was written when the senior author was visiting the Department of Civil Engineering and Geological Sciences of the University of Notre Dame. Financial support for this visit was generously made available by RIVM and the University of Notre Dame.
REFERENCES 1. Bear, J., Dynamics of Flui& in Porous Media. Elsevier, New York, USA, 1972. 2. Ben Salah, M. D., Influence des contrastes de viscositk et de densiti sur le dkplacement en milieu poreux de deux fluides miscibles. Revue de l‘lnstitut Frangais du Petrole, (1965) 1237-55. 3. Bues, M. A., Aachib, M. and Zilliox, L., Incidence of heterogeneities on pollutant transport-density and viscosity contrasts of the liquid phase-structure of solid matrix. In Contaminant Transport in Groundwater, 4. H.E. Kobus & W. Kinzelbach, 1990, pp. 251-7. 4. Bues, M. A. & Zilliox, L., D&placement miscible avec contrastes de densitC et de viscositk en milieu poreux. Identification des parametres de d&placementet stabi& en regime de dispersion mkanique. J. Hydrol., 120 (1990) 125-41. 5. Brigham, W. E., Reed, P. W. & Dew, J. N., Experiments on mixing during miscible displacements in porous media. Sot. Petrol. Engng J., March (1961) 1-8. 6. Hassanizadeh, S. M., Derivation of basic equations of mass transport in porous media. Part 1. Macroscopic balance laws. Adv. Water Resour., 9 (1986) 196-206. 7. Hassanizadeh, S. M., Derivation of basic equations of mass transport in porous media. Part 2. Generalized Darcy’s and Fick’s laws. Adv. Water Resour., 9 (1986) 207-22. 8. Hassanizadeh, S. M., Experimental study of coupled flow and mass transport: A model validation exercise. In ModelCare Modelling,
90: Calibration and Reliability
in Groundwater
IAHS Publ. No. 195, 1990, pp. 241-50. 9. Hassanizadeh, S. M. & Gray, W. G., High velocity flow in porous media. Transport in Porous Media, 2 (1987) 52 1-3 1.
High-concentration-gradient
dispersion in porous media
215
10. Hassanizadeh, S. M. 81 Leijnse, A., On the modeling of brine transport in porous media. Water Resow. Res., 24 (1988) 321-30. 11. Hassanizadeh, S. M., Leijnse, A., De Vries, W. J. & Stapper, R. A. M., Experimental Study of Brine Transport in Porous Media. Report 728514005, RIVM, Bilthoven, The Netherlands, 1990. 12. Kempers, L. J. T. M., Effects of Fluid Properties on
17. Scheidegger, A. E., The Physics of Flow through Porous Media. 3rd edn, Univ. of Toronto Press, Toronto, Canada, 1974. 18. Schelkes, K. & Vogel, P., Paleohydrological information as an important tool for groundwater modeling of Gorleben site. In Paleohydrological Methods and Their
Hydrodynamic Dispersion: Comparison of Analytical Models to Numerical Simulations. Proc. 1st SPE/IMA Europ. Conf Mathematics of Oil Recovery. Cambridge,
19. Schincariol, R. A. & Schwartz, F. W., An experimental investigation of variable density flow and mixing in homogeneous and heterogeneous media. Water Resour. Res., 26 (1990) 2317-29. 20. Slobod, R. L. & Howlett, W. E., Effects of gravity segregation in laboratory studies of miscible displacement in vertical unconsolidated porous media. Sot. Petrol. Engng J., March (1964) 1-8. 21. Spitz, K., Dispersion in Porosen Medien: Einth_$ von Inhomogeniteiten und Dichteunterschieden. PhD Thesis, University of Stuttgart, Germany, 1985. 22. Welty, C. & Gelhar, L. W., Stochastic analysis of the effects of fluid density and viscosity variability on macrodispersion in heterogeneous porous media. Water Resour. Res., 27 (1991) 2061-75.
UK, 25-27 July 1989. 13. Kempers, L. J. T. M., Dispersive Mixing in Stable and Unstable Miscible Displacements. PhD Thesis, Delft University of Technology, The Netherlands, 199 1. 14. Krupp, H. K. & Elrick, D. E., Density effects in miscible displacement experiments. Soil Sci., 107 (1969) 372-80. 15. Matheron, G. & De Marsily, G., Is transport in porous media always diffusive? A counter example. Water Resour. Res., 16 (1980) 901-17. 16. Newberg, M. A. & Foh, S. E., Measurement of Longitudinal Dispersion Coefficients of Gas Flowing Through Porous Media. SPE Paper 17731, 1988.
Applications for Radioactive Waste Disposal. Proc. OECD/ NEA Wkshop. OECD, Paris, France, 1992, pp. 237-50.