» »

4μA 1 2π 2π x 1 1 2π b

(ρ1 ’ ρ2 )gw = ’ + ’

cos tanh

» » »

sinh 2πb cosh 2πb sinh 2πb cosh 2π b

b 2π b 2π b

» » » »

(7.10)

is evaluated at the unperturbed position of the Solving Eq. (7.10) for A 1 and substituting the re-

interface (i.e. y = 0). The value of v 1 at y = 0 is sult into Eq. (7.6) gives

obtained by differentiating Eq. (7.4) with respect

‚w (ρ1 ’ ρ2 )gb

to x and evaluating the resulting expression at =w

‚t 4μ

y = 0. The result is

2π b ’1

»2

tanh 2π b ’ sinh 2π b cosh

‚w » » »

2π A 1 2π x — 2πb

= 2πb ’1

cos (7.6) »

+ sinh 2π b cosh

‚t » » » »

2πb

(7.11)

In order to constrain the value of the coef-

¬cient A 1 in Eq. (7.6), we must incorporate the Equation (7.11) describes an exponential growth

role of buoyancy on ¬‚uid motion. As the interface process with time scale „ given by

between the two ¬‚uids is disturbed downward

2π b ’1

»

(or upward) the less (more) dense ¬‚uid will dis- + sinh 2πb

cosh

4μ » »

„= 2πb

place the other. The pressure difference between ’1

(ρ1 ’ρ2 )gb 2

»

tanh 2πb ’ sinh 2πb

cosh 2π b

» » »

2πb

the two sides of the interface resulting from this

(7.12)

buoyancy force is given by

Equation (7.12) provides the necessary relation-

2(P 1 ) y=0 = ’(ρ1 ’ ρ2 )gw (7.7)

ship between the time scale of the instability and

The pressure on y = 0 in the upper layer can the wavelength ».

be found by using a force balance equation that Figure 7.3 plots the relationship between

equates the pressure gradient in the ¬‚ow with the dimensionless growth rate of a perturba-

the viscous stresses: tion (i.e. 1/„ ) and its dimensionless wavelength

‚P ‚ 2u ‚ 2u 2πb/». The fastest growing wavelength can be

’μ +2 =0 (7.8)

determined by differentiating Eq. (7.12) as a

‚x ‚ x2 ‚y

164 INSTABILITIES

ratio of meander wavelength and channel width

100

is a constant equal to about ten for these diverse

lmax= 2.568b

phenomena (Leopold et al., 1964). This similarity

suggests a common mechanism for meandering.

10

(r1 ’ r2)gbt

In this section we present a linear stability

analysis which predicts channel migration rates

4

and the proportionality between meander wave-

10 length and channel width using a very simple

geometric approach that applies to both river

channels and lava channels. Most models of allu-

vial channel meandering stress the importance of

10

102

100 101

10 the cross-sectional circulation or secondary ¬‚ow

2pb/l

in meander initiation (Callander, 1978; Rhoads

Fig 7.3 Plot of the inverse of the dimensionless time scale „ and Welford, 1991). This model does not replace

as a function of the dimensionless wavelength for the

those more sophisticated models, but rather es-

Rayleigh“Taylor instability (Eq. (7.12)).

tablishes suf¬cient conditions for meandering to

occur in a wide variety of channel ¬‚ows. Mean-

function of „ and setting the result equal to zero. ders of the correct scaling form when the bank

shear stress (in the case of sediment transport)

Turcotte and Schubert (2002) give the result of

or temperature (in the case of thermal erosion)

that equation as

is proportional to the curvature of the bank and

»max = 2.568b (7.13)

cross-sectional variations in shear stress or tem-

perature have a linear form. For alluvial chan-

where »max is the fastest-growing wavelength. Per-

nels, this curvature dependence arises from the

turbations of this wavelength have a growth time

centripetal force of ¬‚uid rounding a bend, and

scale given by

in thermal erosion it arises through the excess

13.04μ

„= latent heat of a curved bank. A linear stability

(7.14)

(ρ1 ’ ρ2 )gb analysis predicts that channels are unstable to

growth at all wavelengths larger than a critical

Equation (7.12) predicts that the fastest growing

wavelength given by » = 2π w, where w is chan-

wavelength is equal to about two and half times

nel width, and stable below it (provided that

the layer thickness. The Rayleigh--Taylor instabil-

some channel ¬‚ows exceed the critical condition

ity provides a classic example of a linear stability

required for erosion). The wavelength of fastest

analysis. We will use the linear stability analysis

growth predicted by this analysis is equal to

procedure again in Section 7.5, to explore the in- √

2 3πw or 10.88w, which is very similar to ob-

stability mechanism that creates oscillating arid

servations from alluvial and lava channels.

alluvial channels.

In alluvial rivers, the dependence of channel

migration on the local curvature of the bed has

been stressed in measurements of channel mi-

7.3 A simple model for river

gration by Nanson and Hickin (1983) and in the-

meandering oretical studies by Begin (1981) and Howard and

Knutson (1984). Nanson and Hickin (1983) found

Alluvial rivers, supraglacial meltwater streams, that channel migration rates were a maximum

within a narrow interval of R /w near 3, where

Gulf Stream meanders, and lava channels on

Earth, the Moon, and Venus (Komatsu and Baker, R is the radius of curvature of the meander,

1994) all exhibit meandering with similar propor- and w is the channel width. We will compare

tionality between meander wavelength and chan- our theoretical results to their data. Begin (1981)

nel width (Leopold et al., 1964; Stommel, 1965; showed, by computing the centripetal force nec-

Parker, 1975; Komatsu and Bakerm, 1994). The essary to accelerate ¬‚ow around a bend, that the

7.3 A SIMPLE MODEL FOR RIVER MEANDERING 165

v

t Fig 7.4 Schematic diagram of the

w2 meander instability model

y calculation.

w

0

l

shear stress at the bed is proportional to the cur- Eqs. (7.15) and (7.17) for the cross-sectional and

vature. Howard and Knutson (1984) developed a along-channel variation in shear stress gives

simulation model of meandering where the mi- 2

„ (x, y) = „b,s (1 + w K ) 1 ’ y + „c (7.18)

gration rates were proportional to local curvature w

up to some maximum value. They obtained real-

The shape of the curved-channel bank within

istic meandering channels with their model.

this linear-stability approach is given by yb =

In Figure 7.4 we present the geometry to be

sin(2π x/»). The curvature K is given by

considered in our model. We assume symmetry

d yb /dx 2 . Substituting the sinusoidal bank shape

2

within our 2D model, i.e. that the forces and ero-

into Eq. (7.18) gives

sion on one bank are equal in magnitude and

2

opposite in direction to those on the other bank. 2 2π

„ (x, y) = „b,s 1’ y 1+w

»

First we consider the distribution of boundary w

shear stress in a straight channel of width w. The 2π x

— sin + „c (7.19)

fundamental assumption of this analysis is that »

the cross-sectional shear stress pro¬le is linear,

We further assume that spatial variations in bed

decreasing from a maximum value at the bank

shear stress produce erosion and accretion of the

(y = w/2) to a minimum value at the channel

bank proportional to the variation in shear stress.

centerline:

In other words, the bank will erode in locations

2

„ (y) = „b,s 1 ’ y + „c where the variation in shear stress downriver is

(7.15)

w positive, and a point bar will prograde in places

where „b,s is the maximum shear stress at the where variations in shear stress are negative. The

bank for a straight channel and „c is the mini- rate of migration, therefore, is proportional to

the derivative of „ (x, y) with respect to x, evalu-

mum shear stress at the centerline. In a curving

ated at y = w/2:

channel, empirical data indicate that the bank

shear stress is equal to (Richardson, 2002) ‚„ d 2π x

= sin

‚x »

dt

w 2

y=w/2

„b,c = „b,s +1 (7.16)

2R 2

w 2π 2 2π x

= „b,s ’ sin

where „b,c refers to the bank shear stress for a » »

2 w

curved channel, and R is the radius of curvature (7.20)

of the channel centerline. For incipient mean-

The normalized growth rate of perturbations of

ders, the radius of curvature of the centerline

the bed is given by

is much larger than the channel width. As such,

we can use the approximation (1 + )2 ≈ 1 + 2 2

™ 2π„b,s 2 w 2π

= ’ (7.21)

to rewrite Eq. (7.16) as » »

w 2

w

„b,c ≈ „b,s 1 + = „b,s (1 + w K ) (7.17) which can be simpli¬ed to

R

2

™ 1 2π w

where K is the planform curvature of the cen-

∝ 1’ (7.22)

» »

terline (equal to the inverse of R ). Combining

166 INSTABILITIES

is dif¬cult, however, to compare migration rates

0.8

from many different bends without ¬rst quantify-

lmax = 10.88w

ing the effects of bank-material texture. Lumping

migration rate (m/yr)

0.6

data from meander bends with different erodi-

bilities would increase the scatter compared to

the Nanson and Hickin data, even if individual

0.4

bends closely follow a universal growth curve. Mi-

gration rate data from within individual bends

would help to resolve this question.

0.2

The analysis that led to Eq. (7.22) also ap-

plies to channels carved by thermal erosion.

0.0 Supraglacial meltwater streams (Parker, 1975)

40

30

10 20

0

l/w and lava channels on the Earth, the Moon, and

Venus (Hulme, 1973) form by thermal erosion. In

Fig 7.5 Plot of the normalized growth rate of perturbations the case of lava channels, melting of the channel

as a function of wavelength » (Eq. (7.22)). Also shown (dots) bed by hot lava (and solidi¬cation on the opposite

are the bank-migration data of Nanson and Hickin (1983)

bank) causes channel migration to occur. Curved

along the Beatton River for comparison.

banks have a melting point that depends on the

curvature of the bank. The melting temperature

of the bank is given by

Equation (7.22) is plotted in Figure 7.5 as the solid

curve. The maximum growth rate occurs where

T = Tm + c K (7.23)

the ¬rst derivative of ™/√ with respect to » is

zero. This occurs at » = 2 3πw = 10.88w, which

where T m is the melting point for a straight chan-

is close to the observed scaling relation (Leopold

nel, K is the planform curvature of the bank, and

et al., 1964).

c is a constant that depends on the latent heat of

For comparison, Figure 7.5 shows observed mi-

the bank material. The temperature required to

gration rate data for the Beatton River from Nan-

melt a curved bank is elevated because eroding

son and Hickin (1983). These authors employed

the bank by a given amount requires melting a

dendrochronological techniques to estimate the

larger volume of material per unit surface area

migration rates for sixteen bends of that river.

than for a straight channel. As such, the tempera-

They expressed their data as a function of the

ture boundary condition in lava channels has the

average radius of curvature of the bend. To com-

same curvature dependence as the shear stress

pare their results to our theoretical growth curve,

boundary condition in alluvial channels. In ad-

we used the relation » = π R which is applica-

dition, turbulent eddies will transport heat away

ble for an ideal shape of fully developed me-

from the hot center of the ¬‚ow towards the cool

anders: the sine-generated curve (Langbein and

sides of the ¬‚ow in the same way that turbu-

Leopold, 1966). Their results show the same de-

lent eddies transport shear stress through an al-

pendence with meander wavelength predicted

luvial channel, resulting in an approximately lin-

by our linear stability analysis. The migration

ear temperature pro¬le. These underlying simi-

rate rises rapidly for small »/w, reaches a maxi-

larities between alluvial and lava channels sug-

mum value, and then decreases. Migration rates

gests that the same geometrical instability mech-

decrease rapidly for tightly curved bends with

anism is at work in both cases.

»/w < 2π , consistent with our model™s prediction

that channel bends become stable below that

wavelength. It should be noted that some stud- 7.4 Werner™s model for eolian

ies following Nanson and Hickin (1983) showed

dunes

that migration rates are far less systematic that

the results obtained for the Beatton River. Most

of those studies have concluded that the Nanson Eolian ripples, dunes, and megadunes obey a

and Hickin results are not generally applicable. It striking periodicity. Many attempts have been

7.4 WERNER™S MODEL FOR EOLIAN DUNES 167

made to understanding the controls on eolian Folk, 1976). Wilson (1972) has proposed the same

bedform spacing in terms of both microscale pro- mechanism for ripples and megadunes with the

cesses (e.g. the trajectories of individual grains bedforms on these scales created by different

in saltation and reptation) and macroscale pro- scales of atmospheric circulation. However, care-

cesses (e.g. the interaction between a growing ful ¬eld studies have not found evidence for the

bedform, the wind ¬‚ow above it, and the result- existence of these persistent atmospheric circula-

ing pattern of erosion and deposition). Despite tions (Livingstone, 1986).

the efforts of many scientists, the mechanisms In 1995, Werner made a signi¬cant advance in

responsible for the formation of eolian ripples, eolian geomorphology in a paper that described

dunes, and megadunes and the factors responsi- a very simple but powerful model of eolian dune

ble for their size and spacing are not fully un- formation (Werner, 1995). Werner™s model begins

derstood. Perhaps the most well studied type of with a set of discrete sand slabs on a rectangular

bedform is the eolian ripple. Bagnold (1941) per- grid. Slabs may represent individual sand grains

formed the pioneering work on the formation of or a collection of grains; the number of slabs on

eolian ripples. He proposed that ripples form by the grid is limited for computational reasons, but

a geometrical instability in which a perturbation the model is not sensitive to the size of each slab.

exposes the windward side of the perturbation During each time step of the model, a sand slab is

to more impacts and faster surface creep than picked up from a randomly chosen pixel. If there

the leeward side, resulting in more grains being is no sand at that pixel, nothing happens. If a

ejected into saltation on the windward slope. In slab is present, it is picked up and moved down-

his theory, the spacing of ripples is related to a wind a distance of l grid point and deposited

characteristic saltation distance that is constant back on the surface with a probability that de-

in time. Sharp (1963), however, observed that in- pends on the presence or absence of sand at the

cipient ripples increase in spacing over time until new location. If sand is present at the new loca-

a steady-state spacing is achieved. In addition, ex- tion, the slab is deposited with probability ps . If

perimental studies of grain impacts suggest that no sand is present, it is deposited with a lower

a wide distribution of energies are imparted to probability pns . In nature, sand is more likely to

grains on the bed during each impact, resulting be deposited on a patch of sand than on a bare

in a wide distribution of saltation path lengths desert surface because the boundary layer above

(Mitha et al., 1986). As such, Bagnold™s hypothesis a sandy surface is usually rougher than above a

is not well supported by available data, and the bare, smooth surface, and because energy is ab-

controls on ripple spacing are still not well un- sorbed from the saltating grain by the granular

derstood. Anderson and Bunas (1992) developed bed. If the slab is not deposited, it is transported

a model that produced realistic ripples based on another distance l and again deposited with a

probability given by ps or pns . One additional rule

the microscale processes of grain movement and

controls the probability of deposition: if a sand

bed impact, but the spacing of ripples in their

slab lands in a ˜˜shadow zone,™™ it is deposited with

model was controlled by an ad hoc ˜˜ceiling™™ in

100% probability. A shadow zone is de¬ned as the

the model domain that does not occur in nature.

domain downwind of any topographic high that

The mechanisms controlling the formation

lies below a plane with an angle of 15—¦ to the

and spacing of eolian dunes are also incom-

horizontal. To map the shadow zone, a path is

pletely understood. Eolian dunes represent a dis-

traced from every local ridge down a 15—¦ slope

tinct bedform type from eolian ripples since

in the downwind direction until the plane in-

they occur at a larger scale with no transitional

tersects the surface. The 15—¦ angle of the shadow

bedforms. It has been proposed that dunes re-

zone can be varied, but its value should be signif-

sult from the presence of pre-existing wave-like

icantly less than the angle of repose. Physically,

motions or secondary circulations in the atmo-

the shadow zone represents the zone of recircu-

spheric boundary layer. Periodic motions intrin-

lation on the lee side of a growing dune. In this

sic to the ¬‚ow or resulting from the response

zone, air¬‚ow velocities are greatly reduced and

of the ¬‚ow to a topographic disturbance upwind

any sand in saltation will likely be deposited.

may produce periodic bedforms (Wilson, 1972;

168 INSTABILITIES

wind Fig 7.6 A schematic diagram of

l Werner™s model for eolian dune

shadow

formation and evolution. After

zone

Werner (1995).