Transcript Document

Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
MORPHODYNAMICS OF SAND-BED RIVERS ENDING IN DELTAS
Wax Lake Delta,
Louisiana
Delta from Iron
Ore Mine into Lake
Wabush, Labrador
Delta of Eau Claire River at
Lake Altoona, Wisconsin
1
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
DELTA SHAPE
Most deltas show a 2D spread in the lateral
direction (fan-delta). The channel(s)
migrate and avulse to deposit over the
entire surface of the fan-delta as they
prograde into standing water.
But some deltas forming in
canyons can be approximated as
1D. The delta of the Colorado
River at Lake Mead was confined
within a canyon until recently.
2D progradation
1D progradation
Mangoky River, Malagasy
2
Colorado River at Lake Mead, USA
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
AN EXAMPLE OF A 1D DELTA
View of the Colorado River at
the upstream end of Lake Mead.
Image from NASA
https://zulu.ssc.nasa.gov/mrsid/mrsid.pl
Hoover Dam was closed in
1936. Backwater from the
dam created Lake Mead.
Initially backwater extended
well into the Grand Canyon.
For much of the history of
Lake Mead, the delta at the
upstream end has been so
confined by the canyon that is
has propagated downstream
as a 1D delta. As is seen in
the image, the delta is now
spreading laterally into Lake
Mead, forming a 2D fan-delta.
3
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
HISTORY OF SEDIMENTATION IN LAKE MEAD, 1936 - 1948
Image based on an original from
Smith et al. (1960)
4
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
WHERE DOES THE SAND GO? WHERE DOES THE MUD GO?
Image based on an original from
Smith et al. (1960)
5
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
STRUCTURE OF A DELTA: TOPSET, FORESET AND BOTTOMSET
A typical delta deposit can be divided into a topset, foreset and bottomset. The
topset is coarse-grained (sand or sand and gravel), and is emplaced by fluvial
deposition. The foreset is also coarse-grained, and is emplaced by avalanching.
The bottomset is fine-grained (mud, e.g. silt and clay) and is emplaced by either
plunging turbidity currents are rain from sediment-laden surface plumes.
topset
standing water
gravel/sand
foreset
bottomset
antecedent bed profile
mud
(silt/clay)
6
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
7
(Kostic and Parker, 2003a,b)
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
SIMPLIFICATION: TOPSET AND FORESET ONLY
Here the problem is simplified by considering a topset and foreset only. That is, the
effect of the mud is ignored. It is not difficult to include mud: see Kostic and Parker
(2003a,b).
topset
gravel/sand
standing water
foreset
antecedent bed profile
8
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
DERIVATION OF THE 1D EXNER EQUATION OF SEDIMENT
CONSERVATION
The channel has a constant width. Let x = streamwise distance, t = time, qt = the
volume sediment transport rate per unit width and p = bed porosity (fraction of
bed volume that is pores rather than sediment). The mass sediment transport
rate per unit width is then sqt, where s is the material density of sediment.
Mass conservation within a control volume with length x and a unit width
(Exner, 1920, 1925) requires that:
/t (sediment mass in bed of
control volume) = mass
sediment inflow rate – mass
sediment outflow rate
H
x
or



s (1   p )x  1  s qt x  qt
t
or
qt

(1   p )
t
x
x  x
qt
qt
 1

1
control
volume
9
x
x+x
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
BUT HOW ARE FLOODS ACCOUNTED FOR IN THE EXNER EQUATION OF
SEDIMENT CONTINUITY?
Daily Discharge, Mississsippi River at Tarbert Landing
45000
Fly River, PNG in flood
40000
35000
Q (m3/s)
30000
25000
20000
15000
10000
5000
0
01-Oct- 27-Jun- 24-Mar- 18-Dec- 13-Sep- 10-Jun- 06-Mar- 30-Nov- 27-Aug- 23-May-
10
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
SIMPLE ADAPTATION TO ACCOUNT FOR FLOODING
Rivers move the great majority of their sediment, and are morphodynamically
active during floods. Paola et al. (1992) represent this in terms of a flood
intermittency If. They characterize floods in terms of bankfull flow, which carries
a volume bed material transport rate per unit width qt for whatever fraction of
time is necessary is necessary to carry the mean annual load of the river. That
is, where Gma is the mean annual bed material load of the river, Bbf is bankfull
width and Ta is the time of one year, If is adjusted so that
Gma  sIfBbf qtTa  (1 R)IfBbf qtTa
where  = water density. The river is assumed to be morphodynamically inactive
at other times. Wright and Parker (2005a,b) offer a specific methodology to
estimate If.
The Exner equation is thus modified to

qt
(1  p )
 If
t
x
transport
rate at
bankfull
fraction
of
timeflow
flow 11is
in flood
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
KEY FEATURES OF THE MORPHODYNAMICS OF THE SELF-EVOLUTION OF 1D
SAND-BED DELTAS UNDER THE INFLUENCE WITH BACKWATER
Water discharge per unit width qw is conserved, and is given by the relation
qw  UH
and shear stress is related to flow velocity using a Chezy (recall Cf-1/2 = C) or
Manning-Strickler formulation;
1/ 6
b  Cf U2
, Cf  const (Chezy)
or
Cf 1/ 2
H
 r  
 kc 
(Manning  Strickler )
In low-slope sand-bed streams boundary shear stress cannot be computed from from
the depth-slope product, but instead must be obtained from the full backwater
equation;
or thus
b  gHS

q2w
S  Cf
H
gH3

q2w
x
1
gH3
,
U
U
H
 b
U
 g
g 
t
x
x
x H
q2w
Cf q2w

b  CfU  Cf 2 ,  
H
RgDH 2
2
12
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
DOWNSTREAM BOUNDARY CONDITION FOR THE BACKWATER
FORMULATION
Base level is specified in terms of a prescribed downstream water surface
elevation d(t) = (L, t) +H(L, t) rather than downstream bed elevation (L, t).
The base level of the
Athabasca River, Canada
is controlled by the water
surface elevation of Lake
Athabasca.
Delta of the Athabasca River
at Lake Athabasca, Canada.
Image from NASA
https://zulu.ssc.nasa.gov/mrsid/mrsid.pl
13
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
THE MORPHODYNAMIC PROBLEM
The formulation below assumes a single characteristic bed grain size D, a
constant Chezy resistance coefficient Cz = Cf-1/2 and the Engelund-Hansen (1967)
total bed material load formulation as an example.
(1   p )
q

 -If t
t
x
qt  RgD D
Exner equation
0.05  Cf q 


2 
Cf  RgDH 
2
w
5/2
Bed material load equation

q2w

 Cf
H
x
gH3

q2w
x
1
gH3
Backwater equation
( x, t ) t 0  I ( x)
Specified initial bed profile
qt ( x, t ) x0  qtf (t )
Specified upstream bed material feed rate
(x, t) xL  H(x, t) xL  d (t)
Downstream boundary condition, or base level set
14
in terms of specified water surface elevation.
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
NUMERICAL SOLUTION TO THE BACKWATER FORMULATION OF
MORPHODYNAMICS
The case of subcritical flow is considered here.
At any given time t the bed profile (x, t) is known. Solve the backwater equation
H  
q2w  
q2w 
 1  3 
  
 Cf
3 
gH   gH 
x  x
upstream from x = L over this bed subject to the boundary condition
H(L, t)  d (t)  (L, t )
Evaluate the Shields number and the bed material transport rate from the relations
2
C
q
f w
 
RgDH 2
Find the new bed at time t + t
 t  t
Repeat using the bed at t + t
qt  RgD D
 
0.05 

Cr
1 qt
 t 
t
(1  p ) x t
5/ 2
15
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
IN MORPHODYNAMICS THE FLOW
AND THE BED TALK TO AND
INTERACT WITH EACH OTHER
Honey, could you
scratch my back,
it itches in a place
I can’t reach.
qt

(1   p )
  If
t
x
dH S  S f (H)

dx 1  Fr 2 (H)
Sure, sweetie, but could
you cut my toenails for me
afterward? I can’t reach
‘em very well either.
16
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
THE PROBLEM OF IMPULSIVELY RAISED WATER SURFACE ELEVATION
(BASE LEVEL) AT t = 0
water surface elevation (base
level) is raised at t = 0 by e.g.
installation of a dam
sediment supply
remains constant
at qqsa
ta
M1 backwater curve

Note: the M1 backwater curve was
introduced in the lecture on hydraulics
and sediment transport
antecedent equilibrium bed
profile established with load qta
sa
before raising base level
17
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
RESPONSE TO IMPULSIVELY RAISED WATER SURFACE ELEVATION:
A PROGRADING DELTA THAT FILLS THE SPACE CREATED BY BACKWATER
Ultimate water surface
Initial water surface
Ultimate bed
Initial bed

transient bed profile
(prograding delta)
See Hotchkiss and Parker (1991) for more details.
18
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
FLOW OVER ANTECEDENT BED
Before the water surface is raised, the is assume to be at normal, mobile-bed
equilibrium with antecedent bed slope Sa and water dischage per unit width during
floods qw. It is useful to compute the characteristics of the normal flow that would
prevail with the specified flow over the specified bed. Let Hna and qtna denote the
flow depth and total volume bed material load per unit width that prevail at
antecedent normal mobile-bed equilibrium with flood water discharge per unit width
qw and bed slope Sa. From Slides 19 and 21 of the lecture on hydraulics and
sediment transport, these are given as
Hna
 Cf q
Hna  
 gSI
2
w
qtna
qw
1/ 3



 Sa
qta
Initial bed
0.05  Cf q
 RgD D

Cf  g
2
w
S 


 RD 
2/3
a
5/2
19
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
CHANGES IMPOSED AT t = 0
The reach has length L and the bed has antecedent bed slope Sa. It is assumed for
simplicity that the antecedent bed elevation at the downsteam end of the reach (x =
L) is da = 0, so that antecedent water surface elevation da is given as
da  Hna
At time t = 0 the water surface elevation (base level) at the downstream end d is
impulsively raised to a value higher than da and maintained indefinitely at the new
level. In addition, the sediment feed rate from the upstream end is changed from the
antecedent value qta to a new feed value qtf.
qtf
Hna
qw
 Sa
Initial bed
qta
d
20
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
ULTIMATE EQUILIBRIUM
Eventually the backwater zone (e.g. reservoir zone behind a dam) fills with sediment
and the river established a new normal, mobile-bed equilibrium in consonance with
the flood, water discharge per unit width qw, the sediment supply rate qtf (which
becomes equal to the transport rate qtu at ultimate equilibrium) and the water surface
elevation d. The bed slope Su and flow depth Hnu at this ultimate normal equilibrium
can be determined by solving the two equations below for them.
 Cf q
Hnu  
 gSu
2
w
1/ 3



Hnu
Final bed

qtu = qtf
qw
Su
0.05  Cf q2w  Su2 / 3 

qtu  qtf  RgD D


Cf  g  RD 
5/ 2
A morphodynamic numerical model can then be used to describe the evolution 21
from the antecedent equilibrium to the ultimate equilibrium
d
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
NUMERICAL MODEL: INITIAL AND BOUNDARY CONDITIONS
The channel is assumed to have uniform grain size D and some constant ambient
slope Sa (before changing conditions at t = 0) which is in equilibrium with an
ambient transport rate qta. The reach of interest has length L. The antecedent
bed profile (which serves as the initial condition for the calculation) is then
(x, t) t0  da  Sa (L  x)
where here da can be set equal to zero. The boundary condition at the upstream
end is the changed feed rate qtf for t > 0, i.e.
qt ( x, t ) x0  qtf (t )
where qtf(t) is a specified function (but here taken as a constant). The downstream
boundary condition, however, differs from that used in the normal flow calculation,
and takes the form
(x, t)  H(x, t) xL  d(t)
where d(t) is in general a specified function, but is here taken to be a constant.
Note that downstream bed elevation (L,t) is not specified, and is free to vary22
during morphodynamic evolution.
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
NUMERICAL MODEL: DISCRETIZATION AND BACKWATER CURVE
The reach L is discretized into M intervals of length x bounded by M+1 nodes. In
addition, sediment is fed in at a “ghost” node where bed elevation is not tracked.
L
Feed sediment here!
x 
xi  (i  1)x , i  1..M  1
M
x
i=1
ghost
2
3
M -1
M
i = M+1
L
The backwater calculation over a given bed proceeds as in the lecture on
hydraulics and sediment transport:
HM1  d (t )  M1
Hp  Hi1 
1
Fback (Hi1 )x
2
, Hi  Hi1 
1

Fback (Hi1 )  Fback (Hp )x
2

i  i1
q2w  
q2w 



Fback (H)   S  Cf
1

,
S

gH3   gH3 
x

where i proceeds downward from M to 1.
23
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
NUMERICAL MODEL: SEDIMENT TRANSPORT AND EXNER
 
q2w
,  
, i  1..M  1
2
RgDH
1 qt,i
 i t 
If t , i  1..M  1
1 p x
0.05 
qt,i  RgD D
i
Cf
i t  t
5/ 2

i
qt,i1  qt,i
 qt,i  qt,i1
a

(
1

a
)
, i  1..M
u
qt,i  u x
x

qt,i  qt,i1
x 
, i  M1

x
Note that the difference scheme used to compute qt,i/x is a central
difference scheme only for au= 0.5. With a backwater formulation takes an
advectional-diffusional form and a value of au greater than 0.5 (upwinding)
is necessary for numerical stability. The difference qt,1 is computed using
the sediment feed rate at the ghost node:
qt,g  qt,0  qtf (t)
24
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
INTRODUCTION TO RTe-bookAgDegBWChezy.xls
The worksheet RTe-bookAgDegBWChezy.xls provides both a grapical user
interface and code in Visual Basic for Applications (VBA). A tutorial on VBA is
provided in the workbook Rte-bookIntroVBA; it introduces the concept of modules.
The code for the morphodynamic model is contained in Module 1 of RtebookAgDegBWChezy.xls. It can be seen by clicking “Tools”, “Macros”, “Visual
Basic Editor” from Excel, and then double-clicking “Module1” in the VBA Project
Window at the upper left of the Screen. The Security Level (“Tools”, “Macro”,
“Security”) must be set to no higher than “medium” in order to run the code.
Most of the input is specified in worksheet “Calculator”. The first set of input
includes: water discharge per unit width qw at flood, flood intermittency If, grain
size D, reach length L, Chezy resistance coefficient Cz, antecedent bed slope Sa
and volume total bed material feed rate per unit width during floods qtf. The
specified numbers are then used to compute the normal flow depth Hna at
antecedent conditions, the final ultimate equilibrium bed slope Su and the final
ultimate equilibrium normal flow depth Hnu.
The user then specifies a downstream water surface elevation d. This
value should be > the larger of either Hna or Hnu to cause delta formation.
25
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
INTRODUCTION TO RTe-bookAgDegBWChezy.xls contd.
The following input parameters are then specified on worksheet “Calculator” by the
user: reach length L, time step t, the number of time steps until data is generated
for output (by printing it onto another worksheet in the workbook) Ntoprint, the
number of times data is generated Nprint, number of spatial intervals M and
upwinding parameter au. The total duration of the calculation is thus equal to t x
Ntoprint x Nprint, and the spatial step length x = equal to L/M.
The parameter R is specified in worksheet “AuxiliaryParameter”.
Once all the input parameters are specified, the code is executed by clicking the
button “Do the Calculation” in worksheet “Calculator”.
The numerical output is printed onto worksheet “ResultsofCalc”. The output
consists of the position x, bed elevation , water surface elevation  and flow
depth H at every node for time t = 0 and Nprint subsequent times. The bed
elevations and final water surface elevations are plotted on worksheet
“PlottheData”.
26
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
INTRODUCTION TO RTe-bookAgDegBWChezy.xls contd.
In worksheet “Calculator” the flow discharge qw (m2/s) and bed material feed rate
at flood flow qtf (m2/s) are specified per unit channel width.
In worksheet “MeanAnnualFeedRate” the user can specify a channel width Bf at
flood flow (e.g. bankfull width). The flood discharge Qf = qf Bf in m3/s and the
mean annual bed material feed rate Gma are then computed directly on the
worksheet.
The input for all the cases (Cases A ~ G) illustrated subsequently in this
presentation is given in worksheet “WorkedCases”.
As noted in Slide 25, the code itself can be viewed by clicking “Tools”, “Macros”
and “Visual Basic Editor”, and then double-clicking Module 1 in the VBA Project
Window in the upper left of the screen. Each unit of the code is termed a “Sub” or
a “Function” in VBA. Three of these units are illustrated in Slides 29, 30 and 31.
27
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
GRAPHICAL USER INTERFACE
The graphical user interface in worksheet “Calculator” is shown below.
Calculation of Morphodynamic Response of River to Imposed Backwater
(qw)
(Inter)
(D)
(lamp)
(Cz)
(SI)
(qtf)
(xid)
(Dt)
au
Calculation of ambient river conditions (before imposed change)
Assumed parameters
2
qw
10.00 m /s
Flood discharge
The orange boxes:
If
0.20
Intermittency
indicate the parameters you must specify.
D
0.25 mm
Grain Size
The numbers in the blue boxes:
p
0.40
Bed Porosity
are computed for you
Cz
14.00
Chezy resistance coefficient
Sa
Antecedent bed slope (initial downstream bed elevation LI = 0)
0.0000739
2
qtf
0.000500 m /s
Bed material feed rate during floods
2
qtna
0.000500 m /s
Antecedent bed material transport at normal equilibrium
Hna
Depth of normal flow over initial bed slope SI with flow discharge qw
8.89 m
Su
Ultimate bed slope in equilibrium with sediment feed rate qtf and flow discharge qw
0.0000739
Hnu
Ultimate depth of normal flow in equilibrium with sediment feed rate q tf and flow dicharge qw
8.89 m
d
L
t
Ntoprint
Nprint
M
au
8.89 m
100000
0.1
1000
6
50
1
Imposed downstream water surface elevation
Do not make d less than the lower of HnI or Hnu
m
length of reach
Click the button to perform a calculation
year
time step
Number of time steps to printout
Do the Calculation
Number of printouts
28
Intervals
Here 1 = full upwind scheme, 0.5 = central difference scheme
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
BACKWATER CALCULATION: Sub Do_Fluvial_Backwater
Sub Do_Fluvial_Backwater()
Dim Hpred As Double: Dim fr2p As Double: Dim fr2 As Double:
Dim fnp As Double: Dim fn As Double: Dim Cf As Double
Dim i As Integer
H(M + 1) = xio - eta(M + 1)
For i = 1 To M
fr2p = qw ^ 2 / g / H(M + 2 - i) ^ 3
Cf = (1 / alr ^ 2) * (H(M + 2 - i) / kc) ^ (-1 / 3)
fnp = (eta(M + 1 - i) - eta(M + 2 - i) - Cf * fr2p * dx) / (1 - fr2p)
Hpred = H(M + 2 - i) - fnp
fr2 = qw ^ 2 / g / Hpred ^ 3
fn = (eta(M + 1 - i) - eta(M + 2 - i) - Cf * fr2 * dx) / (1 - fr2)
H(M + 1 - i) = H(M + 2 - i) - 0.5 * (fnp + fn)
Next i
For i = 1 To M
xi(i) = eta(i) + H(i)
Next i
End Sub
29
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
LOAD CALCULATION: Sub Find_Shields_Stress_and_Load
Sub Find_Shields_Stress_and_Load()
Dim i As Integer
Dim taux As Double: Dim qstarx As Double: Dim Cfx As Double
For i = 1 To M + 1
taux = Cfx * (qw / H(i)) ^ 2 / (Rr * g * D)
qstarx = alt/cFS *taux ^ 2.5
qt(i) = ((Rr * g * D) ^ 0.5) * D * qstarx
Next i
End Sub
30
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
IMPLEMENTATION OF EXNER: Sub Find_New_eta
Sub Find_New_eta()
Dim i As Integer
Dim qtback As Double: Dim qtit As Double: Dim qtfrnt As Double: Dim qtdif As Double
For i = 1 To M
If i = 1 Then
qtback = qtf
Else
qtback = qt(i - 1)
End If
qtit = qt(i)
qtfrnt = qt(i + 1)
qtdif = au * (qtback - qtit) + (1 - au) * (qtit - qtfrnt)
eta(i) = eta(i) + dt / (1 - lamp) / dx * qtdif * Inter
Next i
qtdif = qt(M) - qt(M + 1)
eta(M + 1) = eta(M + 1) + dt / (1 - lamp) / dx * qtdif * Inter
time = time + dt
End Sub
31
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
CASE A
This is a case for which a) base level d is unaltered from the
antecedent value, and b) the bed material feed rate qtf is equal to
the antecedent normal equilibrium value qtna. Condition a) is
ensured by setting d equal to the antecedent normal equilibrium
depth Hna, and condition b) is ensured by setting the bed material
feed rate qtf so that the ultimate equilibrium bed slope Su is equal
to the antecedent bed slope Sa. The result is the intended one:
nothing happens over the 600-year duration of the calculation.
Bed evolution (+ Water Surface at End of Run)
18
qw
If
D
p
Cz
Sa
qtf
qtna
Hna
Su
Hnu
d
16
10.00
0.20
0.25
0.40
14.00
0.0000739
0.000500
0.000500
8.89
0.0000739
8.89
2
m /s
mm
2
m /s
2
m /s
m
m
8.89 m
Elevation in m
14
bed 0 yr
bed 100 yr
bed 200 yr
bed 300 yr
bed 400 yr
bed 500 yr
bed 600 yr
ws 600 yr
12
10
8
6
4
2
L
t
Ntoprint
Nprint
M
au
100000
0.1
1000
6
50
1
0
32
-2
0
20000
40000
60000
Distance in m
80000
100000
m
year
Number of
Number of
Intervals
Here 1 = fu
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
CASE B
This case has the same input as Case A, except that the
downstream water surface elevation is raised from 8.89 m to 20
m. The duration of the calculation is 6 years. A delta starts to
form at the upstream end!
Bed evolution (+ Water Surface at End of Run)
25
Elevation in m
20
bed 0 yr
bed 1 yr
bed 2 yr
bed 3 yr
bed 4 yr
bed 5 yr
bed 6 yr
ws 6 yr
15
delta!
10
5
qw
If
D
p
Cz
Sa
qtf
qtna
Hna
Su
Hnu
d
L
t
Ntoprint
Nprint
M
au
10.00
0.20
0.25
0.40
14.00
0.0000739
0.000500
0.000500
8.89
0.0000739
8.89
20000
40000
60000
Distance in m
80000
mm
2
m /s
2
m /s
m
m
20.00 m
100000
0.1
50
6
50
1
0
0
2
m /s
100000
33
m
year
Number of
Number of
Intervals
Here 1 = fu
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
CASE C
Same conditions as Case B, but the time duration is 30 years.
Bed evolution (+ Water Surface at End of Run)
25
Elevation in m
20
bed 0 yr
bed 5 yr
bed 10 yr
bed 15 yr
bed 20 yr
bed 25 yr
bed 30 yr
ws 30 yr
15
10
5
0
0
20000
40000
60000
Distance in m
80000
100000
34
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
CASE D
Same conditions as Case B, but the time duration is 60 years.
Bed evolution (+ Water Surface at End of Run)
25
Elevation in m
20
bed 0 yr
bed 10 yr
bed 20 yr
bed 30 yr
bed 40 yr
bed 50 yr
bed 60 yr
ws 60 yr
15
10
5
0
0
20000
40000
60000
Distance in m
80000
100000
35
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
CASE E
Same conditions as Case B, but the time duration is 120 years.
Bed evolution (+ Water Surface at End of Run)
30
Elevation in m
25
bed 0 yr
bed 20 yr
bed 40 yr
bed 60 yr
bed 80 yr
bed 100 yr
bed 120 yr
ws 120 yr
20
15
10
5
0
0
20000
40000
60000
Distance in m
80000
100000
36
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
CASE F
Same conditions as Case B, but the time duration is 600 years.
Bed evolution (+ Water Surface at End of Run)
The “dam” is filled and overtopped, and ultimate
normal mobile-bed equilibrium is achieved.
30
Elevation in m
25
bed 0 yr
bed 60 yr
bed 120 yr
bed 180 yr
bed 240 yr
bed 300 yr
bed 360 yr
ws 360 yr
20
15
10
5
0
0
20000
40000
60000
Distance in m
80000
100000
37
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
CASE G
Same conditions as Case D, but the upstream feed rate qtf is tripled compared to
the antecedent normal value qtna, so ensuring that Su is greater than Sa
Bed evolution (+ Water Surface at End of Run)
35
This slope is steeper than this slope
Elevation in m
30
bed 0 yr
bed 10 yr
bed 20 yr
bed 30 yr
bed 40 yr
bed 50 yr
bed 60 yr
ws 60 yr
25
20
15
10
5
0
0
20000
40000
60000
80000
100000
Distance in m
Note that the delta front has prograded out much farther than
Case D because the sediment feed rate is three times higher.
qw
If
D
p
Cz
Sa
qtf
qtna
Hna
Su
Hnu
d
L
t
Ntoprint
Nprint
M
au
10.00
0.20
0.25
0.40
14.00
0.0000739
0.001500
0.000500
8.89
0.0001429
7.14
2
m /s
mm
2
m /s
2
m /s
m
m
20.00 m
100000
0.1
100
6
50
1
m
year
Number of tim
Number of pr
Intervals
Here 1 = full
38
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
GENERALIZATIONS
The basic formulation can be generalized to include
• a self-formed channel with varying width,
• channel sinuosity,
• a channel-floodplain complex in which mud as well as sand can deposit,
• a foreset of specified slope and
• a 2D geometry that yields a fan shape to the delta.
These generalizations are implemented for the Wax Lake Delta shown below in the
paper “Large scale river morphodynamics: application to the Mississippi Delta”
(Parker et al,. 2006) included on the CD for this course as file “WaxLake.pdf”.
39
Contribution from the National Center for Earth-surface Dynamics
for the Short Course
Environmental Fluid Mechanics: Theory, Experiments and Applications
Karlsruhe, Germany, June 12-23, 2006
REFERENCES
Exner, F. M., 1920, Zur Physik der Dunen, Sitzber. Akad. Wiss Wien, Part IIa, Bd. 129 (in
German).
Exner, F. M., 1925, Uber die Wechselwirkung zwischen Wasser und Geschiebe in Flussen,
Sitzber. Akad. Wiss Wien, Part IIa, Bd. 134 (in German).
Hotchkiss, R. H. and Parker, G., 1991, Shock fitting of aggradational profiles due to backwater,
Journal of Hydraulic Engineering, 117(9): 1129-1144.
Kostic, S. and Parker, G., 2003a, Progradational sand-mud deltas in lakes and reservoirs. Part
1. Theory and numerical modeling, Journal of Hydraulic Research, 41(2), pp. 127-140.
Kostic, S. and Parker, G.. 2003b, Progradational sand-mud deltas in lakes and reservoirs. Part
2. Experiment and numerical simulation, Journal of Hydraulic Research, 41(2), pp. 141-152.
Paola, C., Heller, P. L. & Angevine, C. L., 1992, The large-scale dynamics of grain-size variation
in alluvial basins. I: Theory, Basin Research, 4, 73-90.
Parker, G., Sequeiros, O. and River Morphodynamics Class of Spring 2006, 2006, Large scale
river morphodynamics: application to the Mississippi Delta, Proceedings, River Flow 2006
Conference, Lisbon, Portugal, September 6 – 8, 2006, Balkema.
Smith, W. O., Vetter, C.P. and Cummings, G. B., 1960, Comprehensive survey of Lake Mead,
1948-1949: Professional Paper 295, U.S. Geological Survey, 254 p.
Wright, S. and Parker, G., 2005a, Modeling downstream fining in sand-bed rivers. I: Formulation,
Journal of Hydraulic Research, 43(6), 612-619.
Wright, S. and Parker, G., 2005b, Modeling downstream fining in sand-bed rivers. II: Application,
40
Journal of Hydraulic Research, 43(6), 620-630.