General rights
Copyright and moral rights for the publications made accessible in the public portal are retained by the authors and/or other copyright
owners and it is a condition of accessing publications that users recognise and abide by the legal requirements associated with these rights.
Users may download and print one copy of any publication from the public portal for the purpose of private study or research.
You may not further distribute the material or use it for any profit-making activity or commercial gain
You may freely distribute the URL identifying the publication in the public portal
If you believe that this document breaches copyright please contact us providing details, and we will remove access to the work immediately
and investigate your claim.
Downloaded from orbit.dtu.dk on: Aug 10, 2022
Pattern dynamics of vortex ripples in sand: Nonlinear modeling and experimental
validation
Andersen, Ken Haste; Abel, M.; Krug, J.; Ellegaard, C.; Søndergaard, L. R.; Udesen, J.
Published in:
Physical Review Letters
Link to article, DOI:
10.1103/PhysRevLett.88.234302
Publication date:
2002
Document Version
Publisher's PDF, also known as Version of record
Link back to DTU Orbit
Citation (APA):
Andersen, K. H., Abel, M., Krug, J., Ellegaard, C., Søndergaard, L. R., & Udesen, J. (2002). Pattern dynamics of
vortex ripples in sand: Nonlinear modeling and experimental validation. Physical Review Letters, 88(23),
234302. https://doi.org/10.1103/PhysRevLett.88.234302
VOLUME 88, N
UMBER 23 PHYSICAL REVIEW LETTERS 10J
UNE 2002
Pattern Dynamics of Vortex Ripples in Sand: Nonlinear Modeling and Experimental Validation
K. H. Andersen,
1
M. Abel,
2
J. Krug,
3
C. Ellegaard,
4
L. R. Søndergaard,
4,5
and J. Udesen
4,5
1
Department of Mechanical Engineering, Technical University of Denmark, DK-2800 Kgs. Lyngby, Denmark
2
Institut für Physik, Universität Potsdam, D-14415 Potsdam, Germany
3
Fachbereich Physik, Universität Essen, D-45117 Essen, Germany
4
Niels Bohr Institute, University of Copenhagen, DK-2100 Copenhagen Ø, Denmark
5
Department of Mathematics and Physics, University of Roskilde, Box 260, DK-4000 Roskilde, Denmark
(Received 30 January 2002; published 22 May 2002)
Vortex ripples in sand are studied experimentally in a one-dimensional setup with periodic boundary
conditions. The nonlinear evolution, far from the onset of instability, is analyzed in the framework
of a simple model developed for homogeneous patterns. The interaction function describing the mass
transport between neighboring ripples is extracted from experimental runs using a recently proposed
method for data analysis, and the predictions of the model are compared to the experiment. An analytic
explanation of the wavelength selection mechanism in the model is provided, and the width of the stable
band of ripples is measured.
DOI: 10.1103/PhysRevLett.88.234302 PACS numbers: 45.70.Qj, 47.20.Lz, 47.54. +r
Ever since the establishment of a conceptual frame-
work for pattern formation [1], the description of patterns
formed in sand by the flow of wind or water has posed
a challenge to the community [2–7]. Despite diverse ef-
forts including coupled map models for ripples and dunes
in air [2], stochastic models for ripples in air [3], or con-
tinuum equations based on the symmetries of the problem
[4], theoretical understanding has remained sparse. For ex-
ample, all models display a coarsening of the ripple/dune
pattern, but the coarsening does not terminate at a final
selected wavelength, as is frequently observed in nature.
Furthermore, the models are heuristic, and it is not pos-
sible to make a quantitative comparison with experiments.
Here we study
vortex ripples [8], which are created by
an oscillatory water flow, such as that generated near the
sand bed by a surface wave. Vortex ripples have attracted
attention as an example of a nonlinear pattern forming
system with a strongly subcritical first bifurcation [6,7,9],
which cannot be described by conventional methods like
amplitude equations [9]. As most other sand patterns, they
display coarsening and saturation at a finite wavelength.
The approach pursued in this Letter combines a simple
model for the fully developed pattern with a sophisticated
data analysis which allows one to extract the key model
ingredient — the interaction function
f共l兲— directly from
the experimental runs. In this way the validity of the
model can be tested, and additional features required for
the description can be identified. As far as we know, this
is the first quantitative comparison between theory and
experiment for a sand pattern. The basic ideas are general,
and we expect that the theoretical formalism combined
with the experimental analysis can be used for related sand
patterns (e.g., the one studied in [5]), or other strongly
nonlinear systems.
One-dimensional vortex ripples can be created in an an-
nular channel [6,10], ensuring the pattern to be subject to
well-defined periodic boundary conditions. Freely grown
ripples are created from a flat bed by a coarsening process,
which eventually saturates at an equilibrium state, where
the ripple length
˜
l is almost independent of the frequency
n of the driving and proportional to the amplitude of the
oscillation of the plate a [8,11]. Here we are mostly con-
cerned with the stability and evolution of the ripple patterns
themselves, and not with the instability of the flat bed,
which has been discussed elsewhere [9,12]. By creating
a homogeneous ripple pattern, where all ripples have the
same length, and changing amplitude and frequency, the
stability of the pattern is probed. In this way it was found
that there is a stable band of ripples, l
min
,
˜
l,l
max
, for
a given set of driving conditions [7].
Recently, a simple model was proposed, which describes
both the coarsening and saturation of ripples, and repro-
duces the existence of a stable band [9]. In its simplest
form, the change of the length l
j
of the ripple j is a func-
tion of l
j
itself and the lengths of the neighboring ripples:
ᠨ
l
j
苷 2f共l
j21
兲 1 2f共l
j
兲 2 f共l
j11
兲 . (1)
The interaction function f共l兲 describes the transfer of
mass, and consequently length, between neighboring
ripples. Behind each ripple, due to the water oscillation,
a separation vortex forms, which causes an exchange
of mass between the ripples. f共l
j
兲 is the amount of
sand which is taken by a ripple of length l
j
from its
downstream neighbor during one-half of the oscillation
(Fig. 1). At the same time ripple j loses an amount of
mass f共l
j21
兲 due to the action of the upstream ripple.
Adding the reverse process for the second half of the
oscillation results in the model (1) [9]. The purpose of
this Letter is a test of this model through comparison with
experimental data.
Experimental setup.—As shown in Fig. 2, an 11 mm
wide, 15 cm high annular Plexiglas channel of diameter
48.6 cm is filled with water, A, and glass beads, B, with
a diameter of 250 6 50 mm. In the middle is a conical
234302-1 0031-9007兾02兾88(23)兾234302(4)$20.00 © 2002 The American Physical Society 234302-1
VOLUME 88, N
UMBER 23 PHYSICAL REVIEW LETTERS 10J
UNE 2002
−0.015
−0.01
−0.005
0
0 0.2 0.4 0.6 0.8 1 1.2 1.4
f(
λ
) [cm/period]
λ
/a
f
l
(λ)
f( )λ
f
c
(λ)
f
r
(λ)
FIG. 1. The interaction function for a 苷 6 cm and n 苷
0.6 Hz. The thick gray line is the smooth average of the three
different functions (mean values subtracted). The right limit
of the plot at l
兾a 苷 1.45 is the limit of stability, where new
ripples are created 共l
max
兲. The inset shows a sketch of the
ripples in the part of the oscillation when the flow is from the
left to the right. The interaction function can be interpreted
as the transport of sand in the trough, across the vertical
dashed line.
mirror, C, which is filmed from above by a CCD camera
with a resolution of 640 3 480 pixels. The channel is
attached to a motor, D, by a 1.5 m long arm E, which
oscillates the channel in an almost sinusoidal fashion.
We create initial conditions with small ripples by oscil-
lating at a small amplitude. The experiment starts when
amplitude and frequency are changed to their desired val-
ues. Throughout this paper, we use only the condition
a 苷 6 cm, n 苷 0.6 Hz. We want to emphasize that there
is nothing special about this condition, the same quali-
tative results are obtained for other parameters. Grow-
ing the ripples from initial conditions with homogeneous
ripples smaller than l
min
results in a final number of
19–21 ripples with a length (averaged over 16 realiza-
tions) of
˜
l兾a 苷 1.25. The uncertainty in the length due to
A
B
C
D
C
B
E
FIG. 2. A sketch of the experimental setup seen from above
(left) and from the side (right). The length of the arm E and the
width of the channel are not to scale. See text.
periodicity is around 5%. The images are transformed back
to Cartesian coordinates, and the lengths l
j
are found by
fitting them to triangles with fixed slope (Fig. 3a). From
series of consecutive l
j
, the temporal change,
ᠨ
l
j
, is es-
tablished, and a space-time plot of the ripple evolution is
constructed. Figure 3b shows the evolution starting from
a small initial wavelength, which leads to the annihilation
of ripples. Around each annihilation the ripples have not
been analyzed, as it becomes inaccurate to fit the triangles
during these events.
Measuring the stable band.—By running the experi-
ment with different values of a until a homogeneous
pattern occurs, we create initial conditions with a different
number N of ripples. From this initial pattern, the experi-
ment is started with the above parameters, and run for
10 000 periods. We have created initial conditions with
N
[ 关17; 24兴 and observed that for N . 22 ripples are
annihilated, while for N , 18 one or more new ripples
are created. We therefore conclude that there exists a stable
band N
[ 关17.5; 22.5兴 (for a 苷 6 cm, n 苷 0.6 Hz),
which corresponds to l
min
兾a 苷 1.13 6 0.03 and
l
max
兾a 苷 1.45 6 0.05. Note that this is a more accurate
measurement of the stable band than that conducted in
a
0
100
200
300
400
500
600
700
5 10 15 20
νt
x/
a
0
100
200
300
400
500
600
700
5 10 15 20
νt
x/
a
(a)
(b)
(c)
FIG. 3. (a) An example of an extracted profile (opaque re-
gion) and the fitted triangles with constant slope. The line is
shown above the profile for clarity. (b) Space-time plot of the
experimental evolution of the position of the ripple crests starting
from ripples with lengths 2.5 cm and evolving with a 苷 6 cm
and n 苷 0.6 Hz. (c) A simulation of the model (1) using the
extracted interaction function and the same initial conditions
as above.
234302-2 234302-2
VOLUME 88, N
UMBER 23 PHYSICAL REVIEW LETTERS 10J
UNE 2002
[6,7], as we are always forcing with the same values of
a and n.
Mass transfer model.—The model (1) was originally
presented in Ref. [9], together with a more refined ver-
sion. A stability analysis of the homogeneous state l
j
⬅
¯
l
shows that ripples are stable (unstable) if the derivative
f
0
共
¯
l兲 , 0 共.0兲 . For a convex interaction function, the
lower stability boundary l
min
therefore lies at l
marg
, where
f
0
共l
marg
兲 苷 0. To investigate nonlinear pattern evolution
within this model, Eq. (1) is supplemented by the rule that
ripples which reach zero length are annihilated and re-
moved from the system of equations, while the remaining
ones are relabeled.
Simulations of the model starting from typical initial
conditions [13] in the unstable band l,l
marg
show that
an equilibrium wavelength is reached which is essentially
independent of the initial wavelength, and depends only
on the interaction function f共l兲. To gain some insight
into the selection mechanism, it is useful to recast the
model into potential form by writing it in terms of the
position x
j
of the troughs between the ripples defined by
l
j
苷 x
j
11
2 x
j
. Then
ᠨ
x
j
苷 2≠V 兾≠x
j
, with the potential
V given by
V 苷 2
N
X
j苷1
Z x
j11
2x
j
0
dl f共l兲 , (2)
where N is the number of ripples. It is then plausible to
conjecture that the equilibrium length
˜
l can be found by
minimizing V for homogeneous ripples, under the con-
straint that the total length L
苷 N
˜
l is conserved. This
implies that
˜
l is determined through the Maxwell construc-
tion applied to f:
Z
˜
l
0
dl f共l兲 苷
˜
lf共
˜
l兲 . (3)
Comparison with numerical simulations shows that (3) sys-
tematically overestimates
˜
l, with a better performance the
steeper the stable branch of f. Our interpretation is that
the deterministic dynamics gets stuck in the multitude of
metastable states of (1); recall that
any homogeneous state
with
¯
l.l
marg
is a stable, stationary solution. Since the
mean ripple length can increase only by annihilations, its
evolution freezes once all ripples are in the stable band.
The wavelength predicted by (3) should therefore be an
upper bound on the actual equilibrium wavelength, which
is true in all cases we have considered.
Data analysis.— Given the time series l
j
共t兲 we want to
(i) evaluate how well the model (1) describes the evolution
of the ripples and (ii) extract the interaction function f共l兲.
For the analysis, we write (1) in the more general form
ᠨ
l
j
苷 2f
l
共l
j21
兲 1 2f
c
共l
j
兲 2 f
r
共l
j11
兲 , (4)
where f
l
, f
c
, and f
r
denote the left neighbor, center,
and right neighbor interaction function, respectively. This
yields an additional degree of freedom as the functions are
not required to be equal.
We want to determine the optimal transformations
f
l
共l
j21
兲, f
c
共l
j
兲, f
r
共l
j11
兲 in the sense that they minimize
the error x
2
苷
P
t
关
ᠨ
l
j
共t兲 1 f
l
共t兲 2 2f
c
共t兲 1 f
r
共t兲兴
2
.
The problem is solved numerically by the alternating
conditional expectation value algorithm [14,15]. The
algorithm works by varying f
l,c,r
in the space of all
measurable functions until convergence to the absolute
minimum of x
2
[14]. The results are nonparametric
functions which are given in numerical form by points,
e.g., 关l
j
, f
c
共l
j
兲兴. An upper bound on the error, e.g., on
f
c
, is given by s
f
c
苷 max
l
共df
c
兾dl兲 ?s
l
, with s
l
being
the measurement error at the point l; analogous estimates
apply to f
l
, f
r
. The error in the points on the l axis is
estimated using the above result and equals roughly the
errors in the given values of l.
The quality of each of the resulting functions is given by
the maximum correlation c
l,c,r
of one of the terms with the
sum of all the others [16]. A value of c 苷 1 implies a per-
fect result, lower values indicate either imperfect modeling
or (measurement) noise or both.
Results.—We have performed the analysis on data from
18 realizations of a coarsening process initiated with ap-
proximately 60 small ripples and run with a 苷 6 cm and
n 苷 0.6 Hz. Figure 1 shows the three functions f
l
, f
c
,
and f
r
averaged over the 18 runs. Clearly there is some
noise, but the functions are very similar as expected from
(1). The maximum correlations are found to be C
l
苷 0.66,
C
c
苷 0.88, C
r
苷 0.72. The uncertainty in estimating l
amounts to about 2 mm, with a maximum slope of f of
around 3 3 10
23
共period兲
21
, the absolute error of f is
6 3 10
24
cm兾period. This means that the result does not
fluctuate due to lack of data, rather the model (4) does not
account for all the variation in the ripple evolution. We
will return to this point later.
The maximum of the function lies around l
marg
苷
1.08a, which is a little larger than what was found by
the numerical simulations in [9]. Concerning the shape
of the function, it is interesting that the slope in the
unstable band l,l
marg
is much larger than the slope
in the stable band l.l
marg
. As the slope is a measure
of the time scale of the dynamics, this implies that the
initial coarsening stage is much faster than the equili-
brating stage.
To use the resulting function to integrate (1) numeri-
cally for evaluating (3), information about the smallest
ripples is needed. With a lack of such information,
we have extrapolated linearly to l 苷 0, Fig. 1, dashed
line. Varying the slope of the extrapolated part by
a factor of 2 in either direction does not change the
final number of ripples. Figure 3c shows a space-time
plot of the numerical integration. The model predicts
a qualitatively similar behavior as the experiment;
however, the instability of the small ripples develops
slower in the model than in the experiment. If the
extrapolated slope for the interaction function is made
steeper, this instability will evolve faster in the model.
234302-3 234302-3
VOLUME 88, N
UMBER 23 PHYSICAL REVIEW LETTERS 10J
UNE 2002
The model also overestimates the final ripple length.
For similar initial conditions and system size, the final
number of ripples is typically 18, whereas the experiment
yields 20.
Evaluating the Maxwell construction (3) using the mea-
sured function and extrapolating to account for the larger
ripples, gives
˜
l 苷 2.2a, which lies far outside of the range
of definition of f共l兲. Evidently, the experimentally deter-
mined f共l兲 belongs to the class of interaction functions
for which the upper bound given by Eq. (3) is not useful.
The upper bound of the stable band, where new ripples
are created, was found to be l
max
苷 1.45a. At this point,
an infinitesimally small ripple, inserted in the trough be-
tween two larger ripples, is able to gain mass from the
neighbors, and thus grows. If the model shall capture this
behavior, Eq. (1) yields that f共l
max
兲 should be smaller
than f共0兲 [9]; otherwise the right-hand side is negative
and a small ripple disappears. The measured interaction
function (cf. Fig. 1), however, suggests the opposite. Even
with the large uncertainty inherent in the extrapolation,
a bend which makes f共1.45a兲 , f共0兲 is hard to imag-
ine. We therefore conclude that the model in the form
(1) is not able to quantitatively predict the creation of new
ripples. The reason for this apparent failure relates to the
model being developed essentially as an expansion around
a homogeneous state [9]. But in the case of a creation, the
state is as inhomogeneous as it can be: a very tiny ripple
flanked by large ripples. To account for this, the interaction
function should be described as a function of the length of
the ripple creating the separation bubble, l
i
, but also of the
length of its neighbors. Ongoing numerical work indicates
that it is most relevant to write f ⬅ f共l
i
, l
i11
兲 where l
i
is the length of the ripples creating the separation bubble,
and l
i11
is the one “touched” by the separation bubble (see
inset of Fig. 1). In fact, the whole coarsening process is
dominated by highly inhomogeneous configurations, and
this might also be why simulating (1), using the measured
interaction function, does not produce exactly the correct
final ripple length.
Conclusions.—We have demonstrated that we can ex-
tract the interaction function f共l兲 from spatiotemporal data
of the evolution of the profile of the sand surface. The
nonlinear data analysis shows its full strength, producing
nonparametric function estimates, where any parametric
approach would have failed. With the interaction function
it is possible to model the evolution of the ripples, in par-
ticular, their coarsening and equilibration.
The prediction from the model of the existence of a
stable band is indeed observed in the experiment. The
lower bound is predicted at l
marg
兾a 苷 1.08 6 0.03, while
the measurement gives l
min
兾a 苷 1.13 6 0.03. In a simi-
lar experiment [6], the lower bound of the stable band was
found to coincide with the final ripple length: l
min
苷
˜
l,
whereas we find that l
min
,
˜
l. The difference between
the two results might be apparent only, as the number of
ripples in the annulus of [6] was 10 or smaller, giving rise
to an uncertainty in the ripple length on the order of the
difference between l
min
and
˜
l.
The dynamics responsible for the evolution of a 1D
ripple pattern are also relevant for 2D ripple patterns. How-
ever, in the latter case, topological defects will be present
[7], and the final length selection might be determined by
the motion of these [9].
It is a pleasure to acknowledge discussions with
T. Bohr, M. van Hecke, and F. Schmidt. J. K. is grateful to
DTU and NBI for gracious hospitality, and to DFG within
SFB237 for support. K. H. A thanks Universität Essen for
hospitality.
[1] M. C. Cross and P. C. Hohenberg, Rev. Mod. Phys. 65,
851
(1993), and references therein.
[2] H. Nishimori and N. Ouchi, Phys. Rev. Lett. 71
, 197
(1993).
[3] B. T. Werner and D. T. Gillespie, Phys. Rev. Lett. 71
,
3230– 3233 (1993).
[4] Z. Csahók, C. Misbah, F. Rioual, and A. Valance, Eur.
Phys. J. E 3
, 71 (2000).
[5] A. Betat, V. Frette, and I. Rehberg, Phys. Rev. Lett. 83
,
88– 91 (1999).
[6] A. Stegner and J. E. Wesfreid, Phys. Rev. E 60
, R3487
(1999).
[7] J. L. Hansen, M. v. Hecke, A. Haaning, C. Ellegaard, K. H.
Andersen, T. Bohr, and T. Sams, Nature (London) 410
, 324
(2001); J. L. Hansen, M. v. Hecke, C. Ellegaard, K. H. An-
dersen, T. Bohr, and T. Sams, Phys. Rev. Lett. 87
, 204301
(2001).
[8] R. A. Bagnold, Proc. R. Soc. London A 187
, 1–15 (1946).
[9] K. H. Andersen, M.-L. Chabanol, and M. v. Hecke, Phys.
Rev. E 63
, 066308 (2001).
[10] M. A. Scherer, F. Melo, and M. Marder, Phys. Fluids 11
,
58 (1999).
[11] P. Nielsen, J. Geophys. Res. 86
, 6467 (1981).
[12] P. Blondeaux, J. Fluid Mech. 218
, 1–17 (1990); K. H. An-
dersen, Phys. Fluids 13
, 58 – 64 (2001).
[13] We consider here initial conditions with some amount of
disorder. Different behavior results in the (unrealistic)
case of a perfectly ordered state with a localized pertur-
bation; see J. Krug, Advances in Complex Systems (to be
published).
[14] L. Breiman and J. H. Friedman, J. Am. Stat. Assoc. 80
,
580 –598 (1985).
[15] H. Voss, P. Kolodner, M. Abel, and J. Kurths, Phys. Rev.
Lett. 83
, 3422– 3425 (1999).
[16] For example, for
f
l
: c
l
苷 具f
l
共t兲 ? S共t兲典兾共具 f
2
l
典具S
2
典兲;
S共t兲 苷
ᠨ
l
i
2 2f
c
1 f
r
has zero mean (without loss of
generality).
234302-4 234302-4