J.
Fl-uid
Mech.
(1987),
vol.
180,
pp.
387---403
Printed
in
Great
Britain
387
A
non-local
description
of
advection-diffusion
with
application
to
dispersion
in
porous
media
By
DONALD
L.
KOCHt
AND
JOHN
F.
BRADY
Chemical
Engineering,
California
Institute
of
Technology,
Pasadena,
CA
91125,
USA
(Received
14
February
1986
and
in
revised
form
6 October
1986)
When
the
lengthscales
and
timescales
on
which
a
transport
process
occur
are
not
much
larger
than
the
scales
of
variations
in
the
velocity
field
experienced
by
a
tracer
particle,
a description
of
the
transport
in
terms
of
a local,
averaged
macroscale
version
of
Fick's
law
is
not
applicable.
Here,
a
non-local
transport
theory
is
developed
in
which
the
average
mass
flux
is
not
simply
proportional
to
the
average
local
concentration
gradient,
but
is
given
by
a convolution
integral
over
space
and
time
of
the
average
concentration
gradient
times
a spatial-
and
temporal-wavelength-
dependent
diffusivity.
The
non-local
theory
is
applied
to
the
transport
of
a passive
tracer
in
the
advective
field
that
arises
in
the
bulk
fluid
of
a porous
medium,
and
the
complete
residence-time
distribution-
space-time
response
to
a
unit
source
input
-
of
the
tracer
is determined.
It
is
also
shown
how
the
method
of
moments
may
be
simply
recovered
as
a special
case
of
the
non-local
theory.
While
developed
in
the
context
of
and
applied
to
tracer
dispersion
in
porous
media,
the
non-local
theory
presented
here
is
applicable
to
the
general
problem
of
determining
the
average
transport
behaviour
in
advection-diffusion-type
systems
in
which
spatial
and
temporal
variations
are
occurring
on
scales
comparable
with
the
scale
of
interest.
1.
Introduction
Many
important
physical
processes
involve
the
transport
of
a
passive
scalar
quantity
in
the
presence
of
a complex
advective
field
that
fluctuates
on
lengthscales
and
timescales
comparable
to
those
on
which
the
transport
process
occurs.
For
example,
the
initial
dissipation
of
a
pollutant
in
the
atmosphere,
in
rivers,
or
in
ground-water
flows
occurs
on
a lengthscale
comparable
to
or
smaller
than
the
scale
of
the
advective
motion.
The
mean
residence
time
of
a species
in
a reservoir
or
chemical
reactor
may
be
smaller
than
the
time
required
for
the
tracer
to
experience
the
full
range
of
velocity
fluctuations
present
in
such
a heterogeneous
medium.
The
oscillatory
flows
in
the
pulmonary
and
cardiovascular
systems
as
well
as
in
tidal
estuaries
give
rise
to
time-dependent
dispersive
effects
that
depend
fundamentally
on
the
inherent
unsteadiness.
When
the
timescales
and
lengthscales
on
which
a
transport
process
occurs
are
not
much
larger
than
the
scales
of
variations
in
the
velocity
field
experienced
by
a
tracer
particle,
a description
of
the
transport
in
terms
of
an
averaged,
macroscale
version
of
Fick's
law
is
not
applicable.
In
this
paper
a more
general,
'non-local'
description
of
transport
is
derived.
In
the
non-local
theory
the
average
mass
flux
is
not
simply
proportional
to
the
average
local
concentration
gradient,
but
is given
by
a convolu-
tion
integral
over
space
and
all
previous
times
of
the
average
concentration
gradient
t
Present
address:
School
of
Chemical
Engineering,
Cornell
University,
Ithaca,
NY
14853,
USA.
388
D.
L.
Koch
and
J.
F.
Brady
times
a temporal-
and
spatial-wavelength-dependent
diffusivity.
The
average
mass
conservation
equation
then
is
a second-order
integrodifferential
equation,
which
takes
on
a
particularly
simple
form
upon
Fourier
transformation
in
time
and
space.
Wavelength-dependent
properties
have
long
been
used
to
treat
non-local
effects
in
physics.
For
example,
the
electric
displacement
in
a dielectric
may
be
related
to
the
electric
field
through
a dielectric
constant,
which
is
in
general
a function
of
temporal
and
spatial
wavelengths
(Jackson
1975;
Landau,
Lifshitz
&
Pitaevskif
1984).
Indeed,
the
overall
linearity
of
Maxwell's
equations
of
electromagnetism
require
that
the
dielectric
displacement
be
a linear functional
of
the
imposed
field,
and
the
non-local
wavelength-dependent
representation
is
an
appropriate
manner
of
expressing
this
functionality.
In
situations
much
closer
to
ours,
non-local
theories
have
been
developed
to
study
wave
propagation
in
random
media
(Keller
1962;
Frisch
1968),
and
a general
approach
to
such
non-local
problems
for
steady
second-order
linear
stochastic
differential
equations
has
been
developed
by
Beran
&
McCoy
(1970).
While
such
non-local
approaches
are
widespread
in
the
physics
and
applied
physics
literature,
they
do
not
seem
to
have
been
applied
to
advection-diffusion
problems.
Previous
analyses
of
transient
effects
in
advectively
enhanced
dispersion
for
Taylor
dispersion
in
a
tube
(Gill
&
Sankarasubramanian
1970;
Chatwin
1970;
Smith
1981),
for
generalized
Taylor
dispersion
theory
(Pagitsas,
Nadim
&
Brenner
1986)
and
for
an
ad
hoc,
one-dimensional
model
of
a
packed
bed
with
adsorption
(Reis
et
al.
1979)
have
been
based
on
an
expansion
of
the
concentration
field
in
moments.
The
resulting
expression
relates
the
mass
flux
to
an
infinite
expansion
in
time-
and
space-derivatives
of
the
concentration
field.
This
procedure
requires
the
determination
of
an
infinite
set
of
phenomenological
coefficients,
which
are
in
general
tensors.
It
is
unclear
how
one
would
obtain
the
infinite
set
of
boundary
conditions
required
to
solve
a
boundary-value
problem
involving
an
infinite-dimensional
mass
conservation
equa-
tion.
In
addition,
it
is difficult
to
calculate
more
than
the
first
few
phenomenological
coefficients,
so
the
method
is
not
applicable
at
short
times
and
short
spatial
wavelengths
when
the
non-Fickian
effects
are
expected
to
be
most
pronounced.
(Smith's
analysis
goes
further
than
a
moments
expansion
by
writing
an
effective
diffusion
coefficient
in
terms
of
a time-delay
integral,
but
the
limitations
of
this
approach
for
short
times
are
still
present.)
In
this
paper
we
recognize
the
inherent
coupling
between
spatial
and
temporal
variations
and
follow
the
non-local
approaches
in
generalizing
the
conventional
local
Fick's
law
expression
for
mass
conservation.
This
approach
represents
a generaliz-
ation
of
the
Eulerian
description
of
transport,
as
opposed
to
the
moments
or
Lagrangian
approaches
discussed
above,
and
will
be
seen
to
have
significant
advantages.
The
non-local
conservation
equation
derived
in
§2
gives
the
ensemble-average
field
of
a scalar
quantity
in
any
advection-diffusion
problem
for
which
the
quantity's
velocity
v
is
known
(at
least
in
a
statistical
sense).
This
conservation
equation
is
applicable
regardless
of
the
lengthscales
and
timescales
over
which
the
scalar
quantity
varies,
and
it
allows
a complete
determination
of
the
concentration
field
for
all
space-time
positions.
Thus,
the
full
residence-time
distribution
of
a
tracer
being
transported
through
the
advective
field
can
be
predicted.
In
§2
we
also
discuss
how
the
ensemble-average
approach
taken
here
reduces
to
an
average
over
a
unit
cell
for
spatially
periodic
systems
and
can
be
modified
to
treat
problems
whose
structure
is
formally
similar
to
Taylor
dispersion
in
a
tube.
We
also
show
how
the
integro-
differential
formulation
can
be
locally
expanded
to
recover
the
method
of
moments
used
in
previous
studies.
A
non-local
description
of
advection-diffusion
389
In
§§3
and
4 we
demonstrate
the
application
of
this
general
theory
to
the
transport
of
a solute
in
the
advective
field
that
arises
in
the
bulk
fluid
of
a porous
medium
consisting
of
fixed
spheres
whose
positions
are
randomly
distributed.
In
a previous
paper
(Koch
&
Brady
1985)
we
showed
that
one
can
capture
the
correct
qualitative
nature
of
this
tr~nsport
process
through
an
asymptotic
analysis
in
small
solids
volume
fraction.
In
this
limit
the
effect
of
the
particles
on
the
fluid
velocity
may
be
approximated
as
that
due
to
point
forces
in
a
medium
that
obeys
Brinkman's
equations
of
motion.
The
residence-time
distributions
(RTDs)-
the
concentration
field
resulting
from
a
point
source
impulse-
predicted
here
show
the
complete
spatial
and
temporal
evolution
of
the
transport
process
and
the
approach
to
the
long-time
limit.
In
the
long-time
limit
(long
times
and
large
lengthscales),
we
recover,
of
course,
the
classical
local
Fickian
behaviour.
The
velocity
field
considered
here
is
a simple
example
of
a stochastic
velocity
field,
and
as
such
it
may
be
viewed
as
a
crude
model
of
turbulence.
The
general
theory
presented
here
has
much
broader
application
than
determining
the
RTDs
resulting
from
the
fluid-induced
dispersion
-
the
so-called
mechanical
dispersion-
in
a porous
medium.
In
particular,
the
dispersion
in
porous
media
caused
by
the
non-mechanical
mechanisms
associated
with
the
retention
of
the
solute
species
in
regions
of
zero
velocity
(Koch
&
Brady
1985)
has
a pronounced
affect
on
the
RTDs
and
the
approach
to
the
long-time
limit.
These
issues
are
addressed
in
a
separate
work
(Koch
1986;
Koch
&
Brady
1987a).
The
theory
presented
here
is
also
the
natural
starting
point
for
examining
transport
in
systems
which
are
not
spatially
homo-
geneous.
An
application
to
dispersion
in
heterogeneous
porous
media
is
given
in
Koch
&
Brady
(1987b).
2.
The
non-local
transport
equations
In
this
section
we
shall
follow
the
non-local
approach
used
in
other
physical
settings
(Keller
1962;
Frisch
1968;
Beran
&
McCoy
1970)
and
derive
the
general
averaged
conservation
equation
for
a scalar
quantity
c
(such
as
a solute
concentration)
that
satisfies,
at
each
point
in
the
medium,
the
detailed
conservation
equation
where
oc
-+V·q
=
S
Cit
'
D
q
=
vc
=
uc-
MV(Mc)
(2.1a)
(2.1b)
is
the
flux,
u
is
the
sum
of
the
fluid
velocity
and
the
velocity
resulting
from
any
external
force
that
may
be
acting
on
the
solute,
Dis
the
molecular
diffusivity,
and
M
is
the
solubility
(or
solute's
activity
coefficient).
The
sourceS
of
solute
is
taken
to
be
independent
of
concentration
and
includes
the
solute
injected
into
the
bed
initially.
In
general
the
diffusivity
and
solubility
are
functions
of
position
and
are
different,
for
example,
in
the
various
phases
of
a multiphase
medium.
Thus,
D
and
Mare
both
generalized
functions,
which
may
undergo
a
step
change
at
interphase
boundaries.
The
writing
of
the
contribution
to
the
flux
from
the
effects
of
diffusion
and
solubility
as
(D/
M)
V(Mc)
has
been
done
for
convenience
to
remove
any
specific
reference
to
boundary
conditions.
Thus,
(2.1)
applies
throughout
the
entire
medium.
In
the
application
to
dispersion
in
porous
media
consisting
of
a fixed
solid
phase
and
390
D.
L.
Koch
and
J.
F.
Brady
a mobile
fluid
phase,
the
requirement
that
this
term
be
non-singular
gives
continuity
of
mass
flux
and
solute
equilibrium
at
the
interphase
boundaries.
A
step
change
in
solubility
Mat
an
interphase
boundary
necessitates
a
step
change
inc
in
order
that
the
solute
activity
Me
is continuous
and
V(Mc)
is non-singular;
this
step
change
in
c
is
the
interphase
partitioning
of
the
solute
at
local
interphase
equilibrium.
A
step
change
in
D
/
M
due
to
different
molecular
diffusivities
in
the
two
phases
in
turn
requires
a
step
change
in
V(Mc),
so
that
the
flux
(D/
M)
V(Mc)
is
continuous
at
the
boundary
and
its
divergence
is non-singular.
Thus,
the
solute velocity
v
defined
by
(2.1
b)
must
be
interpreted
as
a differential
operator,
which
includes
the
effects
of
convection,
molecular
diffusion,
and
the
interphase
boundary
conditions
on
the
solute
transport.
We
are
interested
in
cases
in
which
the
spatial
variations
in
the
velocity
field
are
complex
and/or
the
precise
position
of
the
source
relative
to
the
velocity
field
is
uncertain,
so
that
the
solute
velocity
is
best
described
stochastically.
Thus,
we
define
an
ensemble
average
(
) as
the
average
over
an
ensemble
of
velocity
fields.
The
precise
specification
of
this
ensemble
necessarily
depends
on
the
application.
For
example,
in
a spatially
periodic
medium
this
ensemble
average
reduces
to
an
average
over
a
unit
cell
(Koch
1986).
In
applying
(2.1)
and
(2.2)
below
to
Taylor
dispersion
in
a
channel
or
tube,
the
ensemble
average
would
be
replaced
by
an
average
over
the
channel
or
tube
cross-section.
In
this
context
it
may
also
prove
more
convenient
to
explicitly
use
boundary
conditions
at
the
tube
wall
rather
than
incorporate
them
in
the
differential
equation
as
we
have
done
here.
This
requires
no
change
in
the
development
below,
because
the
interpretation
of
the
solute
velocity
vas
an
operator
implies
the
imposition
of
boundary
conditions.
Averaging
(2.1a)
yields
a~~>+
v·
(q)
=
(S).
(2.2)
In
order
to
obtain
an
equation
relating
the
average
concentration
to
the
source,
we
must
express
the
average
flux
(q)
=
(vc)
in
terms
of
the
average
or
'macroscale'
properties
of
the
medium.
To
this
end,
we
write
the
concentration
as
c
=
Peq(c)+c',
(2.3)
where
Peq•
the
equilibrium
distribution
of
the
solute
in
the
absence
of
an
average
concentration
gradient,
is defined
by
V·(vPeq)
=
v{
uPeq-
~
V(MPeq)J
=
0,
(Peq)
=
1,
(2.4a)
(2.4b)
and
c'
is
the
deviation
of
the
concentration
from
local
equilibrium.
In
a two-phase
medium
(fluid-solid)
the
solubility
requirement
for
local,
microscale
thermodynamic
equilibrium
at
the
interphase
boundaries
takes
the
form
i.e.
M
=
1
in
the
fluid
phase
and
m
in
the
solid
phase.
Peq
gives
the
partitioning
of
the
solute
between
the
phases
on
the
macroscopic
scale,
and
in
the
case
considered
here
of
a solenoidal
fluid
velocity
field
u
in
a
multi
phase
medium,
Peq
is
constant
in
A
non-local
description
of
advection-diffusion
391
each
phase,
but
undergoes
a
step
change
at
interphase
boundaries,
such
that
MPeq
is constant,
i.e.
{
1-cp!cpm
1
peq
=
1
m(1-cp+cpm-
1
)
in
the
fluid,
in
the
solid,
where
cp
is
the
volume
fraction
of
the
solid
phase.
If
the
fluid
were
compressible,
the
fluid
motion
could
cause
the
equilibrium
distribution
to
vary
within
the
fluid
phase.
Substituting
(2.3)
in
the
equation
for
the
average
mass
flux
gives
(q)
=
(vPeq)
(c)+
(vc').
(2.5)
In
interpreting
(2.5),
we
must
recall
that
the
solute
velocity
vis
an
operator.
Thus,
the
first
term
on
the
right-hand
side
of
(2.5),
obtained
by
operating
on
~q(c)
with
the
operator
v
and
then
averaging,
is
(uPeq)(c)-(
(~)V(MPeq(c)))
=
V(c)-(DPeq)
V(c),
where
(2.6)
is defined
to
be
the
solute
average
velocity
at
equilibrium.
Thus,
using
(2.3)
and
(2.4)
we
may
rewrite
(2.5)
as
(q)
=
V(c)-(DPeq)V(c)+(v'c'),
(2.7)
where
v'
=
v-
V.
An
equation
for
the
concentration
disturbance
may
be
obtained
by
substituting
(2.3)
into
(2.1)
and
using
(2.4)
to
give
QC
1
I
o(c)
at+V·(vc)
=
-vPeq·V(c)-PeqTt+Sd(S),
(2.8)
where
S
=
Sd(S),
i.e.
Sd
is
the
distribution
of
the
source
relative
to
its
average
in
the
microstructure.
The
solution
of
(2.8)
is
c'
=
J
dx
1
foo
dt
1
P(x-x
1
,t-t
1
){
-v'~q·V(c)-~q[a~~>
+
V·V(c)
]+Sd(s)},
(2.9a)
where
oP
oP
[D
J
at+
V·
(vP)
=at+
V·
(uP)-
V·
M
V(MP)
=
8(x-x
1
)
8(t-t
1
).
(2.9b)
Here
P,
the
Green
function
for
the
left-hand
side
of
(2.8),
is
the
transition
probability:
the
probability
of
finding
a
tracer
particle
at
x
at
time
t
given
that
it
was
at
x
1
at
time
t
1
.
Using
the
solution
(2.9)
for
the
concentration
disturbance
the
average
mass
flux
(2.7)
may
be
written
as
(q)
=
V(c)-
Jdx
1
r
00
dt
1
{D(x-xvt-t
1
)'V(c(x
1
,t
1
))
[
o(c(x
1
,
t
1
))
V·
(
>]
(S
>}
+w(x-x1,t-tl)
otl
+
vl
c(xl,tl)
-a(x-xl't-tl)
(xl'tl)
'
(2.10a)