Fowler A. Mathematical Geoscience

Подождите немного. Документ загружается.

450 7 Groundwater Flow

Homogenisation The technique of homogenisation is no more than the tech-

nique of averaging in the spatial domain, most often formulated as a multiple scales

method. Whole books have been written about it, for example those by Bensoussan

et al. (1978) and Sanchez-Palencia (1983). For application to porous media, see, for

example, Ene’s article in the book edited by Cushman (1990).

Piping Many dams are built of concrete, and in this case the problems associated

with seepage do not arise, owing to the virtual impermeability of concrete. Earth

and rockﬁll dams do exist, however, and are liable to failure by a mechanism called

piping. The Darcy ﬂow through the porous dam causes channels to form by eroding

away ﬁne particles. The resultant channelisation concentrates the ﬂow, increasing

the force exerted by the ﬂow on the medium and thus increasing the erosion/collapse

rate of the channel wall. We can write Darcy’s law as a force balance on the liquid

phase,

0 =−φ∇p −

φμ

−φρ

gk (7.275)

(k being vertically upwards) and φμv

/k is an interactive drag term; then the cor-

responding force balance for the solid phase is

0 =−(1 −φ)∇p

φμ

−(1 −φ)ρ

gk, (7.276)

where p

is the pressure in the solid. For a granular solid, we can expect grain

motion to occur if the interactive force is large enough to overcome friction and

cohesion; the typical kind of criterion is that the shear stress τ satisﬁes

τ ≥c +p

tanφ, (7.277)

but in view of the large conﬁning pressure and the necessity of dilatancy for soil

deformation, the piping criterion will in practice be satisﬁed at the toe of the dam

(i.e. the front), and piping channels will eat their way back into the dam, in much

the same way that river drainage channels eat their way into a hillslope. A simpler

criterion at the toe then follows from the necessity that the effective pressure on the

grains be positive. A lucid discussion by Bear and Bachmat (1990, p. 153) indicates

that the solid pressure is related to the effective pressure p

which controls grain

deformation by

=(1 −φ)(p

−p), (7.278)

and in this case the piping criterion at the toe is that p

< 0 in the soil there, or

∂p

/∂z > 0. From (7.275), (7.276) and (7.278), this implies piping if

μv

>(ρ

−ρ

)(1 −φ)g, (7.279)

where v is the vertical component of v

. This criterion is given by Bear (1972).

More generally, piping can be expected to occur if p

reaches 0 in the soil interior

(ignoring cohesion). Sellmeijer and Koenders (1991) develop a model for piping.

7.12 Notes and References 451

Taylor Dispersion Taylor dispersion is named after its investigation by Taylor

(1953), who carried out experiments on the dispersal of solute in ﬂow down a tube.

The dispersion is enabled by the combination of differential axial advection by the

down tube velocity, typically a Poiseuille ﬂow, and the rapid cross stream diffusion

which renders the cross-sectional concentration proﬁle radially uniform. The theory

of Taylor is somewhat heuristic; it was later elaborated by Aris (1956). For a formal

derivation using asymptotic methods, see Fowler (1997, p. 222, Exercise 2).

Its application to porous media stems from the conceptual idea that the pore space

consists of a network of narrow tubules connected at pore junctions. If the tubes are

of radius a and length d

, the latter corresponding to grain size, then the Darcy

ﬂux |u|∼πφU, while the pore radius a ∼ d

√

φ. This would suggest a Taylor

dispersion coefﬁcient of

≈

48D

∼

|u|

48π

Dφ

, (7.280)

as opposed to the measured values which more nearly have D

∼|u|. Taylor dis-

persion in porous media has been studied by Saffman (1959), Brenner (1980) and

Rubinstein and Mauri (1986), the latter using the method of homogenisation.

Bioﬁlm Growth Monod kinetics was described by Monod (1949), by way of

analogy with enzyme kinetics, where one considers the uptake of nutrients as oc-

curring through a series of fast intermediary reactions; when two nutrients control

growth, as in respiration, it is usual to take the growth rate as proportional to the

product of two Monod factors (Bader 1978). A variety of enhancements to this

simple model have also been proposed to account for nutrient consumption due to

maintenance, inactivation of cells in adverse conditions, and other observed effects

(Beeftink et al. 1990; Wanner et al. 2006).

Bacteria in soils commonly grow as attached bioﬁlms on soil grains, with a thick-

ness of the order of 100 µ. A variety of models to describe bioﬁlm growth have been

presented, with an ultimate view of being able to parameterise the uptake rate of

contaminant species in soils and other environments (Rittmann and McCarty 1980;

Picioreanu et al. 1998;Eberletal.2001; Dockery and Klapper 2001; Cogan and

Keener 2004).

Remediation Sites The three sites described in Sect. 7.8 are under study by the

Groundwater Restoration and Protection Group at the University of Shefﬁeld, led by

Professor David Lerner. The site at Four Ashes is described by Mayer et al. (2001),

that at Rexco by Hüttmann et al. (2003), and that at St. Alban’s by Wealthall et al.

(2001).

The description of the two species reaction front given by (7.192) is similar to

that for a diffusion ﬂame (Buckmaster and Ludford 1982) in combustion, and also

corrosion in alloys (Hagan et al. 1986). It is not conceptually difﬁcult to extend

this approach to an arbitrary number of reactions, although it may become awkward

when multiple reaction fronts are present (see, for example, Dewynne et al. 1993).

452 7 Groundwater Flow

Diagenesis The ﬁrst order reaction kinetics (7.225) for the smectite–illite transi-

tion was proposed by Eberl and Hower (1976). Information on solubility limits is

given by Aagaard and Helgeson (1983) and Sass et al. (1987). The asymptotic ap-

proximation called here the weak solubility limit is called solid density asymptotics

by Ortoleva (1994). Details of the use of the weak solubility approximation can be

found in Fowler and Yang (2003).

Compaction Interest in compaction is motivated by its occurrence in sedimentary

basins, and also by issues of subsidence due to groundwater or natural gas extraction

(see, for example, Baú et al. 2000). The constitutive law used here for effective

pressure is that of Smith (1971); it mimics the normal consolidation behaviour of

compacting sediments (such as soils), and is further discussed by Audet and Fowler

(1992) and Jones (1994).

Athy’s law comes from the paper by Athy (1930). Smith (1971) advocates the use

of the high exponent m =8in(7.265). Further details of the asymptotic solution of

the compaction proﬁles are given by Fowler and Yang (1998). Freed and Peacor

(1989) show examples of the ﬂattened porosity proﬁles at depth.

Early work on pressure solution in sedimentary basins was by Angevine and Tur-

cotte (1983) and Birchwood and Turcotte (1994). More recently, Fowler and Yang

(1999) derived the viscous compaction law. An extension to viscoelastic compaction

has been studied by Yang (2000).

Seals One process which we have not described is the formation of high pressure

seals. In certain circumstances, pore pressures undergo fairly rapid jumps across

a ‘seal’, typically at depths of 3000 m. Such jumps cannot be predicted within the

conﬁnes of a simple compaction theory, and require a mechanism for pore-blocking.

Mineralisation is one such mechanism, as some seals are found to be mineralised

with calcite and silica (Hunt 1990). In fact, a generalisation of the clay diagenesis

model to allow for calcite precipitation could be used for this purpose. As it stands,

(7.232) predicts a source for φ, but mineralisation would cause a corresponding sink

term. Reduction of φ leads to reduction of diffusive transport, and the feedback is

self-promoting. Problems of this type have been studied by Ortoleva (1994), for

example.

7.13 Exercises

7.1 Show that for a porous medium idealised as a cubical network of tubes, the per-

meability is given (approximately) by k = d

/72π, where d

is the grain

size. How is the result modiﬁed if the pore space is taken to consist of pla-

nar sheets between identical cubical blocks? (The volume ﬂux per unit width

between two parallel plates a distance h apart is −h



/12μ, where p



is the

pressure gradient.)

7.13 Exercises 453

7.2 A sedimentary rock sequence consists of two types of rock with permeabilities

and k

. Show that in a unit with two horizontal layers of thickness d

and

, the effective horizontal permeability (parallel to the bedding plane) is



where f

/(d

), whereas the effective vertical permeability is given

−1

⊥

−1

Show how to generalise this result to a sequence of n layers of thickness

,...,d

Hence show that the effective permeabilities of a thick stratigraphic se-

quence containing a distribution of (thin) layers, with the proportion of layers

having permeabilities in (k, k +dk) being f(k)dk,aregivenby





∞

kf (k ) d k , k

−1

⊥



∞

f(k)dk

7.3 Groundwater ﬂows between an impermeable basement at z =h

(x,y,t) and

a phreatic surface at z = z

(x,y,t). Write down the equations governing the

ﬂow, and by using the Dupuit approximation, show that the saturated depth h

satisﬁes

φh

kρg

∇.[h∇z

where ∇ = (∂/∂x, ∂/∂y). Deduce that a suitable time scale for ﬂows in an

aquifer of typical depth h

and extent l is t

=φμl

/kρgh

I live a kilometer from the river, on top of a layer of sediments 100 m thick

(below which is impermeable basement). What sort of sediments would those

need to be if the river responds to rainfall at my house within a day; within a

year?

7.4 A two-dimensional earth dam with vertical sides at x = 0 and x = l has a

reservoir on one side (x<0) where the water depth is h

, and horizontal dry

land on the other side, in x>l. The dam is underlain by an impermeable

basement at z =0.

Write down the equations describing the saturated groundwater ﬂow, and

show that they can be written in the dimensionless form

u =−p

,ε

w =−(p

+1),

+ε

=0,

and deﬁne the parameter ε. Write down suitable boundary conditions on the

impermeable basement, and on the phreatic surface z =h(x, t).

Assuming ε 1, derive the Dupuit–Forchheimer approximation for h,

=(hh

)

in 0 <x<1.

454 7 Groundwater Flow

Show that a suitable boundary condition for h at x =0 (the dam end) is

h =1atx =0.

Now deﬁne the quantity

U =



pdz,

and show that the horizontal ﬂux

q =



udz=−

∂U

∂x

Hence show that the conditions of hydrostatic pressure at x = 0 and con-

stant (atmospheric) pressure at x =1 (the seepage face) imply that



qdx=

Deduce that, if the Dupuit approximation for the ﬂux is valid all the way to

the toe of the dam at x =1, then h =0atx =1, and show that in the steady

state, the (dimensional) discharge at the seepage face is

kρgh

2μl

Supposing the above description of the solution away from the toe to be

valid, show that a possible boundary layer structure near x = 1 can be de-

scribed by writing

x =1 −ε

X, h =εH, z =εZ, p =εP,

and write down the resulting leading order boundary value problem for P .

7.5 I get my water supply from a well in my garden. The well is of depth h

(relative to the height of the water table a large distance away) and radius r

Show that the Dupuit approximation for the water table height h is

∂h

∂t

kρg

∂

∂r



∂h

∂r



If my well is supplied from a reservoir at r =l, where h =h

, and I withdraw

a constant water ﬂux q

, ﬁnd a steady solution for h, and deduce that my well

will run dry if

πkρgh

μ ln[l/r

]

Use plausible values to estimate the maximum yield (litres per day) I can use

if my well is drilled through sand, silt or clay, respectively.

7.13 Exercises 455

7.6 A volume V of efﬂuent is released into the ground at a point (r =0) at time t.

Use the Dupuit approximation to motivate the model

∂h

∂t

kρg

∂

∂r



∂h

∂r



h =h

at t =0,r >0,



∞

r(h −h

)dr =V/2π, t > 0,

where h

is the initial height of the water table above an impermeable base-

ment. Find suitable similarity solutions in the two cases (i) h

=0, (ii) h

> 0,

h −h

h

, and comment on the differences you ﬁnd.

7.7 Fluid ﬂows through a porous medium in the x direction at a linear velocity U.

At t =0, a contaminant of concentration c

is introduced at x =0. If the longi-

tudinal dispersivity of the medium is D, write down the equation which deter-

mines the concentration c in x>0, together with suitable initial and boundary

conditions. Hence show that c is given by



erfc



x −Ut

√



+exp





erfc



x +Ut

√



where

erfc ξ =

√



∞

−s

ds.

[Hint: you might try Laplace transforms, or else simply verify the result.]

Show that for large ξ, erfcξ = e

−ξ

[

√

πξ

+···], and deduce that if x =

Ut +2

√

Dt η, with η =O(1), then

≈

erfc η +O



√



Hence show that at a ﬁxed station x =X far downstream, the measured proﬁle

is approximately given by

c ≈c



1 −

erfc







1/2



t −



This is called the breakthrough curve, and indicates that dispersion causes

breakthrough to occur over a time interval (at large distance) of order t

(DX/U

)

1/2

.IfD ≈aU, show that the ratio of t

to t

=X/U is t

∼

(a/X)

1/2

7.8 Rain falls steadily at a rate q (volume per unit area per unit time) on a soil

of saturated hydraulic conductivity K

(=k

g/μ, where k

is the saturated

456 7 Groundwater Flow

permeability). By plotting the relative permeability k

and suction character-

istic σψ/d as functions of S (assuming a residual liquid saturation S

), show

that a reasonable form to choose for k

(ψ) is k

−cψ

. If the water table

is at depth h, show that, in a steady state, ψ is given as a function of the di-

mensionless depth z

∗

=z/z

, where z

=σ/ρ

gd (σ is the surface tension, d

the grain size), by

∗

−z

∗

ψ −



sinh



(ln

∗

−cψ)



sinh



∗





where h

∗

=h/z

, providing q

∗

=q/K

< 1. Deduce that if h z

, then ψ ≈

∗

near the surface. What happens if q>K

7.9 Derive the Richards equation

∂S

∂t

=−

∂

∂z



(S)



∂p

∂z

+ρ



for one-dimensional inﬁltration of water into a dry soil, explaining the mean-

ing of the terms, and giving suitable boundary conditions when the surface ﬂux

q is prescribed. Show that if the surface ﬂux is large compared with k

g/μ,

where k

is the saturated permeability, then the Richards equation can be ap-

proximated, in suitable non-dimensional form, by a nonlinear diffusion equa-

tion of the form

∂S

∂t

∂

∂z



∂S

∂z



Show that, if D =S

, a similarity solution exists in the form

S =t

F(η), η=z/t

where α =

m+2

, β =

m+1

m+2

, and F satisﬁes



)



=αF −βηF



=−1atη =0,F→0asη →∞.

Deduce that



=−(α +β)



Fdη−βηF,

where η

(which may be ∞)iswhereF ﬁrst reaches zero. Deduce that F



< 0,

and hence that η

must be ﬁnite, and is determined by



Fdη=

α +β

What happens for t>F(0)

−1/α

7.13 Exercises 457

7.10 Write down the equations describing one-dimensional consolidation of wet

sediments in terms of the variables φ,v

,p,p

, these being the porosity,

solid and liquid (linear) velocities, and the pore and effective pressures. Ne-

glect the effect of gravity.

Saturated sediments of depth h lie on a rigid but permeable (to water) base-

ment, through which a water ﬂux W is removed. Show that

∂p

∂z

−W,

and deduce that φ satisﬁes the equation

∂φ

∂t

∂

∂z



(1 −φ)



∂p

∂z

−W



If the sediments are overlain by water, so that p = constant (take p = 0) at

z =h, and if φ =φ

+p/K, where the compressibility K is large (so φ ≈φ

show that a suitable reduction of the model is

∂p

∂t

−W

∂p

∂z

∂

∂z

where c =K(1−φ

)k/μ, and p =0onz =h, p

=μW/k. Non-dimensionalise

the model using the length scale h, time scale h

/c, and pressure scale

μW h/k. Hence describe the solution if the parameter ε = μWh/k is small,

and ﬁnd the rate of surface subsidence. What has this to do with Venice?

7.11 Write down a model for vertical ﬂow of two immiscible ﬂuids in a porous

medium. Deduce that the saturation S of the wetting phase satisﬁes the equa-

tion

∂S

∂t

∂

∂z



eff



+gρ



=−

∂

∂z



eff

∂p

∂z



where z is a coordinate pointing downwards,

−p

,ρ=ρ

−ρ

−1

eff



−1



q is the total downward ﬂux, and the sufﬁxes w and nw refer to the wetting

and non-wetting ﬂuid, respectively. Deﬁne the phase mobilities M

.Givea

criterion on the capillary suction p

which allows the Buckley–Leverett ap-

proximation to be made, and show that for q =0 and μ

μ

, waves typ-

ically propagate downwards and form shocks. What happens if q =0? Is the

Buckley–Leverett approximation realistic—e.g. for air and water in soil? (As-

sume p

∼2γ/r

, where γ =70 mN m

−1

, and r

is the pore radius: for clay,

silt and sand, take r

=1 µ, 10 µ, 100 µ, respectively.)

7.12 A model for snow-melt run-off is given by the following equations:

u =



∂p

∂z

+ρ



458 7 Groundwater Flow

k =k

∂S

∂t

∂u

∂z

=0,



−S



Explain the meaning of the terms in these equations, and describe the assump-

tions of the model.

The intrinsic permeability k

is given by

=0.077 d

exp[−7.8ρ

/ρ

where ρ

and ρ

are snow and water densities, and d is grain size. Take d =

1 mm, ρ

= 300 kg m

−3

, ρ

= 10

kg m

−3

, p

= 1kPa,φ = 0.4, μ = 1.8 ×

−3

Pa s, g = 10 m s

−2

, and derive a non-dimensional model for melting

of a one metre thick snow pack at a rate (i.e. u at the top surface z = 0) of

−6

−1

. Determine whether capillary effects are small; describe the nature

of the model equation, and ﬁnd an approximate solution for the melting of an

initially dry snowpack. What is the (meltwater ﬂux) run-off curve?

7.13 Consider the following model, which represents the release of a unit quan-

tity of groundwater at t = 0 in an aquifer −∞ <x<∞, when the Dupuit

approximation is used:

=(hh

)

h =0att =0,x=0,



∞

−∞

hdx =1

(i.e., h =δ(x) at t =0). Show that a similarity solution to this problem exists

in the form

h =t

−1/3

g(ξ), ξ =x/t

1/3

and ﬁnd the equation and boundary conditions satisﬁed by g. Show that the

water body spreads at a ﬁnite rate, and calculate what this is.

Formulate the equivalent problem in three dimensions, and write down the

equation satisﬁed by the similarity form of the solution, assuming cylindrical

symmetry. Does this solution have the same properties as the one-dimensional

solution?

7.14 The tensor D

(i, j =1, 2, 3) has three invariants

III

(Summation over repeated indices is implied.) Show that the invariants of the

tensor

=α

⊥

uδ

+(α



−α

⊥

)

7.13 Exercises 459

where u =|u| and δ

is the Kronecker delta (= 1ifi = j , = 0ifi = j ), are

the same as those of the tensor

D =

⎛

⎝



u 00

0 α

⊥

u 0

00α

⊥

⎞

⎠

7.15 Suppose that a doubly porous medium consists of a periodic sequence of

blocks M with boundaries (fractures) ∂M. The concentration of a chemical

reactant c is taken to be a function of the fast space variable X and the slow

space variable x = εX, and we assume that c =¯c(x) + ε

(X), where the

sufﬁx k refers to fractures (f ) or matrix block (m).

Let G(X, Y) be a Green’s function satisfying

∇

G =δ(X −Y) in M, G =0forY ∈∂M,

and suppose that

∇

ψ = 0inM,

ψ = χ on ∂M.

Show that

ψ =



∂M

∂G(X, Y)

∂N

χ(Y)dS(Y),

where

∂G

∂N

=n.∇

G(X, Y), and hence show that

∂ψ

∂N



∂M



∂M

K(X, Y)χ(Y)dS(Y),

where

K(X, Y) =

∂

G(X, Y)

∂N

Now suppose that the ﬂuctuating matrix and fracture concentrations of the

chemical reactant are given by

∇

=Pe



∂ ¯c

∂t

.∇

¯c



+PeΛ ¯c −∇

¯c ≡R

in M,

subject to c

on ∂M, and

∇

−

(1 −φ

)

∂c

∂N



∂M

=Pe



∂ ¯c

∂t

.∇

¯c



+PeΛ ¯c −∇

¯c ≡R

on ∂M,

subject to conditions of periodicity with zero mean.