. 1
( 14)


This page intentionally left blank
Quantitative Modeling of Earth Surface Processes
Assuming some knowledge of calculus and
Geomorphology is undergoing a renaissance made
basic programming experience, this quantitative
possible by new techniques in numerical modeling,
textbook is designed for advanced geomorphology
geochronology and remote sensing. Earth surface
courses and as a reference book for professional
processes are complex and richly varied, but analyt-
researchers in Earth and planetary sciences look-
ical and numerical modeling techniques are pow-
ing for a quantitative approach to earth surface
erful tools for interpreting these systems and the
processes. Exercises at the end of each chapter
landforms they create.
begin with simple calculations and then progress
This textbook describes some of the most effec-
to more sophisticated problems that require com-
tive and straightforward quantitative techniques for
puter programming. All the necessary computer
modeling earth surface processes. By emphasizing a
codes are available online at www.cambridge.org/
core set of equations and solution techniques, the
book presents state-of-the-art models currently em-
ployed in earth surface process research, as well
Jon Pelletier was awarded a Ph.D. in geological sci-
as a set of simple but practical tools that can be
ences from Cornell University in 1997. Following
used to tackle unsolved research problems. Detailed
two years at the California Institute of Technology
case studies demonstrate application of the meth-
as the O.K. Earl Prize Postdoctoral Scholar, he was
ods to a wide variety of processes including hill-
made an associate professor of geosciences at the
slope, ¬‚uvial, eolian, glacial, tectonic, and climatic
University of Arizona where he teaches geomorphol-
systems. The computer programming codes used in
ogy. Dr Pelletier™s research involves mathematical
the case studies are also presented in a set of ap-
modeling of a wide range of surface processes on
pendices so that readers can readily utilize these
Earth and other planets, including the evolution
methods in their own work. Additional references
of mountain belts, the transport and deposition of
are also provided for readers who wish to ¬ne-
dust in arid environments, and ¬‚uvial and glacial
tune their models or pursue more sophisticated
processes on Mars.
Quantitative Modeling of
Earth Surface Processes
Jon D. Pelletier
University of Arizona
Cambridge, New York, Melbourne, Madrid, Cape Town, Singapore, São Paulo

Cambridge University Press
The Edinburgh Building, Cambridge CB2 8RU, UK
Published in the United States of America by Cambridge University Press, New York
Information on this title: www.cambridge.org/9780521855976

© J. D. Pelletier 2008

This publication is in copyright. Subject to statutory exception and to the provision of
relevant collective licensing agreements, no reproduction of any part may take place
without the written permission of Cambridge University Press.
First published in print format 2008

ISBN-13 978-0-511-42310-9 eBook (EBL)

ISBN-13 978-0-521-85597-6 hardback

Cambridge University Press has no responsibility for the persistence or accuracy of urls
for external or third-party internet websites referred to in this publication, and does not
guarantee that any content on such websites is, or will remain, accurate or appropriate.

Preface page ix

Chapter 1 Introduction 1

1.1 A tour of the ¬‚uvial system 1
1.2 A tour of the eolian system 12
1.3 A tour of the glacial system 19
1.4 Conclusions 29

Chapter 2 The diffusion equation 30

2.1 Introduction 30
2.2 Analytic methods and applications 34
2.3 Numerical techniques and applications 57
Exercises 63

Chapter 3 Flow routing 66

3.1 Introduction 66
3.2 Algorithms 66
3.3 ˜˜Cleaning up™™ US Geological Survey DEMs 70
3.4 Application of ¬‚ow-routing algorithms to estimate
¬‚ood hazards 72
3.5 Contaminant transport in channel bed sediments 74
Exercises 85

Chapter 4 The advection/wave equation 87

4.1 Introduction 87
4.2 Analytic methods 88
4.3 Numerical methods 90
4.4 Modeling the ¬‚uvial-geomorphic response of the southern
Sierra Nevada to uplift 93
4.5 The erosional decay of ancient orogens 101
Exercises 107

Chapter 5 Flexural isostasy 109

5.1 Introduction 109
5.2 Methods for 1D problems 111
5.3 Methods for 2D problems 113
5.4 Modeling of foreland basin geometry 116
5.5 Flexural-isostatic response to glacial erosion in the
western US 120
Exercises 123

Chapter 6 Non-Newtonian ¬‚ow equations 125

6.1 Introduction 125
6.2 Modeling non-Newtonian and perfectly plastic ¬‚ows 125
6.3 Modeling ¬‚ows with temperature-dependent viscosity 130
6.4 Modeling of threshold-sliding ice sheets and glaciers over
complex 3D topography 132
6.5 Thrust sheet mechanics 147
6.6 Glacial erosion beneath ice sheets 149
Exercises 160

Chapter 7 Instabilities 161

7.1 Introduction 161
7.2 An introductory example: the Rayleigh--Taylor instability 162
7.3 A simple model for river meandering 164
7.4 Werner™s model for eolian dunes 166
7.5 Oscillations in arid alluvial channels 169
7.6 How are drumlins formed? 174
7.7 Spiral troughs on the Martian polar ice caps 183
Exercise 187

Chapter 8 Stochastic processes 188

8.1 Introduction 188
8.2 Time series analysis and fractional Gaussian noises 188
8.3 Langevin equations 191
8.4 Random walks 193
8.5 Unsteady erosion and deposition in eolian environments 194
8.6 Stochastic trees and diffusion-limited aggregation 196
8.7 Estimating total ¬‚ux based on a statistical
distribution of events: dust emission from playas 199
8.8 The frequency-size distribution of landslides 205
8.9 Coherence resonance and the timing of ice ages 210
Exercises 221

Appendix 1 Codes for solving the diffusion
equation 222

Appendix 2 Codes for ¬‚ow routing 235

Appendix 3 Codes for solving the advection
equation 242

Appendix 4 Codes for solving the ¬‚exure equation 256

Appendix 5 Codes for modeling non-Newtonian
¬‚ows 263

Appendix 6 Codes for modeling instabilities 267

Appendix 7 Codes for modeling stochastic
processes 274

References 278
Index 290

The colour plates are to be found between pages 108 and 109.

Geomorphology is undergoing a renaissance gists, on the other hand, are adept at reading
made possible by new techniques in numeri- the geomorphic and sedimentary records to ad-
cal modeling, geochronology, and remote sens- dress big questions, but those records often can-
ing. Advances in numerical modeling make it not be fully interpreted using ¬eld observations
possible to model surface processes and their and geochronology alone. A fourth challenge is
feedbacks with climate and tectonics over a geomorphic prediction. In order for applied geo-
wide range of spatial and temporal scales. The morphology to realize its full potential, geomor-
Shuttle Radar Topography Mission (SRTM) has phologists must be able to predict where geo-
mapped most of Earth™s topography at much morphic hazards are most likely to occur, tak-
higher spatial resolution and accuracy than ever ing into account the full complexity of processes,
before. Cosmogenic dating and other geochrono- feedbacks, and the multi-scale heterogeneity of
logic techniques have provided vast new data on Earth™s surface.
surface-process rates and landform ages. Model- Analytical and numerical modeling are pow-
ing, geochronology, and remote sensing are also erful tools for addressing these challenges. First,
revolutionizing natural-hazard assessment and modeling is useful for establishing relationships
mitigation, enabling society to assess the hazards between process and form. Using quantitative
posed by ¬‚oods, landslides, windblown dust, soil models for different processes, modeling allows
erosion, and other geomorphic hazards. us to determine the signatures of those processes
The complexity of geomorphic systems poses in the landscape. In some cases, the histories of
several challenges, however. First, the relation- external forcing mechanisms (e.g. climate and
ship between process and form is often dif¬cult tectonics) can also be inferred. These linkages are
to determine uniquely. Many geomorphic pro- important because there is generally no direct
cesses cannot be readily quanti¬ed, and it is often observation of landforms that enables us to un-
unclear which processes are most important in derstand how they evolved. Modeling has played
controlling a particular geomorphic system, and a signi¬cant role in recent contributions to our
how those processes interact to form the geomor- understanding of many classic landform types.
phic and sedimentary records we see today. Ter- In this book, I will present analytic and numer-
races and sedimentary deposits on alluvial fans, ical models for many different landform types,
for example, are controlled by climate, tectonics, including drumlins, sand dunes, alluvial fans,
and internal drainage adjustments in a way that and bedrock drainage networks, just to name a
geomorphologists have not been able to fully un- few. Modeling is particularly useful for explor-
ravel. Second, surface processes are strongly in¬‚u- ing the feedbacks between different components
enced by ¬‚uid motions, and most classic geomor- of geomorphic systems (i.e. hillslopes and chan-
phic techniques (e.g. ¬eld mapping) are not well nel networks). Channel aggradation and incision,
suited to quantifying ¬‚uid dynamics and their for example, controls the base level of hillslopes,
interactions with the surface. Third, the geomor- and hillslopes supply primary sediment ¬‚ux to
phic community must bridge the gap between channels. Yet many geomorphic textbooks treat
process-based geomorphology and Quaternary ge- these as essentially independent systems. In some
ology. Process-based geomorphologists have made cases, the boundary conditions posed by one as-
great strides in quantifying transport and erosion pect of a geomorphic system can be considered
laws for geomorphic systems, but this approach to be ¬xed. In many of the most interesting un-
has not yet led to major advances in big geologic solved problems, however, they cannot. An un-
questions, such as how quickly the Grand Canyon derstanding of ¬‚uvial-system response to climate
was carved, or how mountain belts respond to change, for example, cannot be achieved with-
glacial erosion, for example. Quaternary geolo- out a quantitative understanding of the coupled

evolution and feedbacks between hillslopes and behavior and solution methods for diffusion, the
channels. reader will be poised to understand diffusive phe-
Numerical modeling is also useful in geomor- nomena as they arise in different contexts. Sim-
phology because of the central role that ¬‚uid dy- ilarly, boundary-layer and non-Newtonian ¬‚uid
namics play in landform evolution. Geomorphol- ¬‚ows are also examples of equation sets that
ogy can be roughly de¬ned as the evolution of have broad applicability in Earth surface pro-
landforms by the ¬‚uid ¬‚ow of wind, liquid wa- cesses. Non-Newtonian ¬‚uid ¬‚ows are the basis
ter, and ice above the surface. Yet, ¬‚uid mechan- for understanding and modeling lava ¬‚ows, de-
ics generally plays a minor role in many aspects bris ¬‚ows, and alpine glaciers. In the course of
of geomorphic research. Many landform types developing my quantitative skills and applying
evolve primarily by a dynamic feedback mech- them to a wide range of research problems, I
anism in which the topography in¬‚uences the have been continually amazed at how often a
¬‚uid ¬‚ow (and hence the shear stresses) above the concept or technique from one area can pro-
topography, which, in turn, controls how the to- vide the missing piece required to solve a long-
pography evolves by erosion and deposition. Nu- standing problem in another ¬eld. It is this cross-
merical modeling is one of the most successful fertilization of ideas and methods that I want
ways to quantify ¬‚uid ¬‚ow in a complex environ- to foster and share with readers. By emphasizing
ment, and hence it can and should play a central a core set of equations and solution techniques,
role in nearly all geomorphic research. In this readers will come away with powerful tools that
book, numerical models of ¬‚uid ¬‚ow serve as the can be used to tackle unsolved problems, as well
basis for many of the book™s numerical landform as speci¬c knowledge of state-of-the-art models
evolution models. Despite the power of numeri- currently used in Earth surface process research.
cal models, they cannot be used in a vacuum. All The book is designed for use as a textbook for an
numerical modeling studies must integrate ¬eld advanced geomorphology course and as a ref-
observations, digital geospatial data, geochronol- erence book for professional geomorphologists.
ogy, and small-scale experiments to motivate, Mathematically, I assume that readers have mas-
constrain, and validate the numerical work. tered multivariable calculus and have had some
Earth surface processes are complex and experience with partial differential equations.
richly varied. Most books approach the inher- The exercises at the end of each chapter begin
ent variety of geomorphic systems by serially with simple problems that require only the main
cataloging the processes and landforms charac- concepts and a few calculations, then progress to
teristic of each environment (hillslope, ¬‚uvial, more sophisticated problems that require com-
eolian, glacial, etc.). The disadvantage of this ap- puter programming.
proach is that each geomorphic environment is The purpose of this book is not to provide
presented as being essentially unique, and com- an exhaustive survey of all analytic or numeri-
mon dynamical behaviors are not emphasized. cal methods for a given problem, but rather to
This book follows a different path, taking advan- focus on the most powerful and straightforward
tage of a common mathematical framework to methods. For example, advective equations can
emphasize universal concepts. This framework fo- be solved using the upwind-differencing, Lax--
cuses on linear and nonlinear diffusion, advec- Wendroff method, staggered-leapfrog, pseudo-
tion, and boundary-value problems that quantify spectral, and semi-Lagrangian methods, just to
the stresses in the atmosphere and lithosphere. name a few. Rather than cover all available meth-
The diffusion equation, for example, can be used ods, this book focuses primarily on the simpler
to describe hillslope evolution, channel-bed evo- methods for readers who want to get started
lution, delta progradation, hydrodynamic disper- quickly or who need to solve problems of mod-
sion in groundwater aquifers, turbulent disper- est computational size. In this sense, the book
sion in the atmosphere, and heat conduction will follow the model of Numerical Recipes (Press
in soils and the Earth™s crust. By ¬rst provid- et al., 1992) by providing tools that work for
ing the reader with a solid foundation in the most applications, with additional references for

readers who want to ¬ne-tune their models. The 1D and when we model the evolution of surfaces
appendices provide the reader with computer we will refer to the model as 2D. This conven-
code to illustrate technique application in real- tion may seem strange at ¬rst, but it makes the
world research problems. Hence, many of the ap- book more consistent overall. For example, if we
plications I cover in this book necessarily come model heat ¬‚ow in a thin layer using the diffu-
from my own research. In focusing on my own sion equation, everyone would agree that we are
work, I don™t mean to imply that my work is solving a 2D problem for T (x, y, t). If we use the
the only or best approach to a given problem. diffusion equation to model the evolution of a
The book is not intended to be a complete survey hillslope described by h(x, y, t), that, too, is math-
of geomorphology, and I knowingly have left out ematically a 2D problem even though it represents
many important contributions in favor of a more a 3D landform.
focused, case-study approach. I wish to acknowledge my colleagues at the
Modelers are not always consistent in the way University of Arizona for taking a chance on
that they use the terms one- (1D), two- (2D), an unconventional geomorphologist and for cre-
and three-dimensional (3D) when referring to a ating such a collegial working environment. I
model. Usually when physicists and mathemati- also wish to thank my graduate students Leslie
cians use the term 1D they are referring to a Hsu, Jason Barnes, James Morrison, Steve De-
model that has one independent spatial variable. Long, Michael Cline, Joe Cook, Maria Banks, Joan
This can be confusing when applied to geomor- Blainey, Jennifer Boerner, and Amy Rice for help-
phic problems, however, because a model for the ful conversations and collaborations. I hope they
evolution of a topographic pro¬le, h(x, t), for ex- have learned as much from me as I have from
ample, would be called 1D even though it repre- them. Funding for my research has come from
sents a 2D pro¬le. Similarly, a model for the evo- the National Science Foundation, the Army Re-
lution of a topographic surface h(x, y, t) would be search Of¬ce, the US Geological Survey, the NASA
described as 2D even though the surface itself is Of¬ce of Space Science, the Department of En-
three-dimensional. In this book we will use the ergy, and the State of Arizona™s Water Sustain-
convention that the dimensionality of the prob- ability Program. I gratefully acknowledge support
lem refers to the number of independent spa- from these agencies and institutions. Finally, I
tial variables. Therefore, when we model 2D to- wish to thank my wife Pamela for her patience
pographic pro¬les we will classify the model as and support.
Chapter 1


and Jones, 1999). The structural style of exten-
1.1 A tour of the ¬‚uvial system sional faulting varies greatly from place to place
in the western US, but Figure 1.2 illustrates one
simple model of extension that can help us un-
The ¬‚uvial system is classically divided into
derstand the topography of the modern Basin
erosional, transportational, and depositional
and Range. The weight of the high topography
regimes (Schumm, 1977). In the Basin and Range
in the region created a vertical compressive force
province of the western US, Cenozoic tectonic
on the lower crust. That compressive force was
extension has produced a semi-periodic topog-
accompanied by an extensional force in the hori-
raphy with high ranges (dominated by erosion)
zontal direction. Mechanically, all rocks respond
and low valleys (dominated by deposition). In this
to a compressive force in one direction with an
region, all three ¬‚uvial-system regimes can be
extensional force in the other direction in or-
found within distances of 10--20 km. As an in-
der to preserve volume. This combination of com-
troduction to the scienti¬c questions addressed
pressive vertical forces and extensional horizon-
in this book, we start with a tour of the pro-
tal forces created faults at angles of ≈ 5--30—¦ to
cess zones of the ¬‚uvial system, using Hanau-
the horizontal (with smaller angles at greater
pah Canyon, Death Valley, California, as a type
example (Figure 1.1).
Crustal blocks ¬‚oat on the mantle, with the
1.1.1 Large-scale topography of the basin weight of the crust block supported by the buoy-
and range province ancy force that results from the bottom of the
crustal block displacing higher-density mantle
The large-scale geomorphology of the Basin and
rock. The shape of each crustal block partially
Range is a consequence of the geometry of faults
determines how high above Earth™s surface the
that develop during extension and the subse-
block will stand. Normal faulting creates a series
quent isostatic adjustment of each crustal block.
of trapezoidal crustal blocks that, while mechan-
In the late Cretaceous, the Basin and Range was
ically coupled by faults, can also rise and fall ac-
an extensive high-elevation plateau, broadly simi-
cording to their shape. A crustal block with a
lar to the Altiplano-Eastern Cordillera of the cen-
wide bottom displaces a relatively large volume
tral Andes today. When subduction of the Far-
of high-density mantle. The resulting buoyancy
allon plate ceased beneath the western US, this
requires that the block stand higher above Earth™s
region became a predominantly strike-slip plate
geoid in order to be in isostatic equilibrium. Con-
boundary and the horizontal compressive force
versely, a crustal block with a relatively narrow
that supported the high topography and over-
bottom will displace a smaller volume of mantle
thickened crust of the region could no longer
and will hence stand lower relative to the geoid.
withstand the weight of the overlying topogra-
This kind of adjustment is the fundamental
phy. The result was regional extension (Sonder

Fig 1.1 Hanaupah Canyon and its alluvial fan, Death Valley,

Fig 1.3 Oblique perspective image of the high peaks of
Hanaupah Canyon. Reproduced with permission of
‚h = 0

divide h

Fig 1.2 Schematic diagram of a simple model for Basin and
h regolith
Range extension and isostatic adjustment. In (a), horizontal
extension and vertical compression cause the creation of
trapezoidal crustal blocks bounded by a series of normal
faults. Crustal blocks with wide bottoms (ranges) rise to
maintain isostatic balance, while blocks with narrow bottoms
(basins) fall relative to Earth™s geoid.

reason why extension causes the semi-periodic to-
pography of the Basin and Range Province. channel head
Fig 1.4 Schematic diagram of a hillslope pro¬le from divide
1.1.2 The hillslope system to channel head.
Hillslopes in the high elevations of Hanaupah
Canyon are characterized by steep, planar topo-
graphic pro¬les and abrupt, knife-edge divides because of this change in drainage direction. At
(Figure 1.3). West of the divide, water and sedi- the divide itself, the water and sediment ¬‚ux is
ment drains towards Panamint Valley, while on zero, since this is the transition point from west-
the east side they drain towards Death Valley. If ern drainage to eastern drainage. Sediment ¬‚ux
we quantify the ¬‚ux of water and sediment along on hillslopes is directly related to the topographic
a pro¬le that includes the drainage divide, the gradient. Therefore, since the sediment ¬‚ux is
zero, the topographic gradient, ‚h/‚ x, must also
¬‚ux of water and sediment must be negative on
one side of the divide and positive on the other be zero (Figure 1.4).

Hillslopes in bedrock-dominated landscapes knife-edge drainage divides of Hanaupah Canyon
are composed of a system of two interacting sur- are a signature of landslide-dominated, nonlin-
faces: the topographic surface, with elevations ear transport on hillslopes (Roering et al., 1999).
given by h(x, y), and the underlying weathering Equation (1.2) states that bedrock lowering
front, given by b(x, y) (Figure 1.4). The difference is a maximum for bare bedrock slopes and de-
between these two surfaces is the regolith depth, creases exponentially with increasing regolith
·(x, y). The topographic and weathering-front sur- thickness. This relationship has been inferred
faces are strongly coupled because the shape from cosmogenic isotope analyses on hillslopes
of the topography controls erosion and deposi- (Heimsath et al., 1997). Conceptually, the in-
tion, which, in turn, changes the values of ·(x, y) verse relationship between weathering/regolith-
(Furbish and Fagherazzi, 2001). The values of production rate and regolith thickness is a con-
·(x, y), in turn, control weathering rates. The sim- sequence of the fact that regolith acts as a
plest system of equations that describes this feed- buffer for the underlying bedrock, protecting it
back relationship is given by: from diurnal temperature changes and in¬ltrat-
ing runoff that act as drivers for physical and
‚· ρb ‚b ‚ 2h
=’ ’κ 2 (1.1) chemical weathering in the subsurface. In Figure
‚t ρs ‚t ‚x
1.4, the weathering front is shown as an abrupt
and transition from bedrock to regolith, but in nature
‚b this boundary is usually gradual.
= ’P 0 e’·/·0 (1.2)
‚t Equations (1.1) and (1.2) can be solved for the
where t is time, ρb is the bedrock density, ρs steady-state case in which regolith thickness is
is the sediment density, κ is the hillslope diffu- independent of time at all points along the hill-
sivity, x is the distance along the hillslope pro- slope pro¬le:
¬le, P 0 is the regolith-production rate for bare
ρb P 0 1
bedrock, and ·0 is a characteristic regolith depth. · = ·0 ln (1.3)
ρs κ ’ ‚ 2 h
Equation (1.1) states that the rate of change ‚ x2
of regolith thickness with time is the differ- Equation (1.3) implies that, in steady state, re-
ence between a ˜˜source™™ term equal to the rate golith thickness decreases as the negative (down-
of bedrock lowering multiplied by the ratio of ward) curvature increases, becoming zero (i.e.
bedrock to sediment density, and a ˜˜sink™™ term bare bedrock) where the curvature reaches a crit-
equal to the curvature of the topographic pro- ical value of
¬le. This curvature-based erosion model is the
‚ 2h ρb P 0
classic diffusion model of hillslopes, ¬rst pro- (1.4)
‚ x2 ρs κ
bar e
posed by Culling (1960). The diffusion model of
hillslope evolution, discussed in Chapter 2, is a Equation (1.2) is not universally applicable.
consequence of conservation of mass along hill- As regolith thickness decreases below a critical
slope pro¬les, and the fact that sediment ¬‚ux is value in arid regions, for example, the land-
proportional to topographic gradient if certain scape is unable to store enough water to sup-
conditions are met. The diffusion model of hill- port signi¬cant plant life. Plants act as weath-
slope evolution does not apply to the steep hill- ering agents (e.g. root growth can fracture rock,
slopes of Hanaupah Canyon, however. In steep canopy cover can decrease evaporation, etc.).
landscapes, sediment ¬‚ux increases nonlinearly Therefore, in some arid environments weather-
with the topographic gradient as the angle of sta- ing rates actually increase as regolith thickness
bility is approached and landsliding becomes the increases, rather than decreasing with thickness
predominant mode of sediment transport. The ef- as the exponential term in Eq. (1.2) predicts
fects of mass movements on hillslope evolution (Figure 1.5). Using a vegetation-limited model of
can be captured in a nonlinear hillslope diffu- weathering in landscape evolution models re-
sion model that will also be discussed in Chap- sults in landscapes with a bimodal distribution of
ter 2. The steep, planar hillslopes and abrupt, slopes (i.e. cliffs and talus slopes) similar to many

channel bed and channel bank, and the channel
P0 bed and shoreline, just to name a few.
P0e’h/h0 The transition between the hillslope and ¬‚u-
vial channel system occurs where the shear stress
‚h of overland ¬‚ow is suf¬cient to entrain hillslope
‚t material. In addition, the rate of sediment exca-
vation from the channel head by overland ¬‚ow
must be greater than the rate of sediment in¬ll-
ing by creep and other hillslope processes. Empir-
ically, this transition occurs where the product
h of the topographic gradient and the square root
of contributing area is greater than a threshold
Fig 1.5 Models for the relationship between regolith
value (Montgomery and Dietrich, 1999):
production and regolith thickness, illustrating the exponential
model of Heimsath et al. (1997) and the alternative 1
S A 2 = X ’1 (1.5)
vegetation-limited model of Anderson and Humphrey (1990).
where S is the topographic gradient or slope, A is
The density contrast between the bedrock and sediment is
assumed to be zero for simplicity. the contributing area, and X is the drainage den-
sity (i.e. equal to the ratio of the total length of all
the channels in the basin to the basin area). The
value of X depends on the texture and perme-
arid-region hillslopes (Anderson and Humphrey,
ability of the regolith, hillslope vegetation type
1990; Strudley et al., 2006). It should also be
noted that the parameters P 0 and κ may be ap- and density, and on the relative importance of
different hillslope processes.
proximately uniform along some hillslope pro-
¬les, but in general the hillslope and regolith-
1.1.3 Bedrock channels
thickness pro¬les will in¬‚uence the hillslope hy-
drology (e.g. the ratio of surface to subsurface Channels are divided into alluvial and bedrock
¬‚ow), which will, in turn, modify the values of P 0 channel types depending on whether or not al-
and κ through time. Thinner regolith, for exam- luvium is stored on the channel bed. In bedrock
ple, can be expected to increase surface runoff, channels, the transport capacity of the channel
thereby increasing κ and decreasing P 0 in a pos- is greater than the upstream sediment ¬‚ux. Sed-
itive feedback. iment storage is rare or nonexistent in these
Equations (1.1) and (1.2) suggest that even cases. In alluvial channels, upstream sediment
in well-studied geomorphic systems such as hill- ¬‚ux is greater than the transport capacity of the
slopes, quantifying the behavior of even the most channel and alluvium ¬lls the channel bed as a
basic elements is a signi¬cant challenge. This result. The distinction between bedrock and al-
challenge provides an opportunity for mathemat- luvial channels is important for understanding
ical modeling to play an essential role in ge- how they evolve. In order for bedrock channels
omorphic research, however, because most ge- to erode their beds they must pluck or abrade
omorphic systems are too complex to be fully rock from the bed. Once rock material is eroded
understood and interpreted with ¬eld observa- from a bedrock river it is usually transported
tions, measurements, and geochronologic tech- far from the site of erosion because the ¬‚ow ve-
niques alone. The coupled hillslope evolution locity needed to transport material is typically
model of Eqs. (1.1) and (1.2) also illustrates a much lower than the velocity required for pluck-
larger point: many geomorphic systems of great- ing or abrasion. Bedrock channels erode their
est interest involve the coupling of different pro- beds by a combination of plucking, abrasion, and
cess domains (in this case the hillslope weather- cavitation. Plucking occurs when the pressure
ing and sediment transport regimes). Other key of fast-¬‚owing water over the top of a jointed
process-domain linkages occur at the juncture be- rock causes suf¬cient Bernoulli lift to dislodge
tween the hillslope and the channel head, the rock from the bed. Abrasion occurs as saltating

erosion is usually linear with slope (i.e. n ≈ 1 and
bedload impacts the bed, chipping away small
m ≈ 1/2 (Kirby and Whipple, 2001)).
pieces of rock. Cavitation occurs when water boils
under conditions of very high pressure and shear In steady state, uplift and erosion are in bal-
ance and ‚h/‚t = 0. Equation (1.6) predicts that
stress. Imploding bubbles in the water create
channel slope, S = |‚h/‚ x|, is inversely propor-
pressures suf¬cient to pulverize the rock. Con-
ditions conducive to cavitation are most likely tional to a power-law function of drainage area,
to occur in rare large ¬‚oods. The relative impor- A, for drainage basins in topographic steady state:
tance of plucking and abrasion depend primar- U ’m/n
S= A (1.7)
ily on lithology. Abrasion is generally considered
to be the dominant process in massive bedrock
Equation (1.7) predicts that, as drainage area in-
lithologies such as granite (Whipple et al., 2000).
creases downstream, the channel slope must de-
In sedimentary rocks, abrasion and plucking are
crease in order to maintain uniform erosion rates
both likely to be important.
across the landscape. This inverse relationship be-
Several mathematical models exist for mod-
tween slope and area is responsible for the con-
eling bedrock channel evolution. In the stream-
cave form of most bedrock channels.
power model, ¬rst proposed by Howard and Kerby
Equation (1.7) can be used to derive an ana-
(1983), the rate of increase or decrease in chan-
lytic expression for the steady-state bedrock chan-
nel bed elevation is equal to the difference be-
nel longitudinal pro¬le. To do this, we must re-
tween the uplift rate and the erosion rate. The
late the drainage area, A, to the distance along
erosion rate is a power function of drainage area
the channel from the channel head, x. The re-
and channel-bed slope:
lationship between area and distance depends
‚h ‚h on the planform (i.e. map-view) shape of the
= U ’ K w Am (1.6)
‚t ‚x basin. For a semi-circular basin, the square root
of drainage area is proportional to the distance
where h is the local elevation, t is time, U is
from the channel head. Adopting that assump-
uplift rate, K w is a constant that depends on
tion together with the empirical observation
bedrock erodibility and climate, A is the drainage
m/n = 0.5, Eq. (1.7) becomes
area, m and n are empirical constants, and x
is the distance along the channel. The general ‚h C
=’ (1.8)
form of the stream-power model also includes ‚x x
an additional constant term that represents a
where ‚h/‚ x has been substituted for S and C is a
threshold shear stress for erosion (Whipple and
constant that combines the effects of uplift rate,
Tucker, 1999). The stream-power model is empiri-
bedrock erodibility, climate, and basin shape. In-
cally based: Howard and Kerby (1983) found that
tegrating Eq. (1.8) gives
Eq. (1.6) successfully reproduced observed erosion
rates at a site in Perth Amboy, New Jersey, mea- h ’ h 0 = ’C ln (1.9)
sured by repeat survey. Their analysis cannot,
however, rule out whether other variables corre- where h 0 is the elevation of the base of the chan-
nel at x = L . Figure 1.6 shows the longitudinal
lated with drainage area and slope are the funda-
mental controlling variables of bedrock erosion. pro¬le of the main channel of Hanaupah Canyon
Nevertheless, the stream power is based on the together with a plot of Eq. (1.9). The model ¬ts the
observed pro¬le quite well using C = 0.9 km.
physically-reasonable assumption that the ero-
sive power of ¬‚oods increases as a function of The evolution of a hypothetical bedrock chan-
drainage area (which controls the volume of wa- nel governed by the stream-power model is illus-
ter routed through the channel) and channel-bed trated qualitatively in Figure 1.7a. A small reach
slope (which controls the velocity of that water of the channel™s longitudinal pro¬le is shown at
for otherwise similar drainage areas and precip- times t1 , t2 , and t3 . If the section is relatively
itation intensities). Calibrations of the stream- small and has no major incoming tributaries,
power model in natural systems suggest that the drainage area can be considered uniform

ing in the propagation of knickpoints upstream
(a) as topographic waves. The rate of knickpoint
propagation is equal to K w A m , i.e. knickpoints
propagate faster in larger channels, wetter cli-
mates, and areas of more easily-eroded bedrock.
The stream-power model will serve as the type
example of the advection/wave equation studied
in Chapter 4.
The Kern River (Figures 1.8 and 1.9) pro-
vides a nice example of knickpoint propaga-
3 tion in action. Two distinct topographic surfaces
have long been recognized in the landscape of
the southern Sierra Nevada (Webb, 1946) (Figure
1.8). The Boreal surface is a high-elevation, low-
h (km)

longitudinal profile
relief plateau that dips to the west at 1—¦ (Fig-
model fit
ure 1.8b). The Chagoopa Plateau is an interme-
1 diate topographic ˜˜bench™™ that is restricted to
the major river canyons and inset into the Bo-
real Plateau (Webb, 1946; Jones, 1987). Figure 1.8b
bedrock alluvial
0 maps the maximum extents of the Chagoopa
0 20
5 10 15
and Boreal Plateaux based on elevation ranges
x (km)
of 1750--2250 m (Chagoopa) and 2250--3500 m a.s.l.
Fig 1.6 (a) Shaded-relief image of the Hanaupah Canyon (Boreal). Associated with each surface are promi-
drainage network and alluvial fan. Location of longitudinal nent knickpoints along major rivers. Knickpoints
pro¬le shown as white curve. (b) Longitudinal pro¬le of main along the North Fork Kern River, for example
Hanaupah Canyon channel, together with best-¬t to Eq. (1.8).
(Figure 1.9b), occur at elevations of 1600--2100 m
and 2500--3300 m a.s.l. The stepped nature of the
Sierra Nevada topography is generally considered
to be the result of two pulses of Cenozoic and/or
bedrock late Cretaceous uplift (Clark et al., 2005; Pelletier,
channel 2007c). According to this model, two major knick-
t3 t2 t1 points were created during uplift, each initiating
a wave of incision that is still propagating head-
(b) ward towards the range crest.
Recent work has highlighted the importance
t1 t
2 t3 of abrasion in controlling bedrock channel evo-
lution. In the abrasion process it is sediment, not
water, that acts as the primary erosional agent.
In the stream-power model, the erosive power is
Fig 1.7 Schematic diagrams of the evolution of (a) bedrock
assumed to be a power function of drainage area.
and (b) alluvial channels through time, illustrating the
advective behavior of bedrock channels and the diffusive Although sediment ¬‚ux increases with drainage
behavior of alluvial channels. area, upstream relief also plays an important role
in controlling sediment ¬‚ux. As such, the stream-
power model does not adequately represent the
throughout the reach. If drainage area is uni- abrasion process. Sklar and Dietrich (2001, 2004)
form and n = 1 is assumed, the erosion rate developed a saltation-abrasion model to quan-
is proportional to the channel slope according tify this process of bedrock channel erosion. In-
to the stream-power model. Accordingly, steeper sights into their model can be gained by re-
portions of the bed will erode faster, result- placing drainage area with sediment ¬‚ux in the

118 W
119 W along-dip profile
(b) 1o
profile in (b)
0 50
x (km)
San Joaquin R.
profile in (b)
37 N

inset in (c)
Kings R.

Kaweah R.

along-strike profile
36 N
Kern R.

100 150
0 50
x (km)
3.0 Boreal
profile in (f) (km)
2.0 Chagoopa
Kern R.
x (km)
0 30
Chagoopa (f)
3.0 knickpoints
Kern R.
profile in (e)
slope map
x (km)
0 100
S (m/m)
0.0 0.5 >1.0
3.0 4.0
1.0 2.0

Fig 1.8 Major geomorphic features of the southern Sierra Nevada. (a) Shaded relief map of topography indicating major rivers
and locations of transects plotted in b. (b) Maximum extents of the Chagoopa and Boreal Plateaux based on elevation ranges of
1750“2250 m and 2250“3500 m a.s.l. Also shown are along-strike and along-dip topographic transects illustrating the three levels
of the range in the along-strike pro¬le (i.e. incised gorges, Chagoopa and Boreal Plateaux) and the westward tilt of the Boreal
Plateau in the along-strike transect. (c) and (d) Grayscale map of topography (c) and slope (d) of the North Fork Kern River,
illustrating the plateau surfaces (e) and their associated river knickpoints (f). For color version, see plate section. Modi¬ed from
Pelletier (2007c). Reproduced with permission of Elsevier Limited.

stantial stream power drain the plateau. In the
(a) sediment-¬‚ux-driven model, however, low erosion
Boreal Plateau
rates on the plateau limit the supply of sediment
that acts as cutting tools to channels draining
the plateau edge. As a result, the sediment-¬‚ux-
driven model leads to much slower rates of knick-
point retreat than those predicted by the stream-
power model in a plateau-type landscape (Gas-
parini et al., 2006).
The sediment-¬‚ux-driven model also implies
Kern R. that the evolution of hillslopes and channels
is more intimately linked than the stream-
Boreal Plateau power model would suggest. In the stream-
power model, hillslope evolution plays no ex-
plicit role in bedrock channel evolution. In the
sediment-¬‚ux-driven model, however, increased
hillslope erosion will supply more cutting tools
to bedrock channels downstream. Downstream
channels will respond to this increased supply
with faster incision, further lowering the base
level for hillslopes in a positive feedback. These
ideas will be further explored in Chapter 4 where
the stream-power and sediment-¬‚ux-driven mod-
els are compared in detail.
Fig 1.9 Virtual oblique aerial photographs of portions of
As one travels through the ¬‚uvial system, the
the Kern River basin. (a) The upper portion of the basin is
transport capacity of the main channel gener-
characterized by the low-relief Boreal Plateau. The mainstem
ally increases with downstream distance as a re-
Kern River is incised 1“2 km into the Boreal Plateau. (b)
sult of increasing contributing area. The effect of
Relief between the Boreal Plateau and Kern River is
declining channel slope, however, partially off-
accommodated by a series of channel knickpoints along the
sets that effect. As a result, transport capacity
Kern River and its tributaries.
increases downstream, but not as rapidly as the
sediment supply that the channel is required to
stream-power model to obtain a sediment-¬‚ux-
transport. The bedrock-alluvial transition occurs
driven model:
where the sediment supply exceeds the transport
‚h ‚h
= U ’ KsQm capacity. In most drainage basins of the west-
‚t ‚x
ern US, the bedrock-alluvial transition occurs up-
where Q s is the sediment ¬‚ux and K s is a new stream of the mountain front (e.g. Figure 1.10).
coef¬cient of erodibility. Only in regions that are tectonically most active
The predictions of the stream-power and does the bedrock-alluvial transition occur at the
sediment-¬‚ux-driven models are similar for cases mountain front itself.
of uniformly uplifted, steady-state mountain
belts. In such cases, erosion is spatially uniform
1.1.4 Alluvial channels
(everywhere balancing uplift) and therefore sed-
iment ¬‚ux Q s is proportional to drainage area Alluvial channels evolve in a fundamentally dif-
A. The predictions of the two models are very ferent style than that of bedrock channels (Figure
different, however, following the uplift of an 1.7). Alluvial channel evolution is governed by a
initially low-relief plateau. In such cases, the conservation of mass relationship which states
stream-power model predicts a relatively rapid that the change in channel-bed elevation is equal
response to uplift because large rivers with sub- to the gradient in bedload sediment ¬‚ux along

load transport. If channel width is assumed to be
uniform along the longitudinal pro¬le, the com-
bination of Eqs. (1.11) and (1.12) gives the classic
diffusion equation (Begin et al., 1981):
‚h ‚ 2h
=κ 2 (1.13)
‚t ‚x
where the diffusivity is given by κ =
(BQb )/(c 0 w 0 ), and w 0 is the uniform channel
width. The diffusive evolution of alluvial chan-
nels of uniform width is illustrated schematically
in Figure 1.7b.

1.1.5 Alluvial fans
Deposition occurs on alluvial fans primarily be-
cause of channel widening near the mountain
Fig 1.10 Oblique perspective image of channels of front. Abrupt widening decreases the transport
Hanaupah Canyon. Smaller channels carve directly into capacity with distance downfan. Under such con-
bedrock, while the largest channels have wider beds ¬lled
ditions, Eq. (1.11) predicts aggradation and fan
with alluvium. Reproduced with permission of DigitalGlobe.
Hanaupah Canyon has one of the most spec-
tacular alluvial fans in the western US. The size
the channel pro¬le:
of the Hanaupah Canyon fan re¬‚ects the large
‚h 1 ‚(wqs )
=’ sediment ¬‚ux draining from Hanaupah Canyon,
‚t c0 ‚ x
which, in turn, is a result of the high relief
where h is the elevation of the channel bed, t is and semi-arid climate of this drainage basin.
time, c 0 is the volumetric concentration of bed The semi-arid climate maximizes runoff intensity
sediment, w is the channel width, qs is the sedi- while minimizing erosion-suppressing vegetation
ment discharge per unit channel width, and x is cover. Many alluvial fans in the western US, in-
the distance downstream. Equation (1.11) is sim- cluding Hanaupah Canyon fan, have a series of
ply a statement of conservation of mass, i.e. the distinct terraces (Figure 1.11) that rise like a ¬‚ight
channel bed must aggrade if the sediment-¬‚ux of stairs from the active channel, with older ter-
gradient is negative (if more sediment enters the races occurring at higher topographic positions
reach from upstream than leaves it downstream) relative to the active channel. Alluvial-fan ter-
and incise if the sediment-¬‚ux gradient is posi- races in many areas of the western US can be cor-
tive (if more sediment leaves the reach than en- related on the basis of their elevation above the
ters it). A number of different relationships exist active channel, their degree of desert pavement
for quantifying bedload sediment ¬‚ux in alluvial and varnish development, and their extent of
channels, but one common approach quanti¬es degradation. Terraces form as a result of changes
sediment discharge as a linear function of chan- in the ratio of sediment supply to transport ca-
nel gradient and a nonlinear function of the dis- pacity through time. During times when the ratio
charge per unit channel width: of sediment supply to transport capacity is high,
aggradation and channel widening operate in a
qs = ’B positive feedback that grows the fan vertically
and radially. During times when the ratio of sedi-
where B is a mobility parameter related to grain ment supply to transport capacity is low, incision
size, Q is water discharge, and b is a constant. and channel narrowing operate in a positive feed-
The value of b is constrained by sediment rating back that causes the channel to entrench into
curves and is generally between 2 and 3 for bed- older fan deposits, leaving an abandoned terrace.

main channels on the fan branch into dozens of
distributary channels separated by horizontal dis-
tances of only a few hundred meters or less.
The triggering mechanisms for alluvial-fan
terrace formation have long been debated, but
climatic changes most likely play a signi¬cant
role in controlling the changes in sediment sup-
ply that trigger fan aggradation and incision.
A growing database of surface and stratigraphic
age estimates suggests that Quaternary geomor-
phic surfaces and underlying deposits of the west-
ern US can be correlated regionally (Christensen
and Purcell, 1983; Bull, 1991; Reheis et al., 1996;
Bull, 1996). Several studies have documented ¬ll
events and/or surface exposure dates between
700--500 ka, 150--120 ka, 70--50 ka, and 15--10 ka
Fig 1.11 Oblique perspective image of the terraces of corresponding to the Q2a, Q2b, Q2c, and Q3 ge-
Hanaupah Canyon fan. Reproduced with permission of omorphic surfaces identi¬ed by Bull (1991) (sur-
DigitalGlobe. Younger terraces are lighter in color in this
face exposure ages correspond to youngest ages
image, representing the relatively limited time available for
in these ranges). Similarity of ages regionally has
desert pavement and varnish formation on younger surfaces.
provided preliminary support for the hypothesis
that Quaternary alluvial-fan terraces are gener-
Over time, episodes of aggradation, incision, and ated by climatic changes.
lateral reworking produce a nested sequence of Changes in sediment supply due to climatic
terraces. changes can result from several factors, but vari-
Multiple cut and ¬ll cycles on alluvial fans ations in drainage density are likely to play a
create a spatially-complex, distributary channel signi¬cant role in controlling the temporal vari-
network that presents a challenge to ¬‚oodplain ations in sediment ¬‚ux from drainage basins
managers in the western US. The Tortolita Moun- in the western US. Terrace formation during
tains fan northwest of Tucson, Arizona is a clas- the Pleistocene--Holocene transition is associated
sic example. Figure 1.12 presents four views of with a ten-fold increase in sediment supply (e.g.
this topographic complexity using a shaded re- Weldon, 1980). It is unlikely that a change in pre-
lief image, aerial photo, grayscale map of a nu- cipitation alone could account for such large in-
merical model of ¬‚ow depth during a recent ex- creases in sediment supply. During Pleistocene
treme ¬‚ood, and a sur¬cial geologic map. The sur- climates, vegetation densities were higher at
¬cial geologic map (Figure 1.12d) was constructed most elevations across the western US. Higher
by integrating soil development and other indi- vegetation density results in a lower drainage
cators of terrace age to group the terraces into density. During times of lower drainage density,
distinct age ranges (Gile et al., 1981; McFadden accommodation space is created in hollows for
et al., 1989). The surface age represents the ap- the deposition of sediment eroded from higher
proximate time since deep ¬‚ooding occurred on up on the hillslope, thereby lowering sediment
the terrace because soils would be stripped from ¬‚uxes from the basin relative to long-term geo-
a terrace subjected to deep scour and buried on a logic averages. During humid-to-arid transitions,
surface subjected to signi¬cant ¬‚uvial deposition. drainage densities increase, removing sediment
Sur¬cial geologic mapping indicates that ¬‚ood stored as colluvium in hollows during the pre-
risk (which is inversely correlated with surface vious humid interval. Sediment ¬‚uxes decrease
age) varies greatly even at scales less than 1 km. when the drainage density reaches a new maxi-
The modeled ¬‚ow depths also illustrate the spa- mum in equilibrium with the drier climate. Ac-
tial complexity of ¬‚ooding. As Figure 1.12c, the cording to this model, it is the change in climate,

111.10° W
111.15° W
32.50° N
1 km
Wild Burro


32.45° N

aerial photo
shaded relief

(c) Cochie Map units
(d) Qy1


surficial geology
model flow depth

Fig 1.12 (a) Shaded-relief image, (b) aerial photo, (c) ¬‚ow
phase with climate because it is responding to
depth model prediction, and (d) sur¬cial-geologic map of the
the rate of change of climate, not its absolute
Tortolita Mountains fan northwest of Tucson. Multiple cycles
state. The duration of the climate oscillation also
of cutting and ¬lling on this fan have resulted in a complex
plays a role in landscape response. Short-term cli-
distributary channel network in which ¬‚ooding from con¬ned
matic oscillations may not allow suf¬cient time
channels upstream branches into dozens of active channels
downstream on the fan. for sediment to be stored on the landscape during
times of lower drainage density. Long-term cli-
mate changes, however, allow time for more sed-
not the absolute state of climate (wet or dry), iment to be stored during low-drainage-density
that most strongly controls sediment ¬‚uxes from humid intervals, resulting in a larger response
drainage basins. when the climate shifts to more arid conditions.
Figure 1.13 schematically illustrates the rela- In the western US, the elevation zone of
tionships between precipitation, vegetation den- greatest vegetation change is well constrained
sity, and sediment ¬‚ux in the western US during through the Last Glacial Maximum (LGM) by
¬‚uctuating climates according to the ¬‚uctuating packrat middens. Today, the modern tree line (i.e.
drainage density model. Precipitation and vegeta- the location where shrub vegetation transitions
tion changes oscillate in phase over geologic time to mature woody species) occurs at approximately
scales. Sediment ¬‚ux, however, oscillates out of 2000 m a.s.l. (with variations of a few hundred

high frequency
low frequency





Fig 1.13 Schematic diagram of the relationships between
precipitation, vegetation density, and sediment ¬‚ux for
Fig 1.14 A HiRISE Camera image of erosional gullies on
low-frequency (long-duration) and high-frequency (short
Mars and their adjacent alluvial fans.
duration) climate changes.

scientists are currently grappling with many of
meters due to slope aspect). At the LGM, the tree
the same science questions that geomorpholo-
line was approximately 500 m a.s.l. In the western
gists puzzle over on Earth. Figure 1.14 presents
US, therefore, the elevational range from 500--
one of many images of the Martian surface that
2000 m has undergone signi¬cant changes in veg-
is tantalizingly similar to ¬‚uvial landforms on
etation type and density many times within the
Earth. Figure 1.14 shows a series of erosional val-
Plio--Quaternary period (Spaulding, 1990). It is
leys in the upper right-hand side of the image,
within this elevation range that we can expect to
transitioning into entrenched depositional fans
see the most signi¬cant impact of climate change
with multiple terrace levels. Whether these ter-
on ¬‚uvial drainage basin processes. Of course, at
races represent the effects of deposition during
much higher elevations some drainage basins in
dry or wet ¬‚ow events is not currently known.
the western US were also glaciated. Glaciation has
Either way, the similarity between many of the
a huge impact on ¬‚uvial processes downstream
¬‚uvial landforms on Mars and those on Earth is
by delivering large volumes of sediment due to
striking, especially given that liquid water is un-
the erosional ef¬ciency of alpine glaciers.
stable on Mars under the temperature and pres-
The ¬‚oor of Death Valley is a saline playa that
sure conditions of today. The wealth of imagery
has been an ephemeral shallow lake many times
and elevation data from Mars provides an excel-
in the past. Pluvial lakes occupied large areas in
lent opportunity to compare and contrast Mar-
the western US during much of the Pleistocene.
tian landforms to those on Earth in order to
During pluvial lake highstands, prominent shore-
learn more about how landforms evolved on both
lines were created by wave-cut action in many
areas of Utah (Bonneville shoreline) and Nevada
(Lahontan shoreline) (Reheis, 1999). Degradation
of these shorelines over time provides natural ex-
1.2 A tour of the eolian system
periments in hillslope evolution that have been
exploited by hillslope geomorphologists. Fluctu-
ations in lake levels also store and release vast The eolian system can be divided into two types:
amounts of windblown dust in arid regions, systems dominated by the transport of silt and
thereby in¬‚uencing soil formation and hydrology clay, and systems dominated by the transport of
on nearby alluvial-fan terraces (Reheis et al., 1995). sand. Both types of systems involve particle en-
One of the most exciting developments in trainment from the ground, atmospheric trans-
surface process research is the wealth of new port, and redeposition by gravitational settling.
data emerging on the surface of Mars. Planetary The primary difference between these two types

hydrological balance/
water table fluctuations
Franklin Playa
moisture-controlled entrainment l
turbulent dispersion 2
u*t u*td
= Penman“Montieth equation
Eagle M. r
u* u*
2 2

ponded water surface (occasional)
ground surface
h z=0
unsaturated zone
unsaturated flow
q q
Amargosa R. q q

Richards equation

Fig 1.15 Oblique perspective image of Franklin Lake Playa
and adjacent Eagle Mountain piedmont. Franklin Lake Playa
water table (moving boundary)
acts as a strong regional source of dust. Terraces on
z = zm
piedmonts such as that of Eagle Mountain act as net sinks for q = qs saturated zone
dust in arid environments. Reproduced with permission of
DigitalGlobe. Fig 1.16 Schematic ¬gure illustrating the coupling between
the hydrologic and eolian systems on playa surfaces.

of systems is the length and time scales involved
in particle motion. Silt and clay transport in- alluvial fan terraces of Eagle Mountain piedmont
volves travel distances on the order of 1 km or are sourced from Franklin Lake Playa. The dust of
greater and residence times on the order of 1 Franklin Lake Playa is, in turn, sourced from the
hour or more in the atmosphere. Sand, in con- Amargosa Valley drainage system, of which Eagle
trast, moves primarily in saltation. Saltation in- Mountain piedmont is a part. As the Amargosa
volves trajectories on the order of 10--1000 cm and River drains into Franklin Lake Playa, its low gra-
residence times on the order of 1 s. In this sec- dient and distributary channel geometry causes
tion, we will focus ¬rst on the transport of dust in most of the sediments it carries to settle out and
arid environments (i.e. silt and clay transport), us- be deposited on the playa. This dust is then sus-
ing Franklin Lake Playa in Amargosa Valley, Cali- ceptible to windblown transport when the playa
fornia, as a type example. Second, we will explore becomes dessicated. Some of this windblown dust
the transport of sand across an alluvial fan using is deposited immediately downwind (e.g. on the
Whiskey Petes™ alluvial fan as an example. Finally, alluvial fan terraces of Eagle Mountain alluvial
we will explore the formation mechanisms of eo- fan), while some is transported thousands of kilo-
lian dunes. meters away to the ocean or to a river that ulti-
mately drains into the ocean. As such, the dust
1.2.1 The dust cycle in arid environments cycle is not a closed system on the continents: at-
The dust cycle in arid environments refers to the mospheric transport results in a net transport of
transport of dust as it moves from sources (e.g. dust from continents to the ocean. Rock weather-
playas and dry channel beds) to sinks (e.g. alluvial ing continually generates new dust for this cycle.
fan terraces) and ultimately back to sources, pri- The amount of dust emitted from playas
marily by the action of ¬‚oods. Franklin Lake Playa varies greatly in space and time. The hydrologic
and the Amargosa River drainage basin (Figure state of a playa plays a signi¬cant role in con-
1.15) provide a classic example of the dust cycle in trolling this spatial and temporal variation in
action. Most of the dust that is deposited on the dust-emitting potential (Fig. 1.16) (Reynolds et al.,

2007). Generally speaking, wet playas (those with of the surface elements (e.g. sand ripples, plants,
shallow groundwater tables and hence signi¬- rocks, etc.). Particles can be entrained by the wind
cant surface moisture) can be expected to emit when the shear velocity is greater than a thresh-
less dust than dry playas under similar wind con- old value, u—td , given by
ditions, because in wet playas the surface mois-
ρs ’ ρa
ture creates a cohesive force between grains that u—td = A gD (1.15)
is not present in dry playas. Franklin Lake Playa,
however, is one of the most active playas in the where A is a constant of proportionality (equal
western US despite its very shallow water table to 0.1 for air), ρs is the density of sediment, ρa
(less than 3 m below the surface in most por- is the density of air, g is the acceleration due to
tions of the playa). The shallow water table be- gravity, and D is grain size diameter.
neath Franklin Lake Playa causes a vapor dis- If the soil is wet, Eq. (1.15) must be modi¬ed
charge that disrupts the formation of surface to include the effects of cohesion. Chepil (1956)
crusts that would otherwise serve to protect the developed an empirical equation to represent the
surface from erosion. The result is a soft puffy effects of moisture on the threshold shear veloc-
surface on Franklin Lake Playa that promotes ity:
dust emission (Czarnecki, 1997). In general, dust
production on playas is relatively low for very wet 2
u—t = +
u2 (1.16)
or very dry playas, and relatively high for playas —td
ρa θ1.5
in an intermediate state.
where θ is the volumetric soil moisture and θ1.5
Temporal variations in dust emissions from
is the soil moisture at a pressure of ’1.5 MPa (i.e.
a single playa are also quite complex. In many
playas, dust emissions are observed to be nega- the wilting point).
tively correlated with antecedent precipitation. Quantifying the dust cycle is important for
In such cases, ¬‚ood events cause a pulse of geomorphic, pedologic, and hydrologic reasons.
recharge that raises the water table, wetting the Dust deposition, for example, controls the perme-
surface and increasing the threshold wind ve- ability of alluvial fan terraces and their rates of
locity. In other playas, however, dust emissions soil development. Windblown dust is also a very
are positively correlated with antecedent rainfall signi¬cant human health hazard. Windblown
(Reheis, 2006). In these playas, dust emissions are particulate matter includes both natural (e.g.
limited by the supply of ¬ne-grained material to mineral dust from playas) and anthropogenic
the playa. In such playas, signi¬cant dust storms (e.g. smoke from ¬res), but natural sources rep-
are produced only when antecedent ¬‚oods bring resent a signi¬cant portion of the total. Studies
suf¬cient sediment to the playa to be mobilized have shown a correlation between daily mortality
when the playa becomes dessicated. and high levels of particulate matter (PM) in the
Particles are entrained by the wind when the atmosphere (Samet et al., 2000). This means that,
Bernoulli lift created by wind ¬‚owing over the in addition to causing long-term respiratory ail-
top of the particle exceeds the weight of the par- ments, high dust concentrations also cause sud-
ticle. The Bernoulli force, in turn, is directly re- den death in some individuals. Figure 1.17 is a
lated to the shear velocity exerted by the ¬‚ow, photograph of a dust storm that hit Las Vegas,
given by (Bagnold, 1941) Nevada on April 15, 2002. Although such dust
κum storms have long been a hazard, rapid popula-
u— = (1.14)
tion growth in the western US may be contribut-
ln z0
ing to more dust storms through groundwater
where κ is the von Karman constant (equal to withdrawal and land disturbance.
0.4), um is the wind velocity measured at a height Atmospheric transport of particulate matter
zm above the ground, and z0 is the aerodynamic can be modeled as a combination of downwind
roughness length of the surface. The roughness advection, turbulent diffusion, and gravitational
length z0 , in turn, is related to the size and shape settling. The steady-state equation describing

whether or not the particles are falling under
gravity. Silt, for example, has a negligible settling
velocity but a ¬nite depositional velocity because
silt particles are deposited as they fall or lodge
between the crags of a rough surface (e.g. clasts
in a desert pavement).
Dust deposition plays an important role in
creating desert pavement, one of the most enig-
matic landforms of the arid environment. Desert
pavement is a type of surface landform formed
on alluvial fan terraces in low-elevation parts of
the western US and other arid regions around
the world. Desert pavement is characterized by
a monolayer of stones at the surface, each stone
Fig 1.17 Photograph of a dust storm in Las Vegas, Nevada,
tightly sutured to the next like pieces in a jig-
on April 15, 2002, when particulate matter (PM)
concentrations reached above 1400 μg/m3 . saw puzzle. Digging into the pavement reveals
that, in most places, the desert pavement sits
atop a silt-rich layer comprised predominantly
these processes is given by
of wind-blown material (Figure 1.19) (McFadden
‚ 2c ‚ 2c ‚ 2c ‚c ‚x et al., 1987; Wells et al., 1995). While desert pave-
+ 2+ 2 ’u +q =0
‚ x2 ‚y ‚z ‚x ‚z ment was originally thought to represent a lag of
(1.17) coarse material left behind as ¬ner material was
eroded out, the observation of a layer of ¬ne ma-
where K is the turbulent diffusivity, c is the
terial underneath the stony monolayer is most
particle concentration, x is the downwind dis-
consistent with a model for pavement formation
tance, y is the crosswind distance, z is the vertical
in which windblown silt and clay is deposited
distance from the ground, u is the mean wind
in the interstitial spaces between surface stones.
velocity, and q is the settling velocity (model
Once an initial layer of silt and clay is deposited
geometry shown in Figure 1.18a). Solutions to
on and immediately underneath the surface, the
Eq. (1.17) are known as Gaussian plumes. This ver-
reorganization of stony clasts can occur during
sion of the advection-diffusion-settling equation
wetting and drying events, freeze--thaw cycles,
assumes that K and u are uniform. More com-
and bioturbation (all of which require, or are en-
plex models are also available in which K and
hanced by, ¬ne-grained sediments). The reorga-
u vary with height to better represent transport
nization of clasts, in turn, helps protect the un-
processes close to the ground (e.g. Huang, 1999).
derlying ¬ne-grained layer from erosion. In this
Deposition from a Gaussian plume is modeled
way, the desert pavement and underlying eolian
by treating the ground as a sink for particles.
deposit coevolve. This coevolution is thought to
Deposition in this model is characterized by a
take several thousand years to initiate, but may
deposition velocity, p, de¬ned as the fraction of
occur faster in areas of higher dust deposition
the particle concentration just above the ground
rates (Pelletier et al., 2008a).
that undergoes deposition per unit time. In this
The texture of the alluvial-fan parent mate-
model framework, deposition at the ground sur-
rial plays an important role in determining the
face is equal to the downward ¬‚ux due to turbu-
speci¬c mode of pavement formation. In gravel-
lent diffusion and particle settling. This balance
rich parent material (i.e. those found near the
provides a ¬‚ux boundary condition at z = 0:
proximal part of an alluvial fan close to the
+ qc(x, y, 0) = pc(x, y, 0) mountain front), clast motion on the surface
K (1.18)
is unlikely to be signi¬cant until the underly-

Physically, the deposition velocity represents the ing eolian layer is suf¬ciently thick to cause
trapping ability of the surface, independent of expansion/contraction of the near-surface layer


z depositional

wind x

2 km
0.1 0.2 0.4 1.0

wind direction

numerical model results
u = 5 m/s, K = 5 m2/s, p = 0.05 m/s

wind direction
through wetting--drying and freeze--thaw cycles.
Fig 1.18 (a) Schematic diagram of model geometry.
In cases of sand-dominated parent material (e.g.
Depositional topography shown in this example is an inclined
the distal portions of many alluvial fans), the ini-
plane located downwind of source (but model can accept any
tial surface has relatively few clasts with which
downwind topography). (b) Grayscale maps of 3D model
to form a pavement. In these cases, pavement
results illustrating the role of variable downwind topography.
In each case, model parameters are u = 5 m/s, K = 5 m2 /s, clasts must ¬rst be pushed to the surface by
and p = 0.05 m/s. Downwind topography is, from left to right, freeze--thaw or wetting--drying cycles, or possibly
a ¬‚at plane, an inclined plane, and a triangular ridge. Modi¬ed made available through progressive fracturing of
from Pelletier and Cook (2005).
larger clasts (McFadden et al., 2005). As the surface
gains clast-material coverage, lateral migration
serves to interlock and suture the clasts as in the


eolian deposits

sand poor

parent material

Fig 1.19 Layered stratigraphy of typical alluvial fan terrace
(Eagle Mountain piedmont), illustrating stony monolayer
comprising the pavement, underlain by silt-rich eolian
deposits sitting atop parent alluvium.

gravel-dominated parent material sand-dominated parent material
bar bar


clast migration
Fig 1.21 (a) Virtual oblique aerial photograph of the margin
of Roach Playa and the adjacent Whiskey Pete™s alluvial fan
along the southern Nevada/California border. (b) Schematic
diagram showing how signi¬cant along-strike relief prevents
continuous saltation in the proximal portion of the fan,
resulting in local storage and accumulation of sand. On the
distal portion of the fan, conversely, continuous saltation
Fig 1.20 Schematic illustration of two modes of enables an eolian corridor to form.
parent-material formation associated with gravel-dominated
and sand-dominated parent material.
Playa is distinctly different from that of silt-
dominated playas such as Franklin Lake Playa. In
gravel-dominated case. These two distinct modes particular, the topographic relief of the alluvial
of pavement formation are illustrated schemati- fan plays a signi¬cant role in controlling the spa-
cally in Figure 1.20. tial pattern of deposition.

. 1
( 14)