angle of repose

Finally, when the sand slab is deposited it may the example of Figure 7.7) and the direction of

create an oversteepened condition in which the transport (e.g. star dunes form if the sand bed

local slope is higher than the angle of repose. If is initially thin and the wind direction is pe-

this occurs, the model moves the slab down the riodically reversed). Werner™s model provides a

direction of steepest descent until a site is found nice illustration of the power of simpli¬ed mod-

that does not lead to an oversteepened condition. els for understanding the behavior of complex

Figure 7.6 illustrates the model schematically. A systems with positive feedback. Werner™s model

code that implements Werner™s model is provided can easily be modi¬ed to account for variable

in Appendix 6. bed topography. Figure 7.8, for example, illus-

Initially, a rough surface develops in Werner™s trates a dune ¬eld migrating through a crater. On

model from an initially ¬‚at surface due to the Mars, the crestlines of dunes are deformed down-

random entrainment of grains from the surface wind of craters, similar to the effect shown in

(i.e. sand from some areas will be entrained more Figure 7.8.

than in others because of the stochastic nature of The height and spacing of dunes in Werner™s

entrainment in the model). Incipient ridges then model increases continuously over time with-

develop small shadow zones that encourage de- out ever reaching a steady-state value. Numeri-

position on their lee sides. This deposition leads cal models that incorporate the interaction be-

to the development of taller dunes with longer tween dunes and the air¬‚ow above them, in

shadow zones in a positive feedback. Using dif- contrast, predict that dunes eventually achieve

ferent parameter values, Werner (1995) was able an optimum dune height and spacing. In na-

to reproduce a variety of dune types by vary- ture, it is unclear which model is correct. In lo-

ing the initial sand thickness (a relatively thin cations where dunes and megadunes are both

initial sand bed leads to barchan dunes, as in present, these two bedform types appear to

wind

100 — L2 timesteps 200 — L2 timesteps 300 — L2 timesteps

Fig 7.7 Results of Werner™s model for a grid of size

L = 500 and model parameters l = 5, pns = 0.4, ps = 0.6, occur at distinctly different spatial scales with

angle of repose = 30—¦ , and angle of shadow zone = 15—¦ for no intermediate bedforms (Wilson, 1972). In

three model durations. The initial condition is uniform sand other words, even in sand seas mature enough

with a thickness of 10 sand parcels on each pixel.

to have formed megadunes, dunes appear to

7.5 OSCILLATIONS IN ARID ALLUVIAL CHANNELS 169

the American West. In their model, channels ag-

wind grade and widen until a threshold slope and/or

width is reached, initiating incision in the over-

steepened reach and aggradation in the down-

stream reach. In this section, we explore a math-

ematical model for Schumm and Hadley™s arid cy-

cle of erosion and compare the model predictions

to the spatial pattern of oscillations commonly

observed in alluvial channels of the southwest-

ern US.

Oscillating alluvial channels in southern Ari-

zona are classi¬ed as one of three basic types

depending on their geomorphic position. Pied-

mont channels fed by montane drainage basins

are continuous channels that alternate between

braided and narrowly entrenched reaches with

wavelengths on the order of 1 km. Wild Burro

Wash on the Tortolita Piedmont is a classic ex-

ample of this type of oscillating channel (Figure

Fig 7.8 Results of Werner™s model modi¬ed to include

variable bed topography. In this example, a dune ¬eld interacts 7.9a). Channel width (measured from 1 m reso-

with a crater. The bed topography causes deformation of lution US Geological Survey Digital Orthopho-

dune crest lines within and downwind of the crater.

toquadrangles (DOQQs)) and bed elevation (mea-

sured from a high-resolution DEM with average

channel slope subtracted) in Wild Burro Wash are

reach a steady-state spacing of approximately 10--

plotted in Figure 7.9a. These data illustrate that

200 m, with stronger winds and larger grain sizes

oscillations in both bed elevation and channel

associated with greater spacing. This observation

width occur in Wild Burro and that the two os-

suggests that the aerodynamic interaction be-

cillations are in phase.

tween dunes and the air¬‚ow above them is essen-

The second type of oscillating channel is

tial for understanding the height and spacing of

the piedmont discontinuous ephemeral stream.

mature eolian dunes.

These channels are fed only by local piedmont

runoff, and they alternate between incised chan-

nels and sheet-¬‚ow-dominated channel fans with-

7.5 Oscillations in arid alluvial

out well-de¬ned banks (Bull, 1997). The transi-

channels tions between channel fans and incised reaches

occur abruptly at one or more steep ˜˜headcuts™™

Most ¬‚uvial-geomorphic studies that have mod- (Figure 7.9b). Dead Mesquite Wash in southeast

eled the evolution of channel longitudinal pro- Tucson is the classic example given by Bull (1997)

¬les have assumed channels that are uniform (Figure 7.9b). Packard (1974) documented in-phase

in width or a prescribed function of drainage oscillations in bed elevation and width in Dead

area (e.g. Begin et al., 1981; Slingerland and Mesquite, similar to the behavior of Wild Burro.

Snow, 1987; Sinha and Parker, 1996; Dade and The third type of oscillating channel is the

Friend, 1998). Actual channels often vary substan- classic southwestern valley-¬‚oor arroyo. These

tially in width, however, and these variations in systems alternate between narrow, deeply incised

width can strongly in¬‚uence channel evolution. channels and broad, shallow channel fans, but

Schumm and Hadley (1957), for example, devel- they typically drain one or more broad basins

oped their conceptual ˜˜arid cycle of erosion™™ on and adjacent ranges and have wavelengths on the

the basis of the episodic cut-and-¬ll histories and order of 10 km. Abrupt headcuts also character-

oscillating geometries observed in the arroyos of ize distributary-to-tributary transitions in these

170 INSTABILITIES

(a) Wild Burro DEM

(c) Vamori

downstream

300 m

B

downstream

photo

I

N B I

B I B I

>100 cm

1988 flood depths

50“10

30“50

10“30

0“10

w/h long./width

h

120/ w profile B

0.3

60/

0

I

0/

1500 x (m)

2500 0

1000 500

2000

(b) Dead Mesquite headcuts headcuts

incised

B

channel fans

300 m

downstream

I

N

5 km B

N

B I B I B I

Fig 7.9 Examples of oscillating channels in southern

two ways. First, a numerical model for the

Arizona. Alternating reaches: B“braided and I“incised.

coupled evolution of the channel longitudinal

(a) Wild Burro Wash data, including (top to bottom) a

high-resolution DEM (digital elevation model) in shaded relief pro¬le and cross section is presented. This model

(PAG, 2000), a color Digital Orthophotoquadrangle (DOQQ), predicts that arid alluvial channels are un-

a 1:200-scale ¬‚ow map corresponding to an extreme ¬‚ood on

stable over a range of wavelengths controlled

July 27, 1988 (House et al., 1991), and a plot of channel width

by channel width, depth, and slope. Second,

w (thick line) and bed elevation h (thin line, extracted from

a database of channel geometries in southern

DEM and with average slope removed). (b) Dead Mesquite

Arizona is constructed that con¬rms the model

Wash shown in a color DOQQ and a detailed map of channel

predictions and places the three channel types

planform geometry (after Packard, 1974). (c) Vamori Wash

within a single continuum from steep, small-

shown in an oblique perspective of a false-color Landsat

image (vegetation [band 4] in red) draped over a DEM. wavelength-oscillating channels to gently slop-

Channel locations in Figure 7.11a. Modi¬ed from Pelletier ing, long-wavelength-oscillating channels.

and Delong (2005). For color version, see plate section. The equation describing conservation of mass

for sediment states that erosion and deposition

systems. In contrast to the other oscillating- in an alluvial-channel bed are proportional to the

channel types, not all arroyos are periodically gradient of total sediment discharge, or

alternating; some are tectonically controlled by

‚h 1 ‚(wqs )

=’

the elevations of basin divides. In southern Ari- (7.24)

‚t C 0 ‚x

zona, Cooke and Reeves (1976) identi¬ed Vamori

Wash as perhaps the best example of an oscil- where h is the elevation of the channel bed, t

lating arroyo system (Figure 7.9c). Entrenched is time, C 0 is the volumetric concentration of

reaches are associated with low vegetation den- bed sediment, w is the channel width, qs is the

sity, whereas channel fans support a broad, dense speci¬c sediment discharge (i.e. the sediment dis-

riparian zone in Vamori Wash. charge per unit channel width), and x is the dis-

The common oscillatory pattern in these tance downstream. In most bedload transport re-

channel types suggests a common underlying lations, speci¬c sediment ¬‚ux is proportional to

mechanism that acts across multiple spatial the 3/2-power of the bed shear stress. This im-

scales. This section explores this hypothesis in plies that speci¬c sediment discharge is a linear

7.5 OSCILLATIONS IN ARID ALLUVIAL CHANNELS 171

function of channel gradient and a nonlinear are assumed to widen as they aggrade and nar-

function of the discharge per unit channel width, row as they incise. Bull™s model states that chan-

or nels incise if the stream power is greater than a

threshold value and aggrade if the stream power

b

‚h

Q

qs = ’B (7.25) is less than that value. Expressing this relation-

‚x

w

ship in terms of channel widening (aggradation)

where B is a mobility parameter related to grain and narrowing (incision) gives

size and Q is water discharge. The value of b is

‚w wqs ’ w 0 qs0

constrained by sediment rating curves and is be- =’ (7.29)

‚t h0 ’ h

tween 2 and 3 for both suspended-load and bed-

load transport. Here we consider the case b = 2. where qs0 is the equilibrium sediment ¬‚ux

Equation (7.25) assumes that the bed shear stress or stream power, and h 0 is the equilibrium

is much larger than the threshold for particle bank height. Equation (7.29) represents a cross-

entrainment. sectional mass balance: sediment removed from

If the channel width is assumed to be uniform the bank contributes to the local sediment-¬‚ux

along the longitudinal pro¬le, the combination de¬cit wqs ’ w 0 qs0 (i.e. the amount of sediment

of Eqs. (7.24) and (7.25) gives the classic diffusion that cannot be transported out of the reach), pro-

equation: moting further aggradation in the reach.

‚h ‚ 2h The nonlinear term in Eq. (7.28) alters the dy-

=κ 2 (7.26)

namics of channels markedly compared to the

‚t ‚x

diffusive behavior expressed in Eq. (7.26). Along a

where the diffusivity is given by κ = (B Q 2 )/

reach with uniform discharge and grain size, the

(C 0 w 0 ), and w 0 is the uniform channel width.

diffusion equation smoothes out curvatures in

If the channel width is not uniform along

the pro¬le over time. The nonlinear term in Eq.

the longitudinal pro¬le, the chain rule must be

(7.28), however, has the opposite effect through

used when differentiating qs in Eq. (7.24). This ap-

a positive feedback between channel width and

proach introduces an additional nonlinear term

slope. In this feedback, spatial variations in chan-

into Eq. (7.26), to give

nel slope generate variations in width via Eq.

‚h 1 ‚ 2h 1 ‚w ‚h (7.29). Large gradients in channel width, in turn,

= κw 0 ’2 (7.27)

‚t w ‚x w ‚x ‚x

2

increase the nonlinear behavior in Eq. (7.28), fur-

ther localizing erosion and deposition to com-

It is mathematically easier to de¬ne h as the

plete the feedback cycle. This cycle is balanced

difference between the local bed elevation and

by the diffusive term, and the balance between

that of a straight, equilibrium channel with uni-

these two terms controls the oscillation wave-

form discharge. Equation (7.27) then becomes

length.

‚h 1 ‚ 2h 1 ‚w ‚h

= κw 0 ’2 S0 ’ Equations (7.28) and (7.29) may be solved by us-

(7.28)

‚t w ‚x w ‚x ‚x

2

ing linear stability analysis and direct numerical

where S0 is the equilibrium channel slope. Equa- solution. Linear stability analysis works by solv-

tion (7.28) can be used to study the evolution of ing the linear approximation to Eqs. (7.28) and

perturbations from the equilibrium geometry. (7.29) for the growth rate of a small-amplitude

Channel widening and narrowing occurs by a oscillation superimposed on the initial channel

complex set of processes including bank retreat geometry (a channel with speci¬ed equilibrium

and bed scouring. Scouring often leads to chan- width w 0 , bank height h 0 , and slope S0 ).

nel narrowing as ¬‚ow is focused into the scour Linear-stability analysis involves solving for

zone and parts of the former channel bed effec- the growth rate of spatially periodic, small-

tively become part of the bank. Bull™s ˜˜threshold- amplitude perturbations from an equilibrium ge-

of-critical-power™™ concept (Bull, 1979) can be used ometry by retaining only linear terms in the

to relate the rate of channel widening and nar- model equations. In this linear stability analysis,

rowing with the excess stream power if channels we assume that the channel has an average width

172 INSTABILITIES

1

(b)

(a) 1/2

l = p(w0 h0/S0)

h

w 0.8

0.6

d 0.4

x

0.2

unstable stable

0

5

0

(c) 3 4 1/2

1 2 6

l ((w0 h0/S0) )

1

incision

0

h backfilling

random

l

h0

instability growth phase propagation phase

downstream downstream

3

w 2

w0

1

0

2 x (km) 3 4 5 2 x (km) 3 4 5

0 1 1

retaining only linear terms gives

Fig 7.10 Model behavior. (a) Schematic diagram of model

geometry. (b) Growth curve for linear-stability analysis,

) = 2π κ 2π ) + S0 w ei ( 2π x )

2π x 2π x

δh 0 ei ( h 0 ei (

’

where δ is the nondimensional growth rate. This analysis » » »

» » w0 0

predicts unstable behavior for small wavelengths, with a

(7.31)

maximum value at » = π (w 0 h 0 /S0 )1/2 . (c) Numerical

solution. Initial phase (left) is characterized by growth of

and

small, random perturbations and an increase in oscillation

)=’ κ

wavelength with time (line thickness increases with time). 2π

2π x +φ+ π 2π x + π

δw 0 ei ( h 0 ei ( ) (7.32)

» »

2 2

Steady-state phase (right) is characterized by solitary-wave »

w0h0

propagation of oscillations in upstream direction. Incision and

Equation (7.31) implies φ = ’π/2 (i.e. h and w

channel narrowing on leading edge of each wave are balanced

are 90—¦ out of phase). Substituting Eq. (7.32) into

by back¬lling and widening on trailing edge. Modi¬ed from

Pelletier and Delong (2005). Eq. (7.31) gives the following equation for growth

rate:

w 0 and slope S0 . The variables h and w refer to 2

κ κ S0 1

2π

δ= 1’ (7.33)

the variations of h and w from their average val- » h0w0 δ

2

w0

ues, and have solutions of the form

The solution to Eq. (7.33) is

2π x

) eδt , w = w ei ( 2π x +φ ) eδt

h = h 0 ei ( (7.30)

» »

0

2π κ π π

where δ is the growth rate and » is the wave- S0

2

δ= ’ ’ (7.34)

w0» » » w0h0

length. These solutions assume a spatially peri-

odic form, and an exponential increase in am-

plitude with time as the instability grows. Sub- Equation (7.34) is plotted in Figure 7.10b.

stituting Eqs. (7.30) into Eqs. (7.28) and (7.29) and This ¬gure shows that the growth rate of

7.5 OSCILLATIONS IN ARID ALLUVIAL CHANNELS 173

perturbations has a maximum at the ampli¬cation of spatial variations in width

within a range of wavelengths close to 1 km. In

w0h0 the latter stage, oscillations achieve a steady-state

»=π (7.35)

S0 amplitude and propagate upstream as a train of

solitary (i.e., nondispersive) topographic waves. At

This result can also be obtained by differentiating

any given instant, the channel geometry is char-

Eq. (7.34) with respect to » and setting the result

acterized by alternating zones of narrow, deeply

equal to zero.

incised reaches and wide, shallow reaches. The

This analysis predicts the growth curve shown

most abrupt slope break is where the channel

in Figure 7.10b. The growth rate is positive

changes from distributary to incised. This break

(i.e., perturbations are unstable) for small wave-

is similar to the headcuts often observed in ar-

lengths, rises to a steep maximum at

royos and discontinuous ephemeral streams. The

basic two-step model evolution is a robust feature

w0h0

»=π (7.36) of the model and was observed for a wide range

S0

of model parameters.

and quickly becomes negative (i.e., perturbations As the instability develops from random ini-

decay to zero) for larger wavelengths. The growth tial conditions, the dominant wavelength in-

rate is a function of κ (i.e., channels with larger creases, and the channel geometry becomes

values of κ develop oscillations more rapidly), but more regularly periodic. The increase in oscilla-

the wavelength corresponding to the maximum tion wavelength indicates that nonlinear, ¬nite-

growth rate (Eq. (7.36)) is independent of κ. This amplitude effects are an important part of the

independence is important for testing Eq. (7.36) model behavior. This result is con¬rmed by Fig-

against ¬eld measurements because κ is not well ure 7.11b, which compares the oscillation wave-

constrained, but the average channel width, bank lengths predicted by the linear-stability analy-

height, and slope can be readily measured for any sis and the direct numerical solution for dif-

ferent values of (w 0 h 0 /S0 )1/2 . The lower solid

channel.

The two-step Lax--Wendroff scheme (e.g., Press line corresponds to the linear-stability predic-

et al., 1992) was used to integrate Eqs. (7.28) and tion (Eq. (7.36)), and the ¬lled squares and up-