Overview of stirring and mixing
December 7, 2000

1

Diﬀusion by discontinous movements

In 1827 Robert Brown, observed that suspended pollen grains are in an
uninterrupted and irregular “swarming” motion. Brown was a botanist andat ﬁrst he believed that only organic materials exhibited this agitation. Butvery soon he extended his observations to particles of inorganic material, suchas a ground-up fragment of the Sphinx. Through the nineteenth century there
was a intermittent discussion concerning the cause of this Brownian motion,and in 1877 Delsaux suggested that the impact of molecules on a macroscopic
particle produces observable displacements. In 1905, after nearly a century
of debate, Einstein deﬁnitively explained this phenomenon [6, 7].

1.1

Einstein’s derivation of the diﬀusion equation

Our interest here is in Einstein’s derivation of the diﬀusion equation, which
is very diﬀerent from that of Fourier. We consider one-dimensional Brownian
motion by projecting the location of the particle onto a straight line which
we call the x-axis.
Einstein’s assumptions are the following: (i) the particles move independently of one another; (ii) we observe particle positions at time intervals τ
which are much greater than the time intervals between molecular collisions.
As a result, the motion in one interval is independent of what happened in
the previous interval.
In the interval τ each particle has a random displacement ∆ along thex-axis. The probability density function (PDF) of ∆ is φ(∆). This meansthat if we observe N  1 particles for a time τ then the number of particles
1

A random walk with 200 steps

Figure 1: Simulated Brownian motion using MATLAB; the routine rand is used
to generate a sequence of 200 random displacements.

which are displaced through a distance which lies between ∆ and ∆ + d∆ is
dN = N φ(∆) d∆.

(1)

The PDF φ(∆) does not change from interval to interval, and φ is symmetricand normalised: ∞φ(∆) = φ(−∆) ,φ(∆) d∆ = 1.
(2)
−∞

2

The average of any function of ∆, f (∆), is ∞¯f (∆)φ(∆) d∆ .
f≡

(3)

−∞

In particular, ∆2 is the mean square displacement in a single step.
If the concentration of particles at time t is denoted by c(x, t), then theevolution of c is determined from the master equation:
 ∞c(x, t + τ ) =c(x − ∆, t)φ(∆) d∆ .
(4)−∞

The integral over ∆ is a sum over the prior locations of the particles which
are at x at time t + τ . Thus, the number of particles in the interval (x −∆, x − ∆ + d∆) is c(x − ∆, t)d∆ and φ(∆) is the fraction of these particles
which jump from x − ∆ onto x.
If the concentration c(x, t) changes on a length scale which is much greaterthan the root mean square displacement, then we can approximate the integral equation (4) by the diﬀusion equation. This assumption that c is slowly
varying means that it it is sensible to use a Taylor series expansion


 ∞
∆2
c(x, t) + τ ct (x, t) ≈φ(∆) c(x, t) − ∆cx (x, t) +cxx (x, t) d∆. (5)2−∞Next, using (2), we reduce (5) toct (x, t) ≈ Dcxx (x, t) ,

D≡

∆2
.
2τ

(6)

This is the diﬀusion equation, and D is the diﬀusivity.
The greatness of Einstein’s contribution to this subject is not the derivation above but rather his formula for the diﬀusivity of a macroscopic particle
D=RT,
6πN νa(7)

where R is the gas constant, T the absolute temperature, N the Avogadro
number, ν the coeﬃcient of viscosity and a the radius of the particle. Coincidentally (7) was also discovered by William Sutherland in Australia, alsoin 1905. This relation enabled Perrin to determine Avogadro’s number byobserving Brownian displacements [7].
3The diﬀusion equation is an approximation of the more exact master equation. As we try to design parameterizations of nonlocal mixing processes, inwhich scale separation assumptions are shaky, we should pay more attention
to this history and consider the possibility of using integral equations suchas (4). Notice also that if the Taylor expansion in (5) is continued to higher
order then one will usually (i.e. for most kernels φ) obtain a hyperdiﬀusive
term such as cxxxx .

1.2

The method of moments

As a check on the derivation of (6), we take a diﬀerent approach using themethod of moments. A moment of the concentration is an integral of theform ∞xn c(x, t) dx.
(8)−∞The zeroth moment, n = 0 in (8), is the total number of particles: ∞c(x, t) dx.
N=

(9)

−∞

The ﬁrst and second moments can be interpreted as the center of mass andmoment of inertia of the concentration proﬁle.
We expect that N is constant, and it is educational to verify this conservation law for both the master equation and the diﬀusion equation by“taking the zeroth moment”. Integrating (4) from x = −∞ to x = +∞, andchanging the order of the integrals on the RHS gives ∞ ∞N (t + τ ) =d∆ φ(∆)
dx c(x − ∆, t).
(10)
−∞

−∞

Changing variables to x = x − ∆ in the inner integral, and using (2), gives
the particle conservation law that N (t + τ ) = N (t). The diﬀusive analog ofparticle conservation is easily obtained by integrating the diﬀusion equation(6) from x = −∞ to x = +∞. Provided that Dcx vanishes at x = ±∞(physically, there is no ﬂux of particles from inﬁnity) one immediately ﬁnds
that Nt = 0.
Extending the procedure above to higher
 n moments, we can make a comparison between the exact results for x c dx and the diﬀusive approximations of these same integrals. To take the ﬁrst moment of the diﬀusion4

equation, multiply (6) by x and integrate from x = −∞ to x = +∞. Onceagain, we use integration by parts and assume that terms such as xcx and cvanish as x → ±∞. Thus we ﬁnd that the center of mass is stationaryd ∞xc(x, t) dx = 0.
(11)
dt −∞
The same result can be obtained by taking the ﬁrst moment of the masterequation. The center of mass is stationary because in (2) we assume that the
PDF of displacements is symmetric.
Continuing, we come to the second moment. For the diﬀusion equation
we obtain
 ∞

d ∞2
x c dx = 2D
c dx,
(12)
dt −∞
−∞
where, as before, the terms which fall outside the integration by parts are
zero because of the rapid decay of c as x → ±∞. The student should showthat from the master equation
 ∞
 ∞
 ∞
2
2
x c(x, t + τ ) dx −x c(x, t) dx =∆2 φ(∆) d∆.
(13)
−∞

−∞

−∞

Recalling the deﬁnition of the diﬀusivity in (6), we see that in the limit
τ → 0 the diﬀerence equation in (13) can be approximated by the diﬀerential
equation in (12).
The law in (12), that the mean square displacement of a cloud of particles
grows linearly with time, is often taken to be the deﬁning characteristic ofdiﬀusion. As we will see later, there are dispersive processes which haveother power-laws, such as x2 c dx ∝ t1/2 . These processes are referred to as“anomalous diﬀusion”.

2

Diﬀusion by continuous movements

2.1

Lagrangian time series

In 1922 Taylor [11] analyzed the diﬀusing power of a velocity ﬁeld. The basic
concept here is that of a Lagrangian time series, such as the x-velocity of atagged ﬂuid particle, u(t), as a function of time. This data is Lagrangian (i.e.,5

A time series of Lagrangian velocity

u(t),

2
0
-2
0

0.5

1

1.5

2

2.5

3

3.5

4

t

Figure 2: A time series with a spectral peak.
following a “ﬂoat”), not Eulerian (i.e, obtained from a “current meter”). The
velocity time series might look like ﬁgure 2. Clearly there is some regularity:
evenly spaced maxima and minima are obvious, and we might guess that
there is a wave which is producing oscillatory displacements. At the sametime, the velocity is not completely predictable, and there is no obvious law
by which we can anticipate all details of the future using observations of thepast.
The simplest assumption we can make to analyze the process in ﬁgure
2 is that the velocity is statistically stationary. This means that average
properties of the velocity, such as the mean square velocity, are not changing
with time. In operational terms, the assumption of stationarity means that
if we take nonoverlapping and well-separated subsamples of the time series
in ﬁgure 2 then the statistical properties of the subsamples are identical.
If the time series is long enough we can chop it into N chunks, each oflength T . We deﬁne an ensemble average by considering each of the N chunks
as a single realization of a random process. This procedure introduces the
additional assumptions that there is a decorrelation time τT , and thattime averages are equivalent to ensemble averages. Thinking of dispersion,
Taylor imagined that each chunk was an independent particle, labeled n =1, 2, ..., N , executing continuous movements. “Continuous” in this context
means that the velocity of particle n, un (t), is a relatively smooth functionof time, at least in comparison with the jittery motion in ﬁgure 1.
We denote the position of particle n by xn (t), so that if all the particles

6begin at x = 0 thendxn
= un (t),dt txn (t) = un (t ) dt .

⇒

(14)

0

We use angular brackets to denote the ensemble average. As an exampleof this notation, the average velocity of the N particles is1 
un (t).
u ≡
N n=1
N

(15)

Because of the stationarity assumption, u is independent of time, and we
can refer all displacements relative to the position of the center of mass by
writing x = x− u t and un = un − u . To save decorating all our subsequent
x’s and u’s with primes we now assume that u = 0.

2.2

Taylor’s formula

The simplest measure of dispersion about the center of mass is the meansquare displacement, x2 . We can calculate the rate of change of this quantity by ﬁrst noting that:
 t
dx2n
dx2n
and (14) ⇒(16)= 2xn un ,= 2 un (t)un (t ) dt .
dt
dt
0
We now ensemble average (16). Because of stationarity, u(t)u(t ) depends
only on the time diﬀerence t − t . Thus, we introduce the correlation functionC(t − t ) ≡ u(t)u(t ) ,(17)and, after a change of variables, write the ensemble average of (16) as t
d x2= 2 C(t ) dt .
(18)
dt
0
Equation (18) is Taylor’s formula, which relates the variance in particle displacement x2 to an integral of the Lagrangian velocity autocorrelation function C(t).
In the simplest situations the correlation function C(t) decreases rapidly
to zero as t → ∞ so that the integral in (18) converges. In this case, the7dispersion of the ensemble at large times is characterized by a diﬀusivity
x2 ∼ 2Dt, where the diﬀusivity D is related to the correlation function by: ∞C(t) dt.
(19)
D=
0

In statistical physics, (19) is known as the Green-Kubo formula.
Taylor did not claim that turbulent dispersion was governed by the diﬀusion equation, (6). We will return to this point later. For the moment noticethat (6) is an approximation valid only for suﬃcently long times that the integral in (18) has converged to the constant D. This restriction is related to
Einstein’s assumption that particle positions are observed at time intervals
τ which are much greater than the decorrelation time.

3

Diﬀusion and anomalous diﬀusion

In the previous sections we emphasized that the diﬀusion equation (6) is
only valid on times long compared to the decorrelation time τ , and only if
the concentration c(x, t) varies on length scales greater than the width ofthe density φ(∆). These assumptions of scale separation in both time andspace are often not satisﬁed in real ﬂows. Thus, dispersion experiments overthe last ten years have revealed behaviours which are much richer than thosesuggested by the arguments of Einstein and Taylor. Experiments often showthat the growth of variance is described by a power law
x2 ∝ tξ .

(20)

But sometimes ξ = 1 (diﬀusion), while in other cases ξ = 1. If ξ = 1 theprocess is referred to as anomalous diﬀusion.

3.1

Rayleigh-Bénard convection

As an example of hydrodynamic diﬀusion (ξ = 1) and transient subdiﬀusion(ξ = 2/3) we mention the experiments of Solomon and Gollub [9, 8] on thedispersion of passive scalar (either methylene blue or uranine dye, or small
latex spheres) along a chain of Rayleigh-Bénard convection cells (see ﬁgure3). We refer to the passive scalar generically as “tracer”.
Following the experimental procedure in ﬁgure 3, suppose that all of the
tracer is initially released in a single cell. The main question is: how many8

cells, N (t), have been invaded by tracer at time t? If this dispersive process
is described by diﬀusion then we expect that N (t) ∝ t1/2 . With certain
interesting restrictions, this t1/2 -law is the experimental result.
The Rayleigh-Bénard ﬂow can be approximately described using a twodimensional and incompressible velocity ﬁeld, (u, v), obtained from the streamfunction
ψ = k −1 A sin [k (x + B sin ωt)] W (z) ,(u, v) = (−ψy , ψx ) .

(21)

The parameter A controls the amplitude of the ﬂow, k = 2π/λ is thewavenumber, and W (z) is a function which satisﬁes the no-slip boundary
conditions at z = 0 and z = H. The term B sin ωt is a simple model ofthe lateral oscillation of the roll pattern which results from an instability
which occurs when the convection is driven suﬃciently strongly. Because
the ﬂow in (21) is simple, highly structured and deterministic, this is not an
example of turbulent dispersion. Nonetheless, the experimental results canbe summarized using the notion of an eﬀective diﬀusivity.
The Péclet number is
A
P ≡
,
(22)
kκ
where κ is the molecular diﬀusivity of the tracer, is a nondimensional parameter which measures the importance of molecular diﬀusivity to advection.
The Péclet number can be considered as the ratio of the time it takes amolecule to orbit around a convection cell to the diﬀusion time across a cell.
In the experiments described here, P is large and molecules make many circuits around a convection cell before Brownian motion jostles them througha distance as large as k −1 .
There are two cases which must be carefully distinguished:
Steady rolls The rolls are steady if either ω = 0 or B = 0 in (21). In either
case, tracer can pass from one roll to a neighbour only via molecular
diﬀusion. But, because molecules are advected through a distance k −1 ,the dye is transported
along the array of cells with an eﬀective diﬀusivity Deﬀ ∝ Aκ/k  κ. Because Deﬀ → 0 if κ → 0, the transport
is limited by molecular diﬀusion.
Unsteady rolls If B and ω are both nonzero then advection (rather thanmolecular diﬀusion) can take particles through the time-averaged position of the cell boundaries. In this case, there is the possibility oftransport unlimited by weak molecular diﬀusion.
9

Figure 3: Transport of uranine dye along an array of convection cells with kB =0.12; time (from the top): 1, 2, 4 and 10 periods of oscillation. (Figure courtesyof Tom Solomon [10].)In the unsteady case, Solomon and Gollub show that trajectories of particles computed with the model streamfunction (21) are similar to the patterns observed experimentally. In both the numerics and the experiments,
provided that ωB = 0, the transport of particles along the array of cells (inthe x–direction) is due to chaotic advection in the neighbourhood of the rollboundaries. This process is strikingly shown in ﬁgure 3.
A rough summary of the results is that in both the steady and the unsteady cases the dye spreads via a one-dimensional diﬀusive process, ξ = 1
in (20), with√ a local eﬀective diﬀusivity Deﬀ . The number of invaded cells
is N (t) ∝ Deﬀ t. In the unsteady case Deﬀ is √independent of the molecular
diﬀusivity κ, while in the steady case Deﬀ ∝ κ. The eﬀective diﬀusivity
in the unsteady case is enhanced by 1 to 3 orders of magnitude over theeﬀective diﬀusivity of the steady case (which in turn is much greater than10the molecular diﬀusivity, κ).
The summary in the previous paragraph omits many interesting details.
One of the more important caveats is that the eﬀective diﬀusivity in thesteady case only describes the dispersion process at very long times:
N (t) ∝ t1/2

when t 

1
.
k2κ

(23)

The time 1/k 2 κ is an estimate of the time taken for molecular diﬀusion totransport tracer through a distance of order k −1 , from the edge of a cell tothe center1 . In this long time limit, the evolution of the tracer is slower than
the intracellular diﬀusion time 1/k 2 κ and consequently the concentrationis uniform within each roll. The concentration changes rapidly at diﬀusive
boundary layers (with thickness proportional to κ1/2 ) which are located at
the roll boundaries. The intercellular ﬂux across these boundary layers is
responsible for the spread of the tracer from one roll to the next.
The scenario described above does not have time to become establisheduntil t  1/k 2 κ. When t
1/k 2 κ there is still a signiﬁcant dispersion of
tracer through many cells which is described by the anomalous diﬀusion law
N (t) ∝ t1/3

when t

1
.
k2κ

(24)

The anomalous process above relies on molecular diﬀusion passing tracer
quickly cross the cell boundaries before there has been time to reach thecenter of newly invaded cells [3, 4, 13]. Thus there is a transient regime ofsubdiﬀusion which preceeds the ﬁnal asymptotic diﬀusive law in (23).

3.2

Anomalous diﬀusion in two-dimensional turbulence

Cardoso et al [2] conducted an experimental study of dispersion in a quasitwo-dimensional turbulent ﬂow. The experimental apparatus is a shallow
pan of ﬂuid, 30cm by 30cm, and 3mm deep. The pan is ﬁlled with salty
water and ﬂow is driven electromagnetically (E × B forcing). The forcing isarranged so that the basic ﬂow is a square lattice of 30 × 30 counter-rotating
vortices. This ﬂow is almost two-dimensional because of the large disparity
between the horizontal dimensions (30cm) and the vertical dimension (3mm).
Although the forcing produces a regular array of vortices, this simple pattern is unstable and a two-dimensional turbulent ﬂow emerges. Visualization
1

We assume that the aspect ratio of the cells is of order unity, kH = O(1).

11Figure 4: The trajectory of a single particle shows a sequence of long ﬂights
interrupted by trapping events in which the particle circles around a vortex. The
vortex trapping events are indicated by the arrows. (Scanned from Cardoso etal. [2])of the turbulence, using tracer particles, shows that in the statistically equilibrated state there is a population of vortices whose size is two or three times
the injection scale of the forcing. Each vortex emerges, moves, merges withother vortices, and eventually disappears.
Cardoso et al [2] injected dye into this vortex mess and observed the twodimensional dispersion of the dye in the horizontal plane. To measure thegrowth of the dye blob Cardoso et al deﬁned


2
2
Rm ≡
x + y c(x, y, t) dx dyc(x, y, t) dxdy ,
(25)

12

and

Rg ≡

(x2 + y 2 ) c(x, y, t) dx dyc(x, y, t) dxdy .

(26)

The experimental scaling law is(Rg , Rm ) ∼ t0.32±0.04 .

(27)

The exponent 0.32 = 1/2 indicates anomalous diﬀusion. Because the dispersion is slower than diﬀusion, this is subdiﬀusion.
By examining typical particle trajectories, such as the one in ﬁgure 4,
Cardoso et al explained the subdiﬀusive growth in terms of an interrupted
random walk. Consider a random walker who pauses between steps. The
length of the pause, τ , is a random variable; in the experiment of Cardoso etal the pause is a trapping event in which a molecule is sequestered in the coreof a stationary vortex. If the average duration of a pause is well deﬁned then
one can simply use Einstein’s formula (6) with τ replaced by the average time
between steps. However, if the pausing times are very broadly distributed
then the average duration of a pause may be inﬁnite and consequently thedispersion is subdiﬀusive. More details follow in the next section.

3.3

Random walk with pauses

Consider a random walk in which the walker pauses for a random time τ
between steps. The various τ ’s have a probability density function W (τ )(the waiting time PDF). This PDF is normalised, ∞W (τ ) dτ = 1 ,(28)
0

and the average waiting time spent between steps is ∞τ W (τ ) dτ .
τ̄ =

(29)

0

Motivated by the experiments of Cardoso et al., we entertain the notion thatτ̄ is inﬁnite because the integral in (29) diverges. For example, suppose thatfor large τ , W (τ ) ∼ τ −µ . Then τ̄ = ∞ if µ ≤ 2.
13

However, if we only observe a ﬁnite number of steps, then we do not
sample the entire density W (τ ). Speciﬁcally, suppose that after N steps, wehave experienced pauses of duration τ1 , τ2 , ..., τN . We want to estimate thelikely value of τmax (N ) ≡ max{τ1 , τ2 , · · · τN }. The quantity τmax (N ) is useful
because we can argue that the structure of W (τ ) for τ > τmax (N ) cannot besigniﬁcant for the displacement after N steps.
To determine τmax (N ), we turn to probability theory. Consider a randomvariable θ uniformly distributed in the interval [0, 1]. That is, the PDF ofθ is P (θ) = 1 if 0 < θ < 1 and P (θ) = 0 otherwise. Suppose we take N
samples, θ1 , · · · , θN and deﬁne θmin (N ) ≡ min{θ1 , · · · , θN }. In this simple
case it is plausible that θmin ∼ N −1 as N → ∞.
Now the trick is to use θ to represent τ : we write θ = τ p , and adjust p sothat the power–law tail of W (τ ) ∼ τ −µ corresponds to the simple structureof P (θ) = 1. In fact,
P (θ) = W (τ )

dτ
,
dθ

⇒

1 ∼ τ 1−µ−p ,(30)

or p = 1 − µ. Because the minimum value of θ maps to the maximum valueof τ , it follows that
τmax (N ) ∼ N 1/(µ−1) .

(31)

Now we return to (29) to estimate the eﬀective average pause time after
N pauses: τmax2−µτ̄eﬀ =
τ W (τ ) dτ ∼ τmax
.
(32)
0

It is also plausible that the total time t spent on this random walk is givenbyt ∼ N τ̄eﬀ .

(33)

Combining (31), (32) and (33) yields the following scaling relationships:
N ∼ tµ−1 ,τ̄eﬀ ∼ t2−µ ,τmax ∼ t .

(34)

The ﬁnal relation is worthy of comment: it implies a form of self-similarityof the random walk.
14

To conclude, the total displacement of our random walk is proportional√
to N . But, with the random pauses, the scaling against time has beenaltered to√
RMS displacement ∝ N ∼ t(µ−1)/2 .
(35)
This theory can be used to interpret the experiment of Cardoso et al: because
the RMS displacement grows as t1/3 it follows that µ ≈ 5/3. Cardoso et al.
successfully tested this prediction by measuring the PDF of trapping times
inside vortices.

4

Stirring and mixing

4.1

Coﬀee and cream

Appealing to the everyday experience of mixing cream into coﬀee, Eckart[5] argued that the homogeneization of two ﬂuids occurs in three stages.
The distinction between the stages is the value of the concentration gradient
averaged over the domain.
Initial : there are distinct interfaces separating globules of cream and coﬀee.
Within each globule, the concentration of cream is nearly constantand the concentration gradient is close to zero. There is a very large
concentration gradient between regions of coﬀee and cream. But theinterfaces between coﬀee and cream are small in number and not ofgreat area, so the average gradient in the coﬀee mug is small.
Stirring : the cream is mechanically swirled and folded, and molecular diffusion is unimportant. During this second stage the concentration gradients increase.
Mixing : the gradients suddenly disappear and the ﬂuid becomes homogeneous; molecular diﬀusion is responsible for the sudden mixing.
In a chemical reaction, molecules of diﬀerent species must come into contact
for the reaction to occur. Thus, when the species are initially separated, thereaction will not begin until the ﬁnal mixing stage is reached. In this sense
there is an important distinction between coarse-grained homogeneization,
occuring solely as a result of stirring, and mixing at the molecular scale.
15

t=1

t=2

t=4

t=8

t=16

t=32

Figure 5: Solution of ct + (1 − r2 )cθ = (8 × 10−4 )∇2 c. The initial condition isc(x, y, 0) = x.

To illustrate these concepts ﬁgure 5 shows a solution of the advection
diﬀusion equation
ct + (1 − r2 )cθ = κ∇2 c,c(r, θ, 0) = r sin θ(36)where κ = 8 × 10−4 . A particle at a distance r from the origin completesa rotation in a time 2π/(1 − r2 ). Thus particles at smaller values of r will
overtake particles at larger values of r and so the concentration is twisted
into spirals by diﬀerential advection (stirring).
The increase in gradient during the stirring phase is evident in the ﬁgure.
But at approximately t = 16, mixing starts to dominate, and diﬀusion rapidly
reduces the average gradient. From the initial condition, an estimate of thetime it would take unassisted diﬀusion to homogenise the ﬂuid is TD ∼ 1/κ =1250. It is only through the initial process of stirring that the concentration
16

gradient is ampliﬁed or, alternatively, that the spirals are stretched out sothat small diﬀusion homogenizes the tracer at t = 32
TD .

4.2

A straining ﬂow

A simple example of a two-dimensional ﬂow which ampliﬁes concentration
gradients is the hyperbolic strain shown in ﬁgure 6. The streamfunction isψ = −αxy and so the advection diﬀusion equation is
ct + αxcx − αycy = κ∇2 c .

(37)

Notice the dimensions here: α−1 has dimensions “time” and κ has dimensions(length)2 /(time). From these two quantities we can build a combination withthe dimensions of (length):
κ
.
α

(≡

(38)

The length ( will appear prominently in the sequel.
We begin our discussion of hyperbolic strain by obtaining a solution inwhich c is independent of both x and t. In this special case the solution of(37) is√y2cy = A exp − 2 , c(x, ±∞, t) = ± 2πA(.
(39)
2l
The concentration proﬁle is the error function shown in ﬁgure 7. The solution√ shows the steady state balance between advection and diﬀusion: with
2πA( = 1, the concentration c changes smoothly between c = +1 as
y → +∞ to c = −1 as y → −∞. The transition occurs in a front ofwidth (.
We can give an intuitive discussion of how the steady state proﬁle inﬁgure 7 is established as the solution of an initial value problem. Suppose
we had started with the initial condition such as c(x, y, 0) = sgn(y) in which
the transition between c = −1 and c = +1 occurs in a distance much√lessthan (. Then the discontinuity in c initially diﬀuses freely, growing like κt.
Once the width of the front becomes comparable to (, that is when√κt ∼ (,⇒t ∼ α−1 ,(40)
17Figure 6: The straining ﬂow described by the streamfunction ψ = −αxy. The
ﬁgure shows how a circular patch of tracer is stretched out along the x–axis by thehyperbolic strain. If κ = 0 the major axis of the ellipse grows as exp(αt) and theminor axis reduces as exp(−αt) so that the area remains constant.

1
0.5

The erf solution

0
-0.5
-1
-2

-1.5

-1

-0.5

0

0.5

1

1.5

2

y

Figure 7: The time independent error function solution to equation (37).

18the spread is arrested and the steady state in ﬁgure 7 is established.
On the other hand, we can also consider an initial condition in which thetransition between c = 1 and c = −1 occurs on a scale L0  (. In this casethe front is initially compressed by the hyperbolic strain so that the widthis reduced exponentially, L = L0 exp(−αt). Because L0  ( the diﬀusion is
unimportant until the exponential reduction in scale reaches (. That is,
L0 e−αt ∼ ( ,t ∼ α−1 ln(L0 /().

⇒

(41)

These considerations illustrate fundamental importance of ( as the scale onwhich advection and diﬀusion come into balance.

4.3

Lagrangian coordinates: a simple example

The hyperbolic strain also provides a painless illustration of some mathematical techniques which can be used in more complicated problems. We begin
by considering the solution of (37) with κ = 0. With no diﬀusion c is tied to
ﬂuid particles. The position of a ﬂuid particle is related to its initial position(a, b), by solving the diﬀerential equations(ẋ, ẏ) = α(x, −y),⇒(x, y) = eαt a, e−αt b .

(42)

The solution of (37) can now be obtained by arguing that the particle which
is at the point (x, y) at time t began at (a, b) = (exp(−αt)x, exp(αt)y) at
t = 0. Because the a particle carries the concentration it follows that thesolution of (37) as an initial value problem isc(x, y, t) = c0 [exp(−αt)x, exp(αt)y] ,(43)

where c0 (x, y) is the initial condition. The philosophy of this method is that
we care where ﬂuid particles come from, but not where they are going to.
The solution above seems to rely crucially on the restriction that κ = 0.
But now look what happens if we use the Lagrangian coordinates (a, b) in(42) as new independent variables in (37). As an accounting device, it iscomforting to deﬁne τ = t and consider that ∂τ as the time derivative with(a, b) ﬁxed. Thus the transformation rules are




∂a ∂a
∂b ∂b
(∂x , ∂y ) =
(44)
,
∂a +
,
∂b = e−ατ ∂a , eατ ∂b .
∂x ∂y
∂x ∂y
19

and
∂t =

∂τ
∂a
∂b
∂τ +
∂a + ∂b = ∂τ − αa∂a + αb∂b .
∂t
∂t
∂t

(45)

The punchline is that∂t + αx∂x − αy∂y = ∂τ ,(46)

which shows that the change to a Lagrangian description makes the convective derivative trivial.
Substituting the transformations above into (37) gives:
ct = κe−2αt caa + κe2αt cbb .

(47)

Naturally, if κ = 0, we recover our earlier solution in (43). But even if κ = 0
it is often easier to solve (47) than the Eulerian form in (37). For example,
Fourier transforming (47), with (∂a , ∂b ) → i(p, q), gives a simple ordinary
diﬀerential equation in time.
It is instructive to use the method above to solve (37) with the initialconditionc(x, y, 0) = δ(x)δ(y) .

(48)

Physically, this is a spot of dye released in a straining ﬂow. √When αt1
the spot spreads diﬀusively, with a diameter which grows as κt. However
when αt ∼ 1 the diameter of the spot becomes comparable to ( ≡ κ/α,and then the spot stops expanding against the compressive direction of thestrain. However the spot continues to stretch along the extensive direction.
Thus, when αt > 1, the spot becomes a ﬁlament with an equilibrium widthof order ( and an exponentially growing length. These intuitive arguments
are supported by the exact solution:


1
x2
y2
c(x, y, t) =exp − 2 − 2 ,
(49)
4πf g
4f
4g
where f (t) and g(t) are
f2 ≡

κ 2αt
(e − 1) ,
2α

g2 ≡

κ
1 − e−2αt .
2αNotice that the peak concentration ultimately decreases like e−2αt .
20

(50)

12

〈 ∇ c ⋅ ∇ c 〉 using
kl=1, 1/2, 1/4, 1/8

10

8

6

4

2

0

0

0.5

1

1.5

αt

2

2.5

3

3.5

Figure 8: The mean square of the concentration gradient. If k

1 then theconcentration gradient grows until t = t∗ in (54) and then decreases precipitiously.
If k ≥ 1/2 then diﬀusion always overpowers strain and the mean square gradientdecreases monotonically to zero.

4.4

An example of sudden mixing

As a ﬁnal look at the hyperbolic straining ﬂow, we note that a solution of(37) isc(x, y, t) = A(t) cos(ke−αt x) cos(keαt y),

(51)



A(t) = exp −(2 k 2 sinh 2αt .

(52)

where

One route to this exact solution is to look for separable solutions of (47),and then transform back to the Eulerian coordinates (e.g., Young, Rhines &Garrett,1982).
The mean value of the square of the concentration gradient varies with21

time as:{∇c·∇c} =k2cosh(2αt) exp −2(2 k 2 sinh(2αt) ,2(53)where {} denotes an average over a large area. {∇c·∇c} is plotted in Figure 8 for various values of k(. Recalling Eckart’s description of stirring as
increasing the concentration gradient, and mixing as decreasing the concentration gradient, we can see the transition between the two phases occurs at
the peaks of the various curves. If k(
1, then the time it takes to reachthis peak is given by t∗ , whereαt∗ ∼ − ln(k().

(54)

Once again, this is the time taken for the exponential factor e−αt to reduce
initial length of the tracer ﬁeld, k −1 , down to the length ( on which strain
and diﬀusion balance.

4.5

A Welander scrapbook

Stirring was beautifully illustrated in a 1955 paper of Welander’s [12]. This
paper is notable also because of its discussion of the importance of coarsegrained averages. Figures 9, 10 and 11 reproduced from Welander (1955)
show that simple velocity ﬁelds produce spectacular distortion of passive
scalars.
In ﬁgures 9, 10 and 11, some dimensions of the scalar blob are stretchedout while other dimensions are contracted. Batchelor (1952) [1] argued thatin turbulent ﬂows random stretching results in an exponential growth of
the separation between two initially adjacent ﬂuid elements. That is, if we
consider two material elements separated by a distance s0 which is much lessthan the scale of the velocity ﬁeld, then Batchelor argues that the separation
grows ass ∼ s0 eγt .

(55)

The time-scale γ −1 is analogous to α−1 in (37), though in ﬁgures 10 and
11 the exponential straining is driven by a random and unsteady velocity,rather than the simple hyperbolic ﬁeld in ﬁgure 6. Note particularly that theexponential law in (55) is valid until the separation s(t) becomes comparable
to the length scale over which the velocity varies.
22

Figure 9: Welander’s numerical solution illustrating diﬀerential advection by asimple velocity ﬁeld. A checkerboard pattern is deformed by a quasigeostrophic
barotropic solution which models atmospheric ﬂow at the 500mb level. The initial
streamline pattern is shown at the top and the subsequent ﬁgures are at 6 hours,
12 hours, 24 hours and 36 hours, respectively. Notice that each square of thecheckerboard maintains constant area as it deforms.

23

Figure 10: Welander’s experimental illustration of the deformation of a small
coloured square element of a ﬂuid surface. To suppress three dimensional turbulence, a vessel of water is brought into solid body rotation. A ﬂoating ﬁlm of
butanol is divided into square elements by means of a metal grid and one of thesesquare elements is coloured with methyl-red. The solid body rotation is disturbed
by stirring the water and the grid is removed.

24

Figure 11: Further deformation of the butanol square in ﬁgure 10. According to
Batchelor [1], the length of the ﬁlaments increases exponentially with t.

25

References[1] G.K. Batchelor. The eﬀect of turbulence on material lines and surfaces.
Proc. Roy. Soc. London A, 213:349–366, 1952.
[2] O. Cardoso, B. Gluckmann, O. Parcollet, and P. Tabeling. Dispersion ina quasi-two-dimensional turbulent ﬂow: an experimental study. Phys.
Fluids, 8(1):209–214, January 1996.
[3] O. Cardoso and P. Tabeling. Anomalous diﬀusion in a linear array ofvortices. Europhys. Lett., 7(3):225–230, 1988.
[4] O. Cardoso and P. Tabeling. Anomalous diﬀusion in a linear system ofvortices. Euro. J. Mech. B/Fluids, 8(6):459–470, 1989.
[5] C. Eckart. An analysis of stirring and mixing processes in incompressible
ﬂuids. J. Mar. Res., 7:265–275, January 1948.
[6] A. Einstein. Investigations on the theory of the Brownian Movement.
Dover, New York, 1956.
[7] A. Pais. ‘Subtle is the Lord...’. Oxford University Press, Oxford and
New York, 1982.
[8] T.H. Solomon and J.P. Gollub. Chaotic particle transport in time dependent Rayleigh-Bénard convection. Phys. Rev. A, 38:6280, 1988.
[9] T.H. Solomon and J.P. Gollub. Passive transport in steady RayleighBénard convection. Phys. Fluids A, 31:1372–???, 1988.
[10] T.H. Solomon, S. Tomas, and J.L. Warner. Chaotic mixing of immiscible
impurities in a two-dimensional ﬂow. Phys. Fluids A, 10:342–350, 1998.
[11] G.I. Taylor. Diﬀusion by continuous movements. Proc. London Math.
Soc., 20:196–212, 1921.
[12] P. Welander. Studies on the general development of motion in a twodimensional, ideal ﬂuid. Tellus, 7:141–156, 1955.
[13] W.R. Young, A. Pumir, and Y. Pomeau. Anomalous diﬀusion of tracer
in convection rolls. Phys. Fluids A, 1:462–469, 1989.

26

[14] W.R. Young, P.B. Rhines, and C. Garrett. Shear-ﬂow dispersion, internal waves and horizontal mixing in the ocean. J. Phys. Ocean., 12:515–
527, 1982.

27