background image

Bull  Eng  Geol  Env  (2000)  59 : 231

–238  7 Q Springer-Verlag

231

Received:  5  July  1999  7 Accepted:  15  April  2000

M.  Zairi  (Y)  7 M.  J.  Rouis
Ecole  Nationale  d’Ingénieurs  de  Sfax,  BP  W,  3038  Sfax,  Tunisia
e-mail:  Moncef.Zairi6enis.rnu.tn
Tel.:  c216-4-274088
Fax:  c216-4-275595

Numerical  and
experimental  simulation
of  pollutants  migration
in  porous  media

M.  Zairi  7 M.  J.  Rouis

Abstract Quantitative 

prediction 

of 

pollutant

migration  is  now  of  great  interest  in  the  field  of
waste 

management 

and 

storage. 

Numerical

modelling of contaminant transport through porous
media  is  an  efficient  tool  to  make  this  prediction
possible.  Site  selection  and  design  of  protective
liners  for  landfills  and  storage  facilities  are  largely
related  to  the  behaviour  of  contaminants  in  such
structures.  This  behaviour  may  be  well  defined  by
numerical  simulations.  This  paper  proposes  a  two-
dimensional  finite  element  model  for  pollutant
migration.  The  implementation  of  the  advection-
dispersion equation in the numerical model and vali-
dation tests by comparison with analytical and semi-
analytical  solutions  are  presented.  Experimental
laboratory tests on cadmium and chloride migration
through  liner  samples  are  also  discussed  and
compared  with  the  results  from  the  numerical
model.

Résumé La  prévision  quantitative  de  la  migration
des  polluants  présente  un  grand  intérêt  dans  le
domaine de la gestion et du stockage des déchets. La
modélisation  numérique  du  transport  des  polluants
à  travers  un  milieu  poreux  est  un  outil  efficace
rendant possible cette prévision. La sélection de sites
et la conception de barrières de protection pour des
installations  de  stockage  dépendent  fortement  du

comportement  des  polluants  dans  de  telles  struc-
tures. Ce comportement peut être bien défini par des
simulations  numériques.  L’article  présente  un
modèle  2D  en  éléments  finis  pour  la  migration  de
polluants.  L’établissement  des  équations  d’advec-
tion-dispersion  dans  le  modèle  numérique  et  les
tests  de  validation  par  comparaison  avec  des  solu-
tions analytiques et semi-analytiques sont présentés.
Des  expériences  en  laboratoire  sur  la  migration
d’ions  cadmium  et  chlorure  à  travers  des  barrières
de protection reconstituées sont également discutées
et  les  résultats  comparés  à  ceux  de  la  modélisation
numérique.

Keywords Environmental protection 7
Groundwater contamination 7 Laboratory studies 7
Models

Mots  clés Protection de l’environnement 7
Pollution des eaux souterraines 7 Études  en
laboratoire 7 Modèles numériques

Introduction

Numerical  modelling  of  pollutant  migration  in  porous
media has recently received a great deal of attention due to
an  increased  interest  in  the  preservation  of  the  quality  of
the environment and particularly the protection of ground-
water  from  various  pollutants.  Mathematical  modelling  of
contaminant  behaviour  in  soils  is  considered  to  be  a
powerful tool for a wide range of problems concerned with
groundwater quality preservation or rehabilitation.
Enormous  efforts  have  been  made  in  recent  years  to
understand  contaminant  behaviour  in  porous  media  and
to  represent  it  by  mathematical  models.  Ogata  (1970)  was
the  first  to  give  a  non-dimensional  analytical  solution  to

background image

M.  Zairi  M.  J.  Rouis

232

Bull  Eng  Geol  Env  (2000)  59 : 231

–238  7 Q Springer-Verlag

the  mass  transport  equation  in  porous  media.  Others
models  were  then  proposed,  e.g.  those  of  Guerghian  et  al.
(1980)  and  Rowe  and  Booker  (1985).  Although  there  are
now many numerical methods that can be used to solve the
mass transport equation, the finite element model is parti-
cularly appropriate for porous media with complex geom-
etry  and  flows.  It  also  permits  a  good  representation  of
porous media heterogeneity.

Numerical  model

Pollutant  behaviour  in  the  soil  is  subject  to  many  pro-
cesses. The water–soil system is dynamic and many chem-
ical,  biological  and  physical  reactions  occur  which  may
affect  pollutant  behaviour.  Contaminant  transport  can  be
reasonably  approximated  by  one-dimensional  solutions,
and  the  two-dimensional  analysis  will  be  adequate  for  the
vast majority of practical applications. In the case of a clay
layer  acting  as  a  liner  below  a  landfill,  for  example,  the
one-dimensional analysis may be adequate when the width
of the landfill is large compared to the thickness of the clay
liner  and  the  direction  of  contaminant  transport  is
predominantly  vertical.  However,  if  the  dimensions  of  the
landfill  are  comparable  to  the  clay  liner  thickness,  hori-
zontal  and  vertical  transport  may  occur  and  the  two-
dimensional  analysis  is  more  appropriate,  especially  when
concentration at points outside the landfill is of concern.
The  governing  equation  of  two-dimensional  pollutant
migration  in  saturated  porous  media,  also  named  the
advection-dispersion  equation,  is  given  by:

C

t

p

D*

x

f

2

C

x

2

c

D*

y

f

2

C

y

2

P

V*

x

C

x

P

V*

y

C

y

(1)

where  C  is  the  concentration  of  the  considered  contami-
nant  at  a  point  of  coordinates  (x,  y)  at  time  t;  D*pD/R,
V*pV/R;  R  is  the  retardation  factor,  defined  by  Bear
(1972,  in  Auriault  and  Lewandowska  1993)  which  is  given
by:

Rp1crK

d

/n

(1a)

where Vx and Vy are pore water flow velocities in X and Y
directions; Dx and Dy are hydrodynamic dispersion coeffi-
cients  in  X  and  Y  directions;  r  is  soil  bulk  density;  K

is

distribution coefficient; and n is porosity.
In  the  case  of  conservative  species  which  have  no  interac-
tion  with  porous  medium  particles,  the  retardation  factor
is  equal  to  1.  For  species  adsorbed  by  soil  particles,  the
value of this factor depends on the type of element studied
and the soil adsorption capacity expressed by the distribu-
tion or partitioning coefficient. It is appreciated that this is
a  very  simple  description  of  the  interaction  between  the
contaminant  and  the  porous  medium  particles,  but  whilst
it  is  approximate,  it  has  often  been  found  to  be  an
adequate approximation.
K

d

, the distribution coefficient which reflects the degree of

retardation by reversible ion exchange, known as the equi-

librium sorption-desorption model of retardation, requires
particularly  careful  use  in  modelling.  At  high  pH  values,
removal of heavy metals from solutions is generally greater
than the exchange capacity of the soil, due to precipitation
as  metal  carbonate  or  hydroxide.  To  overcome  this,  an
apparent  greater  K

value  is  used  to  include  precipitation

and  fixation  in  modelling  retardation  of  heavy  metals.
In the finite element transformation of Eq. (1), the porous
media (or Eq. 1 domain) is discretized in a number (m) of
elements  (e),  which  constitute  the  finite  element  mesh  or
grid.  The  approximated  solution  for  this  equation  may  be
written  for  each  element  in  the  form:

C

e

(x, y)p

n

A

ip1

N

e

i

C

i

(2)

with

C

e

(x,y): approximate concentration of element (e)

N

i

e

: interpolation functions for each node of element (e)

C

i

: nodal concentration (unknown of the problem)

n:  number  of  nodes  of  the  element  (e).

If C in Eq. (1) is substituted by the expression of Eq. (2), an
approximate solution for Eq. (1)  is  obtained.
The porous medium properties (e.g. hydrodynamic disper-
sion,  flow  velocities,  porosity,  distribution  coefficient)  are
considered  constant  in  each  element  but  may  vary  from
one  element  to  another.  The  second  derivatives  of  the
pollutant-migration-governing  equation  are  generally
undefined for a variety of grid elements. They may be eval-
uated if written in the form of first derivatives, by using an
integration  by  parts  of  the  diffusive  terms  of  the  equation
as given by Green’s theorem for integration in two dimen-
sions (Zienkiewicz 1977).
The  integral  on  the  boundary  of  the  domain  or  contour
integral  resulting  from  the  integration  by  parts  is  intro-
duced in the form of a boundary condition. On the domain
boundary,  with  an  imposed  ingressing  mass  flux,  the
global quantity QC

is substituted for the contour integral,

where Q is the water flux and C

its concentration with the

compound  i.  If  the  previous  developments  are  written  for
each  element  of  the  grid,  the  elementary  equation  may  be
written  in  the  matrix  form  as:

[K

e

]

3

C

1

C

n

4

c

[M

e

]

3

iC

1

it

iC

n

it

4

p

0

(3)

with

[K

e

]p[B]

T

[D][B]c[NN]

T

[V][B]

[M

e

]p[N]

T

[N]

where  elements  of  matrices  [B]  and  [B]

are  composed  of

the nodal values of interpolation function first derivatives.
The  elements  of  matrices  [N],  [N]

and  [NN]

are  the

nodal  values  of  interpolation  functions  and  those  of
matrices  [D]  and  [V]  are  the  elementary  hydrodynamic
dispersion coefficients and flow velocities respectively. The

background image

Pollutant  migration  in  porous  media

Bull  Eng  Geol  Env  (2000)  59 : 231

–238  7 Q Springer-Verlag

233

elementary  matrices  [K

e

]  and  [M

e

]  of  Eq. (3)  are  assem-

bled  for  all  domain  elements  to  give  the  following  global
system,  in  which  m  is  the  number  of  domain  elements.

[K] Cc[M]

iC

it

p

0

(4)

with

[K]p

m

A

ip1

[K

e

]

[M]p

m

A

ip1

[M

e

]

The  global  system  of  Eq. (4)  is  a  first-order  differential
equation  in  time.

Definition  of  boundary  conditions

Two  types  of  boundary  conditions  may  exist  on  the
domain  limits.  The  first  is  a  fixed  concentration  limit
which  is  a  Dirchelet-type  boundary  condition.  The  second
is  a  fixed  concentration  flux  limit  and  a  Neumann-type
boundary  condition.  The  first  condition  is  introduced
when  solving  the  global  equations  system  (Eq. 4)  by  a
restructuring  of  the  global  matrix  [K]  to  eliminate  equa-
tions  corresponding  to  an  imposed  concentration  node
and replacing the corresponding concentration variables in
other  equations  by  the  specified  values.  The  imposed  flux
condition  is  represented  by  the  quantity  QCi  which  is
added to the second member of the elementary equation at
the node i. Each node affected by this condition is consid-
ered as a source term.

Solution  implementation

The  finite  element  spatial  discretization  of  the  governing
equation  of  pollutant  migration  in  saturated  homogenous
porous  media  has  given  a  first-order  differential  equation
in  time:

[KCc[M]

iC

it

p

F

for

t 1 t

0

(5)

C(t

o

)pC

o

In  the  previous  equation,  the  right-hand  side  vector  F
corresponds to the boundary conditions. The resolution of
Eq. (5) is in finding a C(t) function set that satisfies Eq. (5)
at each time t and the initial conditions for tpt

0

. Replacing

C  and  its  time  derivatives  by  their  approximate  finite
differences,  writing  the  boundary  conditions  vector  F  for
times  t  and  tcDt  and  using  the  concentration  increment
D

C  gives:

[M]caD[K] DCpDtF

c

D

t

c

1Pa F

t

P

[KC

t

(6)

with

D

CpC

tcDt

P

C

t

is  the  Euler  constant.  For  a 1 0.5,  Eq. (6)  solution  is

unconditionally  stable.  If  0~a~0.5,  the  solution  is  stable
only  for  a  time  step  Dt~Dt

(critical  time  step)  given  by:

D

t

c

p

2

lP2a l

max

(7)

where  l

max 

is  the  greatest  eigenvalue  of  the  product  of

matrices [M]

–1 

and [K] (Dhatt and Tozot 1981). In order to

solve  Eq. (6),  the  Newton-Raphson  method  (as  described
by  Dhatt  and  Tozot  1981)  may  solve  linear  or  non-linear
non-stationary  problems.  This  method  was  generally  used
to  solve  non-linear  problems.  In  the  constructed  algo-
rithm,  for  each  iteration  Eq. (6)  is  solved  by  the  Gauss
elimination  method.
The  above  developments  were  coded  in  a  FORTRAN
program  named  MIGPO  (Zairi  1996),  using  eight  nodes
quadrilateral  isoparametric  elements  for  grid  generation.
More details about these elements, interpolation functions
and  coordinate  transformations  are  given  in  Zienkiewicz
(1977).

Conditions  for  numerical  solution  stability

The discretization is the preliminary stage in transforming
a  system  of  partial  derivatives  equations  into  an  algebraic
one.  In  this  process  derivative  approximation  generates
truncation errors which cause the numerical dispersion of
the solution. Solution stability studies are generally carried
out  by  a  Fourier  series  development  of  the  amplitude  and
the  phase  of  the  numerical  solution;  they  may  then  be
compared to those of the analytical solution. Guerghian et
al.  (1980)  and  Fletcher  (1988)  have  demonstrated  that  the
stability  of  the  advection-dispersion  equation  by  a
Galerkin  finite  element  method  formulation  is  related  to
two  parameters:  the  Peclet  number  (Pe)  and  the  Courant
number (Cr). These two parameters are expressed as:

PepVDL/D

(8)

CrpVDt/DL

(9)

with

V:  flow  velocity
D

L:  characteristic  length  of  the  element  (DX  or  DY)

D: hydrodynamic dispersion coefficient
D

t:  time  step.

The  solution  stability  studies  relative  to  these  parameters
allow  definition  of  the  domain  of  applicability  of  a  given
numerical  scheme.  These  parameters  also  allow  the  deter-
mination  of  the  optimum  grid  and  time  step  for  a  stable
solution of the problem under consideration. In the litera-
ture,  a  value  of  less  than  or  equal  to  2  is  considered  as  a
limit for the Peclet number; the Courant number generally
ranges  between  0  and  1.6  for  a  Galerkin  finite  element
numerical scheme.
In order to determine the stability domain for the solution
scheme used in the MIGPO model construction relative to
the  Peclet  and  Courant  numbers,  a  40-element,  147-node
grid  is  considered.  Many  simulations  were  executed  for  a
large range of these two factors in order to test the conver-
gence and the stability of the MIGPO solution scheme. The
simulation  tests  have  shown  that  convergent  and  stable
solutions  are  observed  for  a  Peclet  number  ranging  from

background image

M.  Zairi  M.  J.  Rouis

234

Bull  Eng  Geol  Env  (2000)  59 : 231

–238  7 Q Springer-Verlag

Fig. 1

Comparative concentration
breakthrough curves calcu-
lated  from  MIGPO  and  Ogata
solutions for time periods of
a  5 years and b  50 years

0.5 to 2 and a Courant number between 0.25 and 2 and for
a  null  value  of  these  two  numbers.
Non-convergent  solutions  are  defined  as  those  not
reaching  the  convergence  criterion  for  high  iteration
numbers. Numerical oscillation is expressed by negative or
very  high  concentrations  relative  to  initial  concentration
(contamination source). Practically, the stability domain is
reached  by  grid  adjustment,  choosing  an  element  size  (Dx
and  Dy)  that  allows  Peclet  and  Courant  numbers  in  their
optimal ranges. It is also possible to vary the time step and
more  accurate  hydrodynamic  dispersion  coefficients  in
order  to  obtain  a  stable  solution.

Numerical  model  validation  tests

The  efficiency  of  a  mathematical  model  depends  on  its
reliability.  This  is  generally  checked  by  comparisons  with
experimental results or analytical solutions or with numer-
ical  models  whose  efficiency  is  proved.  For  efficiency
testing of the “MIGPO”, three validations were undertaken:
(1)  comparison  with  Ogata’s  analytical  solution  (Ogata
1970);  (2)  comparison  with  Rowe’s  semi-analytical  model
(Rowe  and  Booker  1985);  (3)  experimental  validation  by
comparison with the results of laboratory tests.

Validation  by  comparison  with  Ogata’s  analytical
solution

Ogata (1970) gave a one-dimensional analytical solution of
pollutant migration in porous media. In order to compare
the  MIGPO  results  with  those  obtained  by  this  analytical
method,  a  conservative  solute  with  a  concentration  of
2 mg/l  flowing  through  a  1 m thick  clay  layer  having  a
porosity of 0.43 was considered. The flow velocities were of
0.05 m/year  in  the  flow  direction  and  zero  laterally.  The
hydrodynamic dispersion coefficient was 1.25 and 0.46 m

2

/

year in the parallel and lateral flow directions respectively.
The  grid  used  for  this  simulation  incorporated  four
elements  with  a  time  step  of  5 years.
At  xpo,  a  concentration-type  boundary  condition  of  the
form  C(o,t)pCo  is  assumed  for  both  Ogata  and  the  finite
element  model.  For  the  lower  boundary  of  the  layer,  at
xpL,  the  following  condition  is  applied  to  the  two  solu-
tions:

C

x

L, tp0

(10)

The  concentration  evolution  relative  to  the  initial  concen-
tration  C/C

with  depth  is  presented  in  Fig. 1  for  time

periods of 5 and 50 years. Each of the curves plotted for the
MIGPO corresponds to the concentration evolution at grid
nodes which are lined in the flow direction. Figure 1 shows
that  the  three  profiles  obtained  from  the  MIGPO  results
have  the  same  trend  and  evolution  in  time  as  that  of  the
analytical  solution.  Concentration  values  obtained  by  the
two  methods  are  of  the  same  magnitude,  the  greatest
difference  between  the  two  results  being  some  0.01  (for  a
time  of  5 years).

Validation  by  comparison  with  the  Rowe  semi-analytical
model

The  model  POLUTE,  established  by  Rowe  and  Booker
(1985),  was  also  used  to  verify  the  validity  of  the  MIGPO
results.  This  method  involves  a  transformation  of  all  the
governing  equations  by  introducing  a  Laplace  transform;
the  transformed  equations  are  solved  analytically  and  the
solutions  are  numerically  inverted.  The  flow  of  a  2 mg/l
cadmium solute through a 1 m thick saturated layer having
a  permeability  of  10

–10

m/s  and  a  porosity  of  0.4  was

considered  to  compare  the  results  of  the  two  models.  In
this  theoretical  example,  only  molecular  diffusion  was
considered as a transport mechanism. The dispersion coef-
ficient  used  was  cadmium  in  an  aqueous  solution
(7.17!10

–8

m

2

/s;  in  Quigley  et  al.  1987)  multiplied  by  a

tortuosity  coefficient  of  0.4.  The  cadmium  adsorption  by
soil  particles  was  also  modelled  considering  a  retardation
factor of 3.
The  evolution  with  depth  of  the  ratio  of  the  initial
cadmium  concentration  and  that  after  5  and  50 years  is
presented in Fig. 2. The results from the two models differ
only  slightly  for  the  time  period  of  5 years  (B0.02 mg/l)
with  no  difference  discernible  for  the  50-year  period.
Concentration profiles through the layer appear similar for
both  models  and  show  an  increasing  concentration  with
time.

background image

Pollutant  migration  in  porous  media

Bull  Eng  Geol  Env  (2000)  59 : 231

–238  7 Q Springer-Verlag

235

Fig. 2a,b

Comparative concentration
breakthrough curves calcu-
lated from MIGPO and Rowe
solutions for time periods of
a  5 years and b  50 years

Table 1

Chemical composition of kiln dust and fly ash

Oxide  (%)

CaO

MgO

Na

2

O

K

2

O

Fe

2

O

3

MnO

SiO

2

Al

2

O

3

SO

4

CO

3

Cl

Kiln  dust

42.06

1.56

1.45

2.7

1.83

0.1

13.67

4.3

6.8

24.8

0.73

Fly  ash

19.74

3.66

7.71

1.0

6.06

0.13

31.12

18.2

7.76

4.3

0.3

Table 3

Composition of liners 1 and 2

Liner

Kiln  dust
(%)

Fly  ash
(%)

Silica  steam
(%)

Liner  1

95.24

0

4.76

Liner  2

76.19

19.05

4.76

Table 2

Mineral composition of kiln dust and fly ash

Mineral  phase

Main  phase

Secondary  phase

Kiln  dust

Calcite,  lime,  quartz,
anhydrite,  arcanite,  alite

Magnesite

Fly  ash

Amorphous

Quartz,  magnesite

Many  simulation  tests  have  demonstrated  a  good  agree-
ment between MIGPO and POLUTE results both quantita-
tively and qualitatively. This is more evident for long time
periods and narrow grids.

Experimental  study  of  pollutant  migration  through
liners

The  experimental  study  of  pollutant  migration  through
liner samples is part of a research programme on the use of
fly  ash  and  stabilised  kiln  dust  as  major  components  for
the construction of liners. It is hoped that this may provide
a potential solution to the environmental problems caused
by  this  industrial  refuse  storage,  particularly  in  five  Tuni-
sian cities having cement works and a projected landfill in
their  vicinity.
The  aim  of  the  experimental  work  was  to  outline  the
mechanical  and  hydrogeological  characteristics  of  the
liners being developed. Numerous experimental tests were
undertaken  on  a  variety  of  liners,  e.g.  permeability,
leaching  and  longevity  tests.  In  this  paper  the  pollutant
migration tests for two liners are presented.
The  mineralogical  and  chemical  compositions  of  the  kiln
dust  and  fly  ash  used  are  presented  in  Tables 1  and  2.  As
can  be  seen,  both  the  alkali  (Na,  K)  and  sulphate  concen-
trations  were  high.  The  kiln  dust  is  characterised  by  high
calcite and lime concentrations, while the fly ash was a C-
ASTM type with a moderate lime concentration compared
to that of kiln dust and particularly high concentrations of
alkalis, iron and alumina.
A number of liners were constructed using a variety of kiln
dust and fly ash mixtures. The compositions of liners 1 and
2  discussed  in  this  work  are  presented  in  Table 3.  Silica
steam was added to liners 1 and 2 to ensure a silica source

and  to  stabilise  the  mixtures  (Rouis  1991).  The  prepared
mixtures were put into a PVC cylinder, 70 mm in diameter
and 20 mm in height. The mixture was slightly compacted
in three layers and the liners thus formed were matured for
1 week  in  a  wet  room  prior  to  being  left  saturated  in
distilled water for 3 weeks. Four weeks after preparation of
the  mixture, the distilled water  was  replaced by a solution
with  a  known  concentration  of  cadmium  chloride  (CdCl

2

,

2H

2

O). Samples of the floating solution were taken weekly

over the 67 days of the test. Each time a sample was taken,
it was replaced by an equivalent volume of water so that a
150-mm head was maintained above the liner. A sufficient
quantity  of  cadmium  chloride  was  dissolved  in  this  ad-
ditional  water  in  order  to  maintain  the  chloride  (Cl

)

background image

M.  Zairi  M.  J.  Rouis

236

Bull  Eng  Geol  Env  (2000)  59 : 231

–238  7 Q Springer-Verlag

Fig. 3

Chloride  concentration  in  liner  interstitial  water  at  end  of  diffu-
sion  test

Fig. 4

Variation in cadmium concentration in solute over liner 1

Fig. 5

Cadmium  concentration  in  liners  1  and  2  at  end  of  diffusion
test

concentration  at  1000 mg/l  and  the  cadmium  (Cd

2c

)

concentration at 1500 mg/l. At the end of the test, the 100-
mm  liner  sample  was  cut  into  five  20-mm-thick  slices.
Each  slice  was  dried  for  24 h  at  105 7C,  then  crushed  and
the  chloride  concentrations  determined  on  a  solid  sample
dissolved by acid.
During  the  67 days  of  testing,  no  water  flow  took  place  in
the  apparatus,  which  had  a  leachate  collecting  device.  The
pH of the liner varied from 7 on the first day to 11.5 on the
seventh  day  and  then  stabilised  at  between  11.5  and  12.

Chloride  migration

The  chloride  concentration  of  the  interstitial  water  in  the
liners was obtained by reducing the quantities given by the
chemical  analyses  to  those  of  the  water  in  the  initial
sample prior to drying. This water content is considered to
equal  effective  porosity  of  the  liner,  measured  by  the
mercury  porosimeter.  In  this  way,  chloride  concentration
variation  in  the  solute  over  the  liner  can  be  determined
together  with  the  ion  concentration  profile  in  the  liner.
The  chloride  concentration  of  the  source  solution
decreased  with  this  ion  migration  into  the  liner,  from
1000 mg/l  at  the  beginning  of  each  week  (after  the  adjust-
ment  to  the  initial  concentration)  to  750 mg/l  at  the  end,
just  before  water  sampling.  A  chloride  concentration  of
450 mg/l  in  the  interstitial  water  at  the  liner  base  was
reached after 67 days of flow (Fig. 3). This is related to the
absence  of  any  interaction  between  the  chloride  ion  and
the liner constituents.

Cadmium  migration

The  variation  in  cadmium  concentration  of  the  source
solution with time for liner 1 is presented in Fig. 4. A sharp
decrease in cadmium concentration was recorded from the
first day of the experiment, from 1500 to 11.5 mg/l. It then
stabilised at about 0.05 mg/l, even though a cadmium-rich

solute  was  added  after  each  weekly  sampling.  A  similar
behaviour was recorded for all compositions of liner.
This  decrease  is  related  to  the  precipitation  of  cadmium
carbonate as a result of the increase in pH due to portlan-
dite  dissolution  from  the  liner  into  the  solute.  This  is
confirmed by the cadmium concentration profiles through
liners  1  and  2  (Fig. 5)  which  show  a  sharp  decrease  in  the
first  200 mm  of  sample.  In  their  study  of  the  influence  of
pH  on  selectivity  and  retention  of  heavy  metals  in  clay
soils, Yong and Phadungchewit (1993) found that the pres-
ence  of  carbonates  in  the  soil  contributes  significantly  to
their retention capability. At high pH values, precipitation
mechanisms, (for example as hydroxides and/or as carbon-
ates)  dominate  the  cation  exchange  capacity.
It  remains  difficult  to  be  confident  of  the  real  cadmium
concentration  in  the  interstitial  water  in  the  liner.  It  is
likely  that  most  of  the  analysed  cadmium  in  the  liners  is
from the precipitates which cover their pore surfaces. This

background image

Pollutant  migration  in  porous  media

Bull  Eng  Geol  Env  (2000)  59 : 231

–238  7 Q Springer-Verlag

237

Fig. 6a,b

Comparison of liner chloride
concentration profiles from
MIGPO and experimental
tests  for  a  10 elements grid
and b  20 elements grid

leads  to  difficulties  in  calculating  distribution  and  hydro-
dynamic dispersion coefficients.

Numerical  simulation  of  experimental  tests

The  results  of  the  pollutant  migration  tests  on  liners  have
been  used  for  the  experimental  validation  of  the  MIGPO
model.  For  the  numerical  simulation  of  these  tests  the
following procedure was adopted:
1. The  domain  of  the  liner  sample  was  divided  into  ten

elements.

2. Throughout  the  test,  chloride  concentration  was  main-

tained  at  1000 mg/l  at  nodes  on  the  limit  between  the
cadmium chloride solute and the liner.

3. Chloride 

diffusion 

coefficients 

were 

taken 

as

DxpDyp0.02 m

2

/year.

4. The  flow  velocities  Vx  and  Vy  were  taken  as  zero  by

assuming  that  the  transport  mechanisms  are  limited  to
molecular diffusion only.

5. Conservative  values  of  chloride  ions  were  considered  –

the retardation factor equalling 1.

6. Simulations  were  carried  out  for  ten  time  steps  of

6.7 days each.

Simulation  results

The  simulated  and  measured  interstitial  water  chloride
concentrations through the liner sample over a time period
of  67 days  are  shown  in  Fig. 6a.  The  curves  obtained  with
the MIGPO model represent the average concentration for
each  20-mm-thick  slice  of  liner  sample.  As  can  be  seen
from  the  figure,  the  simulated  concentrations  are  close  to
those obtained experimentally, with both profiles having a
similar appearance.
Experiments were also undertaken using the same simula-
tion conditions but with a finer grid of 20 elements and 79
nodes  (Dxp1.75 cm  and  Dyp2 cm).  The  time  interval  in
this case was 2.68 days. Chloride concentration profiles for
the  simulated  liners  are  shown  in  Fig. 6b.  As  can  be  seen
from  the  figure,  there  is  a  closer  correlation  with  the
experimental results than in the first tests, suggesting that
the  grid  refinement  reduces  the  difference  between  the
MIGPO and experimental concentration profiles.

The  inability  of  the  model  to  reproduce  the  Cd  data  is
related  to  the  fact  that  speciation  of  Cd  with  chloride  is
important  in  evaluating  its  mobility.  Results  from  theore-
tical and experimental studies (Yong and Sheremata 1991)
indicated  that  a  better  insight  into  the  characteristics  of
heavy  metal  adsorption  modelling  and  prediction  of
contaminant transport in soils can be obtained if the effect
of  the  other  constituents  in  the  solutions  is  taken  into
account.  Thermodynamic  calculations  may  be  used  in
predicting the effect of the speciation and complexation by
inorganic  and  organic  ligands.  The  model  presented  here
should  be  completed  with  a  thermodynamic  calculation
sequence  to  allow  the  prediction  of  such  contaminant
transport.

Conclusions

This  paper  has  outlined  a  two-dimensional  finite  element
solution  for  the  advection-dispersion  equation.  The  tech-
niques  discussed  have  been  implemented  in  the  MIGPO
program.  In  order  to  optimise  the  use  of  this  model,  a
limited  study  of  the  numerical  dispersion  and  oscillation
criteria was undertaken by an examination of a large range
of Peclet and Courant numbers.
Validation  of  the  MIGPO  results  was  carried  out  by
comparing the results with those obtained using two other
techniques:  the  Ogata  analytical  solution  and  Rowe  semi-
analytical solution, together with experimental results. The
numerical simulation results were in good agreement with
those of analytical and semi-analytical solutions, while the
concentrations  determined  were  of  the  same  order  of
magnitude  and  with  a  similar  profile.
The  experimental  study  of  pollutant  migration  through
liners indicates that the behaviour of contaminants in such
structures  depends  on  the  type  of  species  considered,  the
composition  of  the  liner  and  the  physical  and  chemical
conditions  and  properties  of  the  medium.  The  numerical
simulation test for chloride migration gave results in good
agreement  with  those  measured  during  the  laboratory

background image

M.  Zairi  M.  J.  Rouis

238

Bull  Eng  Geol  Env  (2000)  59 : 231

–238  7 Q Springer-Verlag

tests,  suggesting  that  this  is  an  appropriate,  conservative
model of contaminant migration. Laboratory measurement
of  the  distribution  of  cadmium  in  the  liner  material  indi-
cated that more complex models are required to represent
the  chemical  and  physical  processes  controlling  its  migra-
tion.

References

Auriault  JL,  Lewandowska  J 

(1993) Homogenisation analysis

of  diffusion  and  adsorption  macro  transport  in  porous  media:
macro  transport  in  the  absence  of  advection.  Géotechnique
43(3) : 457–469

Bear  J 

(1972) Dynamics of fluids in porous media. Elsevier, New

York, 764 pp

Dhatt  G,  Tozot  G 

(1981)  Une  présentation  de  la  méthode  des

éléments finis. Maloine, Paris

Fletcher 

CAJ 

(1988)  Computational  techniques  for  fluid

dynamics  –  fundamental  and  general  techniques.  Springer,
Berlin Heidelberg New York

Guerghian  AB,  Ward  DS,  Cleary  RW 

(1980) A finite element

model  for  the  migration  of  leachate  from  a  sanitary  landfill  in
Long Island, New York. Water Resour Bull 16(5) : 900–906

Ogata  A 

(1970)  Theory  of  dispersion  in  a  granular  medium.

Fluid  movement  in  earth  materials.  US  Geol  Surv  Prof  Pap
411-I

Quigley  RM,  Yanful  EK,  Fernandez  F 

(1987)  Ion  transfer  by

diffusion  trough  clayey  barrieres.  Proc  Conf  on  Geotechnical
Practice  for  Waste  Disposal,  ASCE,  Ann  Arbor,  15–17  June,
pp 137–158

Rouis  MJ 

(1991)  Conception  et  caractérisation  de  barrières

hydrogéologiques à partir de poussières de four de cimenterie.
Thèse  de  doctorat  et  sciences  appliquées,  University  of  Sher-
brooke, Canada

Rowe  RK,  Booker  JR 

(1985) 1-D pollutant migration in soils of

finite depth. ASCE J Geotech Eng 111(4) : 137–158

Yong  RN,  Phadungchewit  Y 

(1993) pH influence on selectivity

and  retention  of  heavy  metals  in  some  clay  soils.  Rev  Can
Géotech 30 : 821–833

Yong  RN,  Sheremata  TW 

(1991)  Effect  of  chloride  ion  on

adsorption of cadmium from a landfill leachate Can Geotech J
28 : 378–387

Zairi  M 

(1996)  Modélisation  de  la  migration  des  polluants  en

milieux poreux. Cas des sites de stockage de phosphogypse de
Sfax.  Thèse  de  doctorat  de  spécialité,  Ecole  Nationale  d’Ingé-
nieurs de Sfax, Sfax

Zienkiewicz  OC 

(1977)  The  finite  element  method  in  engi-

neering science. McGraw-Hill, New York, 787 pp