Gradually varied flow

Uniform flow is an ideal reached only far from any disturbance. Near a dam, a lake, a weir, a change of slope or a culvert, the depth departs from normal and changes gradually along the channel, often over hundreds of metres or more. This is gradually varied flow, and computing its water surface profile tells the engineer how high floods will rise at a building, a bridge or an outfall. For Kestrel Bay the question is concrete: when Lake Superior is high and the creek is in flood, how far upstream does the lake raise the water, and does it reach the pond outfall or the culvert?

Assumptions

The theory of gradually varied flow rests on three assumptions. The flow is steady. The streamlines are nearly parallel, so the pressure is hydrostatic. And the head loss at any section is the same as it would be in uniform flow with the same depth and velocity, so that Manning's equation gives the friction slope:

The last assumption cannot be proved, but decades of measurements show that it works well where the depth changes slowly.

The equation of gradually varied flow

Differentiate the total head \(H = z + y + V^2/2g\) with respect to distance \(x\) along the channel. The bed slope is \(S_0 = -dz/dx\), the friction slope is \(S_f = -dH/dx\), and lesson 2 showed that \(d(y + V^2/2g)/dy = 1 - Fr^2\). Combining these gives

This one equation explains every profile. The numerator compares the bed slope with the friction slope. When the depth is greater than normal, the velocity is lower than in uniform flow, so \(S_f < S_0\) and the numerator is positive; when the depth is less than normal, \(S_f > S_0\) and the numerator is negative. The denominator compares the depth with the critical depth: positive for subcritical flow, negative for supercritical flow. The sign of \(dy/dx\) therefore depends only on where the actual depth lies relative to the normal and critical depths.

Classifying the profiles

Channel slopes are classified by comparing the normal depth with the critical depth: mild (M) if \(y_n > y_c\), steep (S) if \(y_n < y_c\), critical (C) if they are equal, horizontal (H) if \(S_0 = 0\) and adverse (A) if the bed rises in the direction of flow. The normal depth line (NDL) and the critical depth line (CDL) divide the space above the bed into three zones: zone 1 above both lines, zone 2 between them and zone 3 below both. Each profile is named by the slope letter and the zone number.

On a mild slope, which covers most rivers and creeks, three profiles occur:

• M1, in zone 1, above the normal depth. The depth increases downstream, approaching the normal depth asymptotically upstream and a horizontal surface downstream. It is the backwater curve found upstream of a dam, a weir or a high lake level, and it is by far the most common profile in flood studies.
• M2, in zone 2, between the normal and critical depths. The depth decreases downstream. It is the drawdown curve upstream of a free overfall or a sudden steepening of the channel, ending at critical depth.
• M3, in zone 3, below critical depth. The flow is supercritical, and the depth increases downstream until a hydraulic jump returns it to subcritical flow. It occurs below a sluice gate or at the outlet of a steep culvert discharging onto a mild channel.

Steep slopes have the corresponding S1, S2 and S3 profiles. The S1 profile, for example, forms when a steep channel runs into a pool: a jump occurs upstream, and the subcritical depth then rises toward the pool level.

Two rules make profile sketching reliable. Profiles approach the normal depth line asymptotically, because \(dy/dx \to 0\) as \(S_f \to S_0\). And they cross the critical depth line at right angles in theory, because \(dy/dx \to \infty\) as \(Fr \to 1\); in reality, the assumptions fail near critical depth, and the flow passes through it at a control or a jump.

Where to start

Lesson 1 noted that subcritical flow is controlled from downstream and supercritical flow from upstream. Profile computations follow this rule: an M1 or M2 profile is computed upstream from a downstream control, such as a lake level or the critical depth at a brink, and an M3 or S2 profile is computed downstream from an upstream control, such as the depth under a gate. Before any calculation, the engineer should sketch the expected profiles, mark the controls and decide the direction of each computation.

The direct step method

For a prismatic channel, one with constant shape and slope, the simplest method of computing a profile is the direct step method. The depth is divided into small steps, and the distance between each pair of depths is found from the energy equation written as

where \(E\) is the specific energy and \(\bar{S}_f\) is the average friction slope of the two sections. With \(x\) measured upstream, the sign convention is arranged so that \(\Delta x\) is positive; it is simplest to write the formula so that it gives the distance upstream directly, as Sofia's script does. Smaller depth steps give better accuracy, especially near the normal depth, where the profile flattens and the distances grow large.

The standard step method

Natural channels are not prismatic. Their shape, slope and roughness change from one surveyed cross section to the next, so the distances are fixed and the depths must be found. The standard step method does this by trial: starting at a section with a known water surface elevation, the engineer guesses the elevation at the next section upstream, computes the friction loss between the two with the average friction slope, adds expansion or contraction losses, and adjusts the guess until the energy equation balances. This is what HEC-RAS does for every cross section in a floodplain model, and it is the basis of most floodplain mapping in Canada and the United States.

The backwater from Lake Superior

The creek enters Lake Superior over a bed about 2.6 m below the 100-year lake level, which the Lakehead Region Conservation Authority sets at 184.0 m along this part of the shore. Anika takes that as the downstream control. Whether a 100-year lake level should be combined with a 100-year creek flood is a question of joint probability: the two are caused by different weather, so their coincidence is far rarer than either alone. For this check she uses the combination anyway as a conservative bound, and notes that the conservation authority's guidance should be consulted for the combination it expects.

Lake Superior's levels are regulated, within natural limits, by outflow control works at its outlet on the St. Marys River under the International Joint Commission, but the lake still varies by a metre or so over the years and responds to wind setup during storms. The design high level includes these effects.

Sofia's script continues the steps until the depth is within 2 % of normal, printing every fourth step, and interpolates the depth at the pond outfall and at the culvert outlet.

# Gradually varied flow: the backwater profile from Lake Superior up Kestrel Creek (direct step method)
from math import sqrt
UNITS = "SI"
SI = UNITS == "SI"
g, k = (9.81, 1.0) if SI else (32.2, 1.49)
b, z, n, S0, Q = (3.0, 2.0, 0.035, 0.004, 18.0) if SI else (10.0, 2.0, 0.035, 0.004, 640.0)
L = "m" if SI else "ft"
yn, yc = (1.615, 1.182) if SI else (5.286, 3.872)          # from lessons 2 and 3
y_lake = 2.60 if SI else 8.50                               # depth at the mouth with a high lake level
spots = [("pond outfall", 250.0 if SI else 820.0), ("culvert outlet", 1200.0 if SI else 3940.0)]   # distances from the mouth
def geo(y):
    A = (b + z * y) * y; P = b + 2 * y * sqrt(1 + z * z); return A, A / P, b + 2 * z * y
def Sf(y):
    A, R, T = geo(y); return (n * Q / (k * A * R ** (2 / 3))) ** 2
def E(y):
    A, R, T = geo(y); return y + Q * Q / (2 * g * A * A)
def kind(y):
    zone = 1 if y > max(yn, yc) else (2 if y > min(yn, yc) else 3)
    return ("M" if yn > yc else "S") + str(zone)
print(f"yn = {yn} {L} > yc = {yc} {L}: mild slope; depth {y_lake} {L} at the mouth gives an {kind(y_lake)} profile")
print(f"{'y':>6} {'E':>7} {'Sf':>9} {'dx':>8} {'x':>8}")
ys = []; y = y_lake
while y > 1.02 * yn:
    ys.append(y); y -= (y_lake - 1.02 * yn) / 24
if ys[-1] - 1.02 * yn > 1e-6: ys.append(1.02 * yn)
else: ys[-1] = 1.02 * yn
x = 0.0; prof = [(ys[0], 0.0)]
print(f"{ys[0]:6.3f} {E(ys[0]):7.3f} {Sf(ys[0]):9.6f} {'':>8} {0.0:8.0f}")
for y1, y2 in zip(ys, ys[1:]):
    dx = (E(y1) - E(y2)) / (S0 - (Sf(y1) + Sf(y2)) / 2)    # distance upstream (positive)
    x += dx; prof.append((y2, x))
    if len(prof) % 4 == 1 or y2 == ys[-1]:
        print(f"{y2:6.3f} {E(y2):7.3f} {Sf(y2):9.6f} {dx:8.1f} {x:8.0f}")
print()
for name, xs in spots:
    yx = yn
    for (ya, xa), (yb, xb) in zip(prof, prof[1:]):
        if xa <= xs <= xb: yx = ya + (yb - ya) * (xs - xa) / (xb - xa)
    print(f"Depth at the {name}, {xs:.0f} {L} upstream: {yx:.3f} {L}" + ("  (beyond the backwater: normal depth)" if xs > prof[-1][1] else ""))
print(f"Backwater reaches within 2 % of normal depth {x:.0f} {L} upstream of the mouth")

yn = 1.615 m > yc = 1.182 m: mild slope; depth 2.6 m at the mouth gives an M1 profile
     y       E        Sf       dx        x
 2.600   2.636  0.000528                 0
 2.441   2.486  0.000696     11.2       44
 2.282   2.338  0.000931     11.8       91
 2.124   2.193  0.001269     12.9      140
 1.965   2.054  0.001766     14.9      196
 1.806   1.922  0.002518     20.2      267
 1.647   1.801  0.003690     59.9      415

Depth at the pond outfall, 250 m upstream: 1.840 m
Depth at the culvert outlet, 1200 m upstream: 1.615 m  (beyond the backwater: normal depth)
Backwater reaches within 2 % of normal depth 415 m upstream of the mouth

Reading the results

The backwater dies out quickly. The depth falls to within 2 % of the normal depth about 415 m upstream of the mouth. At the pond outfall, 250 m upstream, the depth is 1.84 m, about 0.23 m above normal, so the pond outlet pipe will be partly submerged during a coincident flood and high lake; lesson 12 checks the effect on the pond. At the culvert, 1.2 km upstream, the lake has no effect: the tailwater is the creek's own normal depth of 1.615 m.

The length of a backwater depends mainly on the slope. A rough estimate of its length is the extra depth at the control divided by the bed slope, which here gives about one metre of extra depth over a slope of 0.004, or a few hundred metres, consistent with the computation. On a flat river with a slope of 0.0002, the same extra depth would affect the water level many kilometres upstream, which is why dams on large rivers raise flood levels far from the dam.

Profiles in culverts and channels

Gradually varied flow also appears inside long culverts and storm sewers. A culvert on a mild slope with a low tailwater develops an M2 drawdown toward critical depth at the outlet; with a high tailwater, an M1 profile backs up through the barrel and may fill it. A culvert on a steep slope runs supercritical with an S2 profile from critical depth at the entrance. Recognizing which applies is the key to the outlet control analysis in lesson 12.

Practical points

• Sketch the profile, mark the normal and critical depth lines and the controls, and decide the direction of computation before calculating.
• Use small depth steps near the normal depth, where distances grow rapidly.
• Use the standard step method, or HEC-RAS, for natural channels with surveyed cross sections.
• Think about the joint probability of a high downstream level and a high flow; use agency guidance for the combination.
• A backwater on a steep creek is short; on a flat river it can extend for kilometres.