Channel Routing
Channel Routing
Channel Routing
Channel Routing
Simulate the movement of water through a channel
Used to predict the magnitudes, volumes, and
temporal patterns of the flow (often a flood wave)
as it translates down a channel.
2 types of routing : hydrologic and hydraulic.
both of these methods use some form of the
continuity equation.
Continuity equation
Hydrologic Routing
Hydraulic Routing
Momentum Equation
2
Continuity Equation
The change in storage (dS) equals the difference
between inflow (I) and outflow (O) or :
O - I =
dt
dS
For open channel flow, the continuity equation
is also often written as :
q =
x
Q
+
t
A
c
c
c
c
A = the cross-sectional area,
Q = channel flow, and
q = lateral inflow
Continuity equation
Hydrologic Routing
Hydraulic Routing
Momentum Equation
3
Hydrologic Routing
Methods combine the continuity equation with some
relationship between storage, outflow, and possibly
inflow.
These relationships are usually assumed, empirical,
or analytical in nature.
An of example of such a relationship might be a
stage-discharge relationship.
Continuity equation
Hydrologic Routing
Hydraulic Routing
Momentum Equation
4
Use of Manning Equation
Stage is also related to the outflow via a relationship
such as Manning's equation
2
1
3
2 49 . 1
f h
S AR
n
Q =
Continuity equation
Hydrologic Routing
Hydraulic Routing
Momentum Equation
5
Hydraulic Routing
Hydraulic routing methods combine the continuity
equation with some more physical relationship
describing the actual physics of the movement of the
water.
The momentum equation is the common relationship
employed.
In hydraulic routing analysis, it is intended that the
dynamics of the water or flood wave movement be
more accurately described
Continuity equation
Hydrologic Routing
Hydraulic Routing
Momentum Equation
6
Momentum Equation
Expressed by considering the external forces acting on a
control section of water as it moves down a channel
)
S
-
S
g( =
A
vg
+
2x
A) y (
A
g
+
x
v
V +
t
v
f o
c
c
c
c
c
Henderson (1966) expressed the momentum equation as :
t
v
g
1
-
x
v
g
v
-
x
y
-
S
=
S
o f
c
c
c
c
c
c
Continuity equation
Hydrologic Routing
Hydraulic Routing
Momentum Equation
7
Combinations of Equations
Simplified Versions :
t
v
g
1
-
x
v
g
v
-
x
y
-
S
=
S o f
c
c
c
c
c
c
x
v
g
v
-
x
y
-
S
=
S o f
c
c
c
c
x
y
-
S
=
S o f
c
c
S
f
= S
o
Unsteady -Nonuniform
Steady - Nonuniform
Diffusion or noninertial
Kinematic
Continuity equation
Hydrologic Routing
Hydraulic Routing
Momentum Equation
8
Routing Methods
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
9
Modified Puls
The modified puls routing method is probably most often
applied to reservoir routing
The method may also be applied to river routing for
certain channel situations.
The modified puls method is also referred to as the
storage-indication method.
The heart of the modified puls equation is found by
considering the finite difference form of the continuity
equation.
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
10
Modified Puls
t
S
-
S
=
2
O
+
O
(
-
2
I
+
I 1 2 2 1 2 1
A
Continuity Equation
O
+
t
S
2
=
O
-
t
S
2
+
I
+
I 2
2
1
1
2 1
A
|
.
|
\
|
A
Rewritten
The solution to the modified puls method is accomplished by
developing a graph (or table) of O -vs- [2S/t + O]. In order
to do this, a stage-discharge-storage relationship must be
known, assumed, or derived.
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
11
Modified Puls
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
12
Modified Puls Example
Given the following hydrograph and the 2S/At + O curve, find the
outflow hydrograph for the reservoir assuming it to be completely full at
the beginning of the storm.
The following hydrograph is given:
Hydrograph For Modified Puls Example
0
30
60
90
120
150
180
0 2 4 6 8 10
Time (hr)
D
i
s
c
h
a
r
g
e
(
c
f
s
)
13
Modified Puls Example
The following 2S/At + O curve is also given:
14
Modified Puls Example
A table may be created as follows:
Time I
n
I
n
+I
n+1
2S
n
/t - O
n
2S
n
/t + O
n+1
O
n+1
(hr) (cfs) (cfs) (cfs) (cfs) (cfs)
0
1
2
3
4
5
6
7
8
9
10
11
12
15
Modified Puls Example
Hydrograph For Modified Puls Example
0
30
60
90
120
150
180
0 2 4 6 8 10
Time (hr)
D
i
s
c
h
a
r
g
e
(
c
f
s
)
Next, using the hydrograph and interpolation, insert the Inflow
(discharge) values.
For example at 1 hour, the inflow is 30 cfs.
Time I
n
I
n
+I
n+1
2S
n
/t - O
n
2S
n
/t + O
n+1
O
n+1
(hr) (cfs) (cfs) (cfs) (cfs) (cfs)
0 0
1 30
2 60
3 90
4 120
5 150
6 180
7 135
8 90
9 45
10 0
11 0
12 0
16
Modified Puls Example
Time I
n
I
n
+I
n+1
2S
n
/t - O
n
2S
n
/t + O
n+1
O
n+1
(hr) (cfs) (cfs) (cfs) (cfs) (cfs)
0 0 30
1 30
2 60
3 90
4 120
5 150
6 180
7 135
8 90
9 45
10 0
11 0
12 0
The next step is to add the inflow to the inflow in the next time step.
For the first blank the inflow at 0 is added to the inflow at 1 hour to
obtain a value of 30.
17
Modified Puls Example
Time I
n
I
n
+I
n+1
2S
n
/t - O
n
2S
n
/t + O
n+1
O
n+1
(hr) (cfs) (cfs) (cfs) (cfs) (cfs)
0 0 30
1 30 90
2 60 150
3 90 210
4 120 270
5 150 330
6 180 315
7 135 225
8 90 135
9 45 45
10 0 0
11 0 0
12 0 0
This is then repeated for the rest of the values in the column.
18
Modified Puls Example
Time I
n
I
n
+I
n+1
2S
n
/t - O
n
2S
n
/t + O
n+1
O
n+1
(hr) (cfs) (cfs) (cfs) (cfs) (cfs)
0 0 30 0 0
1 30 90 30
2 60 150
3 90 210
4 120 270
5 150 330
6 180 315
7 135 225
8 90 135
9 45 45
10 0 0
11 0 0
12 0 0
O
+
t
S
2
=
O
-
t
S
2
+
I
+
I 2
2
1
1
2 1
A
|
.
|
\
|
A
The 2S
n
/At + O
n+1
column can then be calculated using the following
equation:
Note that 2S
n
/At - O
n
and O
n+1
are set to zero.
30 + 0 = 2S
n
/At + O
n+1
19
Modified Puls Example
Time I
n
I
n
+I
n+1
2S
n
/t - O
n
2S
n
/t + O
n+1
O
n+1
(hr) (cfs) (cfs) (cfs) (cfs) (cfs)
0 0 30 0 0
1 30 90 30 5
2 60 150
3 90 210
4 120 270
5 150 330
6 180 315
7 135 225
8 90 135
9 45 45
10 0 0
11 0 0
12 0 0
Then using the curve provided outflow can be determined.
In this case, since 2S
n
/At + O
n+1
= 30, outflow = 5 based on the graph
provided.
20
Modified Puls Example
Time I
n
I
n
+I
n+1
2S
n
/t - O
n
2S
n
/t + O
n+1
O
n+1
(hr) (cfs) (cfs) (cfs) (cfs) (cfs)
0 0 30 0 0
1 30 90 20 30 5
2 60 150
3 90 210
4 120 270
5 150 330
6 180 315
7 135 225
8 90 135
9 45 45
10 0 0
11 0 0
12 0 0
To obtain the final column, 2S
n
/At - O
n
, two times the outflow is
subtracted from 2S
n
/At + O
n+1
.
In this example 30 - 2*5 = 20
21
Modified Puls Example
Time I
n
I
n
+I
n+1
2S
n
/t - O
n
2S
n
/t + O
n+1
O
n+1
(hr) (cfs) (cfs) (cfs) (cfs) (cfs)
0 0 30 0 0
1 30 90 20 30 5
2 60 150 74 110 18
3 90 210
4 120 270
5 150 330
6 180 315
7 135 225
8 90 135
9 45 45
10 0 0
11 0 0
12 0 0
The same steps are repeated for the next line.
First 90 + 20 = 110.
From the graph, 110 equals an outflow value of 18.
Finally 110 - 2*18 = 74
22
Modified Puls Example
Time I
n
I
n
+I
n+1
2S
n
/t - O
n
2S
n
/t + O
n+1
O
n+1
(hr) (cfs) (cfs) (cfs) (cfs) (cfs)
0 0 30 0 0
1 30 90 20 30 5
2 60 150 74 110 18
3 90 210 160 224 32
4 120 270 284 370 43
5 150 330 450 554 52
6 180 315 664 780 58
7 135 225 853 979 63
8 90 135 948 1078 65
9 45 45 953 1085 65
10 0 0 870 998 64
11 0 0 746 870 62
12 0 0 630 746 58
This process can then be repeated for the rest of the columns.
Now a list of the outflow values have been calculated and the
problem is complete.
23
Kinematic Wave
Kinematic wave channel routing is probably the most
basic form of hydraulic routing.
This method combines the continuity equation with a
very simplified form of the St. Venant equations.
Kinematic wave routing assumes that the friction slope
is equal to the bed slope.
Additionally, the kinematic wave form of the
momentum equation assumes a simple stage-discharge
relationship.
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
24
Kinematic Wave Basic
Equations
q =
x
Q
+
t
A
L
c
c
c
c
Q = A
m
An explicit finite difference scheme in a space-
time grid domain is often used for the solution
of the kinematic wave procedure.
X
t
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
25
Working Equation
2
q + q
= q =
x
A
-
A
2
A
+
A
m +
t
A
-
A 1) - j (i, j) (i, 1) - j 1, - (i 1) - j (i, 1) - j 1, - (i 1) - j (i,
1) - (m
1) - j (i, j) (i,
(
A
(
A
o
X
t
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
26
Wave Speed TOO Fast?
When the average celerity, c, is greater than the
ratio x/t, a conservative form of these
equations is applied. In this conservative form,
the spatial and temporal derivatives are only
estimated at the previous time step and previous
location.
q =
t
A
-
A
+
x
Q - Q
1) - j 1, - (i j) 1, - i j) 1, - (i j) (i,
A A
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
27
Kinematic Wave Assumptions
The method does not explicitly allow for separation of the
main channel and the overbanks.
Strictly speaking, the kinematic method does not allow
for attenuation of a flood wave. Only translation is
accomplished. There is, however, a certain amount of
attenuation which results from the finite difference
approximation used to solve the governing equations.The
hydrostatic pressure distribution is assumed to be
applicable, thus neglecting any vertical accelerations.
No lateral, secondary circulations may be present, i.e. -
the channel is represented by a straight line.
Channel slopes should be 10% or less.
The channel is stable with no lateral migration,
degradation, and aggredation.
Flow resistance may be estimated via Manning's equation
or the Chezy equation.
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
28
Muskingum Method
Sp = K O
Sw = K(I - O)X
Prism Storage
Wedge Storage
Combined S = K[XI + (1-X)O]
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
29
Muskingum, cont...
O
2
= C
0
I
2
+ C
1
I
1
+ C
2
O
1
Substitute storage equation, S into the S in
the continuity equation yields :
S = K[XI + (1-X)O]
O - I =
dt
dS
t 0.5 + Kx - K
t 0.5 - Kx
- =
C0
A
A
t 0.5 + Kx - K
t 0.5 + Kx
=
C1
A
A
t 0.5 + Kx - K
t 0.5 - Kx - K
=
C2
A
A
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
30
Muskingum Notes :
The method assumes a single stage-discharge
relationship.
In other words, for any given discharge, Q, there can be
only one stage height.
This assumption may not be entirely valid for certain
flow situations.
For instance, the friction slope on the rising side of a
hydrograph for a given flow, Q, may be quite different
than for the recession side of the hydrograph for the
same given flow, Q.
This causes an effect known as hysteresis, which can
introduce errors into the storage assumptions of this
method.
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
31
Estimating K
K is estimated to be the travel time through the reach.
This may pose somewhat of a difficulty, as the travel
time will obviously change with flow.
The question may arise as to whether the travel time
should be estimated using the average flow, the peak
flow, or some other flow.
The travel time may be estimated using the kinematic
travel time or a travel time based on Manning's
equation.
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
32
Estimating X
The value of X must be between 0.0 and 0.5.
The parameter X may be thought of as a weighting coefficient
for inflow and outflow.
As inflow becomes less important, the value of X decreases.
The lower limit of X is 0.0 and this would be indicative of a
situation where inflow, I, has little or no effect on the storage.
A reservoir is an example of this situation and it should be
noted that attenuation would be the dominant process
compared to translation.
Values of X = 0.2 to 0.3 are the most common for natural
streams; however, values of 0.4 to 0.5 may be calibrated for
streams with little or no flood plains or storage effects.
A value of X = 0.5 would represent equal weighting between
inflow and outflow and would produce translation with little
or no attenuation.
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
33
More Notes - Muskingum
The Handbook of Hydrology (Maidment, 1992)
includes additional cautions or limitations in the
Muskingum method.
The method may produce negative flows in the initial
portion of the hydrograph.
Additionally, it is recommended that the method be
limited to moderate to slow rising hydrographs being
routed through mild to steep sloping channels.
The method is not applicable to steeply rising
hydrographs such as dam breaks.
Finally, this method also neglects variable backwater
effects such as downstream dams, constrictions, bridges,
and tidal influences.
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
34
Muskingum Example Problem
Time Inflow
C
0
I
2
C
1
I
1
C
2
O
1
Outflow
0 3 3
1 5
2 10
3 8
4 6
5 5
A portion of the inflow hydrograph to a reach of channel is given
below. If the travel time is K=1 unit and the weighting factor is
X=0.30, then find the outflow from the reach for the period shown
below:
35
Muskingum Example Problem
The first step is to determine the coefficients in this problem.
The calculations for each of the coefficients is given below:
t 0.5 + Kx - K
t 0.5 - Kx
- =
C
0
A
A
C
0
= - ((1*0.30) - (0.5*1)) / ((1-(1*0.30) + (0.5*1)) = 0.167
t 0.5 + Kx - K
t 0.5 + Kx
=
C
1
A
A
C
1
= ((1*0.30) + (0.5*1)) / ((1-(1*0.30) + (0.5*1)) = 0.667
36
Muskingum Example Problem
C
2
= (1- (1*0.30) - (0.5*1)) / ((1-(1*0.30) + (0.5*1)) = 0.167
t 0.5 + Kx - K
t 0.5 - Kx - K
=
C
2
A
A
Therefore the coefficients in this problem are:
C
0
= 0.167
C
1
= 0.667
C
2
= 0.167
37
Muskingum Example Problem
Time Inflow
C
0
I
2
C
1
I
1
C
2
O
1
Outflow
0 3 0.835 2.00 0.501 3
1 5
2 10
3 8
4 6
5 5
The three columns now can be calculated.
C
0
I
2
= 0.167 * 5 = 0.835
C
1
I
1
= 0.667 * 3 = 2.00
C
2
O
1
= 0.167 * 3 = 0.501
38
Muskingum Example Problem
Next the three columns are added to determine the outflow at time
equal 1 hour.
0.835 + 2.00 + 0.501 = 3.34
Time Inflow
C
0
I
2
C
1
I
1
C
2
O
1
Outflow
0 3 0.835 2.00 0.501 3
1 5 3.34
2 10
3 8
4 6
5 5
39
Muskingum Example Problem
This can be repeated until the table is complete and the outflow at
each time step is known.
Time Inflow
C
0
I
2
C
1
I
1
C
2
O
1
Outflow
0 3 0.835 2.00 0.501 3
1 5 1.67 3.34 0.557 3.34
2 10 1.34 6.67 0.93 5.57
3 8 1.00 5.34 1.49 8.94
4 6 0.835 4.00 1.31 7.83
5 5 3.34 1.03 6.14
40
Muskingum-Cunge
Muskingum-Cunge formulation is similar to the
Muskingum type formulation
The Muskingum-Cunge derivation begins with the
continuity equation and includes the diffusion form
of the momentum equation.
These equations are combined and linearized,
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
41
Muskingum-Cunge
working equation
where :
Q = discharge
t = time
x = distance along channel
qx = lateral inflow
c = wave celerity
m = hydraulic diffusivity
Lat
cq
x
Q
x
Q
t
Q
+
c
c
=
c
c
+
c
c
2
2
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
42
Muskingum-Cunge, cont...
Method attempts to account for diffusion by taking
into account channel and flow characteristics.
Hydraulic diffusivity is found to be :
O
BS
Q
2
=
The Wave celerity in the x-direction is :
dA
dQ
C =
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
43
Solution of Muskingum-Cunge
Solution of the Muskingum is accomplished by
discretizing the equations on an x-t plane.
Q
C
+ Q
C
+ Q
C
+ Q
C
=
1 + j
1 + n
Q
L
4
n
1 + j
3
1 + n
j
2
n
j
1
x) - 2(1 +
k
t
2x +
k
t
=
C1
A
A
x) - 2(1 +
k
t
2x -
k
t
=
C2
A
A
x) - 2(1 +
k
t
k
t
- x) - 2(1
=
C3
A
A
x) - 2(1 +
k
t
k
t
2
=
C4
A
|
.
|
\
| A
X
t
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
44
Calculation of K & X
c
x
= k
A
|
|
.
|
\
|
A
=
x c BS
Q
X
O
1
2
1
Estimation of K & X is more physically based
and should be able to reflect the changing
conditions - better.
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
45
Muskingum-Cunge - NOTES
Muskingum-Cunge formulation is actually considered an
approximate solution of the convective diffusion equation.
As such it may account for wave attenuation, but not for
reverse flow and backwater effects and not for fast rising
hydrographs.
Properly applied, the method is non-linear in that the flow
properties and routing coefficients are re-calculated at each
time and distance step
Often, an iterative 4-point scheme is used for the solution.
Care should be taken when choosing the computation
interval, as the computation interval may be longer than
the time it takes for the wave to travel the reach distance.
Internal computational times are used to account for the
possibility of this occurring.
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
46
Muskingum-Cunge - NOTES
Muskingum-Cunge may also be used distributed modeling
Currently Lag and K and DWOPER are available to
perform this operation.
Muskingum-Cunge offers a compromise between the
simplicity of Lag and K and complexity of DWOPER.
Muskingum-Cunge is physically-based but not as data
intensive as DWOPER.
The data inputs needed are:
Control parameters
Hydrologic: Inflow hydrographs
Physical system: channel geometry (cross-sections and
channel profile)
Data outputs: Method will sum and route discharge
hydrographs to overall basin outlet.
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
47
Muskingum-Cunge Example
The hydrograph at the upstream end of a
river is given in the following table. The
reach of interest is 18 km long. Using a
subreach length Ax of 6 km, determine the
hydrograph at the end of the reach using the
Muskingum-Cunge method. Assume c =
2m/s, B = 25.3 m, S
o
= 0.001m and no lateral
flow.
Time (hr) Flow (m
3
/s)
0 10
1 12
2 18
3 28.5
4 50
5 78
6 107
7 134.5
8 147
9 150
10 146
11 129
12 105
13 78
14 59
15 45
16 33
17 24
18 17
19 12
20 10
48
Muskingum-Cunge Example
First, K must be determined.
K is equal to :
c
x
K
A
=
seconds 3000
/ 2
/ 1000 6
=
=
s m
km m km
K
Ax = 6 km, while c = 2 m/s
49
Muskingum-Cunge Example
The next step is to determine x.
|
|
.
|
\
|
A
=
x c BS
Q
x
O
1
2
1
All the variables are known, with B =
25.3 m, S
o
= 0.001 and Ax =6000 m,
and the peak Q taken from the table.
253 . 0
/ ) 6000 )( 2 )( 001 . 0 )( 3 . 25 (
/ 150
1
2
1
3
3
=
|
|
.
|
\
|
=
s m
s m
x
Time (hr) Flow (m
3
/s)
0 10
1 12
2 18
3 28.5
4 50
5 78
6 107
7 134.5
8 147
9 150
10 146
11 129
12 105
13 78
14 59
15 45
16 33
17 24
18 17
19 12
20 10
50
Muskingum-Cunge Example
A curve for Ax/cAt is then needed to determine At.
For x = 0.253, Ax/(cAt) < 0.82
51
Muskingum-Cunge Example
Therefore, At can be found.
seconds 7200 Use
seconds 3658
82 . 0 / 2
/ 1000 6
) 82 . 0 )( /(
82 . 0 ) /(
= A
> A
> A
A > A
< A A
t
t
s m
km m km
t
c x t
t c x
52
Muskingum-Cunge Example
The coefficients of the Muskingum-Cunge method can now
be determined.
) 1 ( 2
2
1
x
K
t
x
K
t
C
+
A
+
A
=
7466 . 0
) 253 . 0 1 ( 2
3000
7200
) 253 . 0 ( 2
3000
7200
1
=
+
+
= C
53
Muskingum-Cunge Example
The coefficients of the Muskingum-Cunge method can now
be determined.
) 1 ( 2
2
2
x
K
t
x
K
t
C
+
A
A
=
4863 . 0
) 253 . 0 1 ( 2
3000
7200
) 253 . 0 ( 2
3000
7200
2
=
+
= C
54
Muskingum-Cunge Example
The coefficients of the Muskingum-Cunge method can now
be determined.
) 1 ( 2
) 1 ( 2
3
x
K
t
K
t
x
C
+
A
A
=
232 . 0
) 253 . 0 1 ( 2
3000
7200
3000
7200
) 253 . 0 1 ( 2
3
=
+
= C
55
Muskingum-Cunge Example
The coefficients of the Muskingum-Cunge method can now
be determined.
) 1 ( 2
2
4
x
K
t
K
t
C
+
A
|
.
|
\
|
A
=
233 . 1
) 253 . 0 1 ( 2
3000
7200
3000
7200
2
4
=
+
|
.
|
\
|
= C
56
Muskingum-Cunge Example
Then a simplification of the original formula can be made.
L
n
j
n
j
n
j
n
j
Q C Q C Q C Q C Q
4 1 3
1
2 1
1
1
+ + + =
+
+ +
+
Since there is not lateral flow, Q
L
= 0. The simplified
formula is the following:
n
j
n
j
n
j
n
j
Q C Q C Q C Q
1 3
1
2 1
1
1 +
+ +
+
+ + =
57
Muskingum-Cunge Example
A table can then be created in 2 hour time steps similar to the
one below:
Time (hr) 0 km 6 km 12 km 18 km
0 10
2 18
4 50
6 107
8 147
10 146
12 105
14 59
16 33
18 17
20 10
22 10
24 10
26 10
28 10
58
Muskingum-Cunge Example
It is assumed at time zero, the flow is 10 m
3
/s at each distance.
Time (hr) 0 km 6 km 12 km 18 km
0 10 10 10 10
2 18
4 50
6 107
8 147
10 146
12 105
14 59
16 33
18 17
20 10
22 10
24 10
26 10
28 10
59
Muskingum-Cunge Example
Next, zero is substituted into for each letter to solve the equation.
n
j
n
j
n
j
n
j
Q C Q C Q C Q
1 3
1
2 1
1
1 +
+ +
+
+ + =
0
1 3
1
0 2
0
0 1
1
1
Q C Q C Q C Q + + =
60
Muskingum-Cunge Example
Using the table, the variables can be determined.
0
1 3
1
0 2
0
0 1
1
1
Q C Q C Q C Q + + =
Time (hr) 0 km 6 km 12 km 18 km
0 10 10 10 10
2 18
4 50
6 107
8 147
10 146
12 105
14 59
16 33
18 17
20 10
22 10
24 10
26 10
28 10
=
=
=
0
1
1
0
0
0
Q
Q
Q
10
18
10
61
Muskingum-Cunge Example
Therefore, the equation can be solved.
0
1 3
1
0 2
0
0 1
1
1
Q C Q C Q C Q + + =
s m Q
Q
/ 89 . 13
) 10 )( 2329 . 0 ( ) 18 )( 4863 . 0 ( ) 10 )( 7466 . 0 (
3 1
1
1
1
=
+ + =
Time (hr) 0 km 6 km 12 km 18 km
0 10 10 10 10
2 18 13.89
4 50
6 107
8 147
10 146
12 105
14 59
16 33
18 17
20 10
22 10
24 10
26 10
28 10
62
Muskingum-Cunge Example
Therefore, the equation can be solved.
1
1 3
2
0 2
1
0 1
2
1
Q C Q C Q C Q + + =
s m Q
Q
/ 51 . 34
) 89 . 13 )( 2329 . 0 ( ) 50 )( 4863 . 0 ( ) 18 )( 7466 . 0 (
3 1
1
1
1
=
+ + =
Time (hr) 0 km 6 km 12 km 18 km
0 10 10 10 10
2 18 13.89
4 50 34.51
6 107
8 147
10 146
12 105
14 59
16 33
18 17
20 10
22 10
24 10
26 10
28 10
63
Muskingum-Cunge Example
This is repeated for the rest of the columns and the subsequent
columns to produce the following table. Note that when you
change rows, n changes. When you change columns, j
changes.
Time (hr) 0 km 6 km 12 km 18 km
0 10 10 10 10
2 18 13.89 11.89 10.92
4 50 34.51 24.38 18.19
6 107 81.32 59.63 42.96
8 147 132.44 111.23 88.60
10 146 149.91 145.88 133.35
12 105 125.16 138.82 145.37
14 59 77.93 99.01 117.94
16 33 41.94 55.52 73.45
18 17 23.14 29.63 38.75
20 10 12.17 16.29 21.02
22 10 9.49 9.91 12.09
24 10 10.12 9.70 9.30
26 10 9.97 10.15 10.01
28 10 10.01 9.95 10.08
64
Full Dynamic Wave Equations
The solution of the St. Venant equations is known as
dynamic routing.
Dynamic routing is generally the standard to which
other methods are measured or compared.
The solution of the St. Venant equations is generally
accomplished via one of two methods : 1) the method
of characteristics and 2) direct methods (implicit and
explicit).
It may be fair to say that regardless of the method of
solution, a computer is absolutely necessary as the
solutions are quite time consuming.
J. J. Stoker (1953, 1957) is generally credited for
initially attempting to solve the St. Venant equations
using a high speed computer.
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
65
Dynamic Wave Solutions
Characteristics, Explicit, & Implicit
The most popular method of applying the implicit
technique is to use a four point weighted finite
difference scheme.
Some computer programs utilize a finite element
solution technique; however, these tend to be more
complex in nature and thus a finite difference technique
is most often employed.
It should be noted that most of the models using the
finite difference technique are one-dimensional and that
two and three-dimensional solution schemes often revert
to a finite element solution.
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
66
Dynamic Wave Solutions
Dynamic routing allows for a higher degree of accuracy
when modeling flood situations because it includes
parameters that other methods neglect.
Dynamic routing, when compared to other modeling
techniques, relies less on previous flood data and more
on the physical properties of the storm. This is
extremely important when record rainfalls occur or
other extreme events.
Dynamic routing also provides more hydraulic
information about the event, which can be used to
determine the transportation of sediment along the
waterway.
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
67
Courant Condition?
If the wave or hydrograph can travel through the
subreach (of length x) in a time less than the
computational interval, t, then computational
instabilities may evolve.
The condition to satisfy here is known as the Courant
condition and is expressed as :
c
dx
dt s
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes
68
Some DISadvantages
Geometric simplification - some models are designed to use
very simplistic representations of the cross-sectional
geometry. This may be valid for large dam breaks where very
large flows are encountered and width to depth ratios are
large; however, this may not be applicable to smaller dam
breaks where channel geometry would be more critical.
Model simulation input requirements - dynamic routing
techniques generally require boundary conditions at one or
more locations in the domain, such as the upstream and
downstream sections. These boundary conditions may in the
form of known or constant water surfaces, hydrographs, or
assumed stage-discharge relationships.
Stability - As previously noted, the very complex nature of
these methods often leads to numeric instability. Also,
convergence may be a problem in some solution schemes.
For these reasons as well as others, there tends to be a stability
problem in some programs. Often times it is very difficult to
obtain a "clean" model run in a cost efficient manner.
Modified Puls
Kinematic Wave
Muskingum
Muskingum-Cunge
Dynamic
Modeling Notes