of 27
1.
Introduction
River deltas are densely populated, ecologically diverse, and socioeconomically valuable landscapes
(Gleick,
2003
; Olson & Dinerstein,
1998
; Vörösmarty et al.,
2009
) that are sensitive to changes in climate
(Bianchi & Allison,
2009
; Giosan,
2014
; Syvitski,
2008
). Previous work has identified three primary in
-
fluences of climate on delta growth, namely 1) changes in sea level, 2) changes in sediment supply and
transport capacity, and 3) changes in flood regime (Blum & Törnqvist,
2000
; Hallet et al.,
1996
; Langbein &
Schumm,
1958
; Members,
1988
; Walling & Webb,
1996
). However, little is known about how these changes
might affect the occurrence of channel avulsions, i.e., catastrophic shifts in river course to the shoreline
(Ganti et al.,
2014
; Mohrig et al.,
2000
; Slingerland & Smith,
2004
). River avulsions are a fundamental pro
-
cess responsible for building delta lobes, and tend to occur at a characteristic location and frequency for
each delta (Figures
1a
and
1b
). Understanding where and when avulsions will occur in response to climate
change is crucial for predicting future flood hazards and sustaining land for coastal cities and ecosystems.
At the downstream boundary of deltas, climate-induced sea-level changes impact delta growth and retreat
(Stanley & Warne,
1994
). Field observations show sea-level rise enhances aggradation in the lower reaches
of deltaic rivers (Fisk,
1945
; Powell,
1875
; Schumm,
1993
), a trend supported by laboratory experiments
(Kim et al.,
2006
; Martin et al.,
2009
) and numerical models (Parker et al.,
2008
; Swenson et al.,
2005
). The
avulsion frequency
A
f
has been found to scale inversely with the channel-filling timescale
c
T
(Figure
1c
)
(Jerolmack & Mohrig,
2007
; Reitz & Jerolmack,
2012
),
1/
Ac
fT
(1)
where
/
cca
THv
is the time required to fill the channel of depth
c
H
at aggradation rate
a
v
(Figure
1d
).
Models also show a regime at low sea-level rise rates, where avulsion frequency is insensitive to sea-level
rise because the rate of topset aggradation that drives avulsion is primarily controlled by sediment supply
and delta progradation (Chadwick et al.,
2019
,
2020
; Ratliff et al.,
2018
). At very high rise rates, avulsion
frequency is expected to either reach an upper limit where nearly the entire sediment supply is deposited
Abstract
Coastal rivers that build deltas undergo repeated avulsion events—that is, abrupt changes
in river course—which we need to understand to predict land building and flood hazards in coastal
landscapes. Climate change can impact water discharge, flood frequency, sediment supply, and sea
level, all of which could impact avulsion location and frequency. Here we present results from quasi-2D
morphodynamic simulations of repeated delta-lobe construction and avulsion to explore how avulsion
location and frequency are affected by changes in relative sea level, sediment supply, and flood regime.
Model results indicate that relative sea-level rise drives more frequent avulsions that occur at a distance
from the shoreline set by backwater hydrodynamics. Reducing the sediment supply relative to transport
capacity has little impact on deltaic avulsions, because, despite incision in the upstream trunk channel,
deltas can still aggrade as a result of progradation. However, increasing the sediment supply relative to
transport capacity can shift avulsions upstream of the backwater zone because aggradation in the trunk
channel outpaces progradation-induced delta aggradation. Increasing frequency of overbank floods
causes less frequent avulsions because floods scour the riverbed within the backwater zone, slowing net
aggradation rates. Results provide a framework to assess upstream and downstream controls on avulsion
patterns over glacial-interglacial cycles, and the impact of land use and anthropogenic climate change on
deltas.
CHADWICK AND LAMB
© 2021. American Geophysical Union.
All Rights Reserved.
Climate-Change Controls on River Delta Avulsion
Location and Frequency
A. J. Chadwick
1,2
and M. P. Lamb
1
1
Division of Geological and Planetary Sciences, California Institute of Technology, Pasadena, CA, USA,
2
Department of
Earth and Environmental Sciences & St. Anthony Falls Laboratory, University of Minnesota, Minneapolis, MN, USA
Key Points:
Sea-level rise causes more frequent
river avulsions at a location relative
to the shoreline set by backwater
hydrodynamics
Increases in upstream sediment
supply can trigger avulsions
upstream of the backwater zone
Increases in flood frequency cause
erosion in the backwater zone,
reducing long-term aggradation rate
and avulsion frequency
Supporting Information:
Supporting Information may be found
in the online version of this article.
Correspondence to:
A. J. Chadwick,
achadwic@umn.edu
Citation:
Chadwick, A. J., & Lamb, M. P. (2021).
Climate-change controls on river
delta avulsion location and frequency.
Journal of Geophysical Research: Earth
Surface
,
126
, e2020JF005950.
https://
doi.org/10.1029/2020JF005950
Received 10 NOV 2020
Accepted 18 MAY 2021
10.1029/2020JF005950
RESEARCH ARTICLE
1 of 27
Journal of Geophysical Research: Earth Surface
on the active delta topset (Chadwick et al.,
2020
) or avulsions do not occur because the entire delta is
drowned (Muto,
2001
; Muto et al.,
2007
; Parker et al.,
2008
; Tomer et al.,
2011
). Enhanced aggradation
during sea-level rise has been linked to more frequent avulsions on the Rhine-Meuse delta (Stouthamer &
Berendsen,
2001
; T
ö
rnqvist,
1994
) and in delta laboratory experiments (Martin et al.,
2009
). However, the
Mitchell River delta provides a counterexample where avulsion frequency was reduced during Holocene
sea-level rise (T. I. Lane et al.,
2017
).
Relative sea-level rise can also affect avulsion locations, and thereby set the location of the delta apex. On
fan deltas, avulsions typically occur at a topographic slope break tied to a canyon outlet or bedrock-alluvial
transition (Ganti et al.,
2014
; Muto et al.,
2016
), which may be approximated as geographically fixed in
simplified models (Jerolmack,
2009
). Fan delta size—characterized by the distance
A
L
between the del
-
ta apex, or avulsion node, and the shoreline—approaches the autostratigraphic length-scale,
auto
L
(Muto
et al.,
2007
),
Aauto
LL
(2)
where
/
autos
Lq
,
s
q
is width-averaged sediment supply, and
is relative sea-level rise rate.
auto
L
repre
-
sents the delta size when rise rate and aggradation rate are in equilibrium, assuming the entire sediment
supply is deposited uniformly on the delta topset (i.e.,
 
vqL
asA
/
). This equilibrium is not sustainable
because part of the sediment supply is delivered to the delta foreset, which necessarily grows thicker as
sea-level rises. As a result, the delta shrinks and eventually drowns, in a process known as autoretreat over
timescales of
2
auto
auto
s
LS
T
q
, where
S
is channel-bed slope (Muto,
2001
; Muto et al.,
2007
; Parker et al.,
2008
;
Tomer et al.,
2011
).
CHADWICK AND LAMB
10.1029/2020JF005950
2 of 27
Figure 1.
Mapped delta lobes of (a) the Mississippi River and (b) the Yellow River overlain on Google Earth imagery, showing characteristic avulsion length,
A
L
, and avulsion frequency,
A
f
(after Coleman et al.,
1998
and Pang & Si,
1979
). (c) Scaling relationship between measured avulsion frequency and the inverse
of the channel-filling timescale (Equation
1
) (Ganti et al.,
2014
; Jerolmack & Mohrig,
2007
). (d) Schematic channel cross section showing aggradation at rate
a
v
. Channel avulsion occurs when the river has aggraded to a critical height, comparable to the channel depth
c
H
, that renders it unstable (Ganti et al.,
2014
;
Mohrig et al.,
2000
). (e) Scaling relationship between measured avulsion length and computed backwater length-scale (Equation
3
) (Chatanantavet et al.,
2012
;
Ganti, Chadwick, Hassenruck-Gudipati, Fuller, et al.,
2016
; Ganti et al.,
2014
).
Journal of Geophysical Research: Earth Surface
In contrast to fan deltas, lowland deltas commonly feature avulsion nodes on unconfined plains without
a topographic slope break (Brooke et al.,
2020
; Ganti et al.,
2014
). Avulsion nodes on lowland deltas have
been documented to shift with movement of the shoreline (Ganti, Chadwick, Hassenruck-Gudipati, Fuller,
& Lamb,
2016
; Ganti et al.,
2014
) to maintain a constant avulsion length
A
L
that scales with the backwater
length-scale
b
L
(Figure
1e
) (Chatanantavet et al.,
2012
; Jerolmack & Swenson,
2007
),
Ab
LL
(3)
where
/
bc
LHS
is the ratio of bankfull river channel depth,
c
H
, to channel-bed slope
S
(Lamb et al.,
2012
;
Paola & Mohrig,
1996
). Therefore, deltas with backwater-influenced avulsions might translate upstream,
rather than reduce in size, with rising sea level (Chadwick et al.,
2020
; Moran et al.,
2017
).
Changes in water and sediment supply cause aggradation and incision on deltas, which might also af
-
fect avulsion frequency and location (Blum & Törnqvist,
2000
; Schumm,
1993
). Early experimental efforts
showed higher sediment supply was associated with enhanced aggradation and greater avulsion frequency
(Ashworth et al.,
2004
; Bryant et al.,
1995
). Similarly, the Rhine-Meuse delta had more frequent avulsions
during a period of increased sediment supply in the Holocene (Stouthamer & Berendsen,
2001
). Avulsion
location may also be influenced by sediment supply, with pulses in sediment supply historically linked to
avulsion sites far upstream of deltas, for example, on the Tacquari megafan (Makaske et al.,
2012
) and in
New Zealand (Korup,
2004
).
Climate change can also affect the magnitude and frequency of large flood events on river deltas
(Knox,
2000
; Members,
1988
; Munoz et al.,
2018
). Flood regimes control the aggradation of both chan
-
nel and floodplain (Leopold & Maddock,
1953
; Naito & Parker,
2019
), with potential consequences for
river avulsion (Brizga & Finlayson,
1990
; Hajek & Edmonds,
2014
; Nicholas et al.,
2018
). For example,
on deltas with flashier flood regimes, field evidence indicates that less in-channel aggradation is neces
-
sary to trigger an avulsion, resulting in more frequent avulsions (Ganti et al.,
2014
; Moodie et al.,
2019
)
(Figure
1d
).
Flood regimes can also affect avulsion location. Chatanantavet et al. (
2012
) hypothesized that variable
flood discharges and non-uniform flow in the backwater zone (i.e., the reach within
b
L
of the shoreline)
cause a spatial maximum in long-term aggradation rate that determined the avulsion site. Scaled physical
experiments and subsequent modeling supported the backwater hypothesis, showing deltas produced a
backwater-scaled avulsion node when subjected to variable flood regimes with subcritical Froude num
-
bers (Chadwick et al.,
2019
; Ganti, Chadwick, Hassenruck-Gudipati, Fuller, et al.,
2016
; Ganti, Chadwick,
Hassenruck-Gudipati, & Lamb,
2016
). Numerical simulations have also isolated special cases where rapid
sea-level rise (Chadwick et al.,
2019
; Moran et al.,
2017
; Wu et al.,
2020
) or the assumption of a horizontal
deltaic floodplain (Chadwick et al.,
2019
; Ratliff et al.,
2021
) can give rise to backwater-scaled aggradation
and avulsion without flood variability.
Here, we build on recent advances in modeling backwater-driven avulsions on deltas (Chadwick
et al.,
2020
; Chatanantavet et al.,
2012
; Moodie et al.,
2019
; Moran et al.,
2017
; Ratliff et al.,
2018
) to
better understand the effect of climate change on avulsion frequency and location. We use the model of
Chadwick et al. (
2019
) that includes backwater hydrodynamics, variable flood discharges, and multiple
deltaic lobes that are abandoned and reoccupied due to avulsion. Previous work used the model to ad
-
dress the origin of a backwater-scaled avulsion node (Chadwick et al.,
2019
), and the effect of relative
sea-level rise on avulsion frequency (Chadwick et al.,
2020
). Here we present new results to explore the
effect of sea-level rise on avulsion location, and the roles of sediment supply and flood regime on avulsion
frequency and location. First, we briefly review the model and how climate forcing is parameterized. We
then present results from numerical simulations with systematic variation of climate forcing, exploring
systematic change in avulsion location and frequency as well as autogenic variability. Finally, we discuss
implications for avulsion dynamics over glacial-interglacial cycles and during modern anthropogenic
climate change.
CHADWICK AND LAMB
10.1029/2020JF005950
3 of 27
Journal of Geophysical Research: Earth Surface
2.
Model Framework
2.1.
Governing Equations
The model is the same as that presented in Chadwick et al. (
2019
) and is briefly reviewed here. The model
is designed to capture the dynamics of repeated delta-lobe construction and avulsion in a simplified, qua
-
si-2D framework. We consider a generic delta plain with an imposed number of lobes that are assumed to
form a branching pattern (Figure
2a
). Each lobe is modeled as a coupled river and floodplain of constant
width (
f
B
), channel sinuosity (
Ω
), wash load ratio
), and bed porosity
(
p
) using sediment mass balance
(Parker,
2004
; Parker et al.,
2008
),



 

1 ΛΩ
1
t
pf
Q
tx
B
(4)
CHADWICK AND LAMB
10.1029/2020JF005950
4 of 27
Figure 2.
(a) Conceptual model in plan view. Black solid lines are active channel of width
c
B
in a floodplain/lobe
of width
f
B
. Broken lines are abandoned channels. Shaded regions are deposits created during avulsion cycles 1–4.
At times of avulsion, the active lobe (in this case, lobe 4) is abandoned and the river is rerouted downstream of the
avulsion location (yellow star) to reoccupy the lowest-elevation abandoned lobe (Lobe 1), where it will begin building
a new lobe filling in the space indicated by (5). (b) Conceptual model in long profile, showing variable definitions and
climate-change boundary conditions of sea-level rise
, sediment supply
s
Q
, and flood regime
w
Q
(c) Same as (b), but
now showing riverbed aggradation and floodplain superelevation (
Δ
) of the active lobe (Lobe 4) relative to the lowest
abandoned lobe (Lobe 1). All delta lobes share a single trunk channel.
Journal of Geophysical Research: Earth Surface
where
is riverbed elevation,
t
is time,
x
is downstream distance, and
t
Q
is the flux of total bed-material
load. The river is prescribed an initial topset profile
0
ending at a foreset at
tf
xx
(Figure
2b
) and Equa
-
tion
4
is integrated over time using finite differences to determine topset aggradation. Foreset progradation
is approximated using a moving boundary formulation, following previous work (Hotchkiss & Parker,
1991
;
Kostic & Parker,
2003
; Swenson et al.,
2000
), such that all sediment delivered to the foreset contributes to its
progradation at a constant slope
a
S
(Text
S1
).
At a given time, the river occupies a single lobe where water and sediment are transported in a river channel
of prescribed depth (
)
c
H
and width
(
c
B
). Water is routed using a quasi-2D backwater equation for water mass
and momentum conservation under quasi-steady flow conditions (Chow,
1959
; Sturm & Tuzson,
2001
),

2
2
22
11
f
SCFr
dH
FrHdB
dx
Bdx
FrFr
(5)
where
H
is flow depth,
 
S
x
is riverbed slope,
f
C
is friction factor,

/
w
FrQBHgH
is Froude
number,
w
Q
is water discharge,
g
is gravity, and
B
is flow width. We assume a uniform flow width in the
channel, and an offshore plume with a constant spreading angle beyond the river mouth (Chatanantavet
et al.,
2012
; Lamb et al.,
2012
). The river mouth location (
m
x
) is set by the intersection of the floodplain pro
-
file (

 
fc
H
) and sea level (
sea
) (Figure
2b
; Text
S1
) (Chadwick et al.,
2019
). The river mouth typically
occurs slightly upstream of the topset-foreset break in our simulations, allowing for a dynamic water-sur
-
face elevation in the zone

m
tf
xxx
that is important for reproducing realistic backwater effects during
high flows (Text
S1
) (Chatanantavet et al.,
2012
; Chatanantavet & Lamb,
2014
). Sediment is routed accord
-
ing to the Engelund and Hansen (
1967
) relation for total bed-material load,

3
t
f
QBRgD
C
(6)
where
R
is submerged specific gravity of sediment,
D
is the median grain size of bed material,
2
( / )/
fw
C Q
BH
RgD
is Shields number,
0.05,
and
2.
Climate change affects model behavior through boundary conditions on Equations
4
5
. At the downstream
boundary, flow depth is determined by sea level,
H
xx
tf
seaxx
tf
||
 

(7)
where

 
seab
Ht
is sea level,
b
H
is initial basin depth offshore, and
is the user-specified relative
sea-level rise rate (Figure
2b
). Far upstream of the backwater zone, the river experiences normal flow con
-
ditions



0
dH
dx
and Equation
5
reduces to




1/3
2
2
fw
n
CQ
H
gSB
(8)
where
n
H
is the normal flow depth. Thus, normal flow depth varies according to the flood regime, which
is prescribed using a probability distribution
P
of monthly mean discharges
w
Q
at the upstream end (Fig
-
ure
2b
) (LeBoutillier & Waylen,
1993
; Stedinger,
1993
). Sediment supply enters the model domain via a
ghost node at the upstream end (Kostic & Parker,
2003
) at user-specified time-averaged rate
s
Q
. The instan
-
taneous sediment supply
s
Q
at the ghost node is covaried with water discharge such that the river maintains
transport capacity (
st
QQ
in Equation
6
) at a constant slope regardless of flow regime, which is a good
approximation for lowland rivers downstream of the hydrographic boundary layer (Mackin,
1948
; Wong &
Parker,
2006
),
 
0
ss
QPdQ
(9)
CHADWICK AND LAMB
10.1029/2020JF005950
5 of 27
Journal of Geophysical Research: Earth Surface

 

 





3
0
/
(10)
where
0
S
is the constant slope at the ghost node, obtained by combining and numerically integrating Equa
-
tions
8
10
for a given flow regime (
,
w
PQ
) and sediment supply (
s
Q
).
Avulsions occur where and when the active lobe first aggrades to a critical height relative to the neighboring
lobes (Figure
2c
), which is referred to as superelevation (
Δ
),

Δ
c
xHH
(11)
where
H
is the avulsion threshold, a dimensionless number of order unity (Ganti, Chadwick, Hassen
-
ruck-Gudipati, & Lamb,
2016
; Ganti et al.,
2014
; Mohrig et al.,
2000
). We measure superelevation as the
height of the active lobe floodplain (
f
) relative to the floodplain of the lowest-elevation abandoned lobe
(
,
fabandoned
), evaluated at the same distance downstream from the trunk channel:








,,
,
for
Δ
for
ffabandoned
mabandoned
fsea
mabandoned
x
xxx
x
x
xx
(12)
where
,
m abandoned
x
is the streamwise coordinate of the abandoned-lobe shoreline (Figure
2c
). Seaward of the
abandoned lobe, superelevation is measured relative to sea level (
sea
).
After each avulsion, the river is routed to the lowest abandoned lobe by joining the bed profile of the active
channel upstream of the avulsion site with the bed profile of the new flow path downstream,







123
MIN
,
,
A
new
abandonedabandonedabandoned
A
x
xx
x
xxxxx
(13)
where
new
is the new riverbed profile after avulsion,
is the riverbed profile before avulsion,
A
x
is the
avulsion location, and

12
,,
abandoned
abandoned
and
3
abandoned
are the three abandoned-lobe long profiles. The
MIN
operator here selects the abandoned profile that has the minimum mean elevation downstream of
the avulsion node. Inactive lobe shapes are unchanged when abandoned (Galloway,
1975
) but are partial
-
ly drowned in cases due to relative sea-level rise. After establishing the new flow path, lobe construction
(Equations
4
10
) and avulsion setup (Equation
11
) begin anew.
2.2.
Dimensional Analysis
We non-dimensionalize the model to allow applicability of our results for a broad range of lowland deltas.
Equations
4
13
are scaled by a characteristic channel depth (
c
H
), channel width (
c
B
), backwater length
(
b
L
), and channel-filling timescale
T
H
QLB
c
c
sbf
p


/
1
1

(Chatanantavet & Lamb,
2014
; Reitz & Jerol
-
mack,
2012
) that are assumed constant (Text
S2
). Dimensional analysis reveals that model behavior is sensi
-
tive to six dimensionless input parameters linked to climate change. Delta response to sea-level rise depends
on the ratio of rates of sea-level rise and delta aggradation (Chadwick et al.,
2020
; Ganti et al.,
2019
; Liang
et al.,
2016
; Muto & Steel,
1997
), estimated by the normalized relative sea-level rise rate
HnT
cc
/
(14)
where

 
1 /2
nN
is a constant representing the average number of avulsions before a given lobe is
reoccupied (Ganti et al.,
2019
) and
N
is the imposed number of delta lobes (Figure
2a
). The denominator
of Equation
14
represents a first-order estimate of the maximum possible aggradation rate, i.e., the delta's
capacity to aggrade in the face of sea-level rise, assuming all sediment is deposited uniformly across the
CHADWICK AND LAMB
10.1029/2020JF005950
6 of 27