- Email: [email protected]
Hydrodynamics of storage release during river ice breakup
Accepted Manuscript Hydrodynamics of storage release during river ice breakup
Spyros Beltaos PII: DOI: Reference:
S0165-232X(16)30416-5 doi: 10.1016/j.coldregions.2017.04.009 COLTEC 2388
To appear in:
Cold Regions Science and Technology
Received date: Revised date: Accepted date:
19 December 2016 31 March 2017 29 April 2017
Please cite this article as: Spyros Beltaos , Hydrodynamics of storage release during river ice breakup, Cold Regions Science and Technology (2017), doi: 10.1016/ j.coldregions.2017.04.009
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
ACCEPTED MANUSCRIPT Hydrodynamics of storage release during river ice breakup Spyros Beltaos
IP
T
Watershed Hydrology and Ecology Research Division, Canada Centre for Inland Waters, Environment and Climate Change Canada, Burlington, Ontario, Canada [email protected]
CR
ABSTRACT
AC
KEY WORDS
CE
PT
ED
M
AN
US
The breakup of the ice cover in cold-region rivers is a brief but seminal period of their hydrologic regime, with important ecological and socio-economic implications. The main driver of ice breakup processes is the flow discharge hydrograph. It is generated by runoff from snowmelt and rainfall but can be modified by rapid release of water from storage as the ice cover recedes by ablation and mechanical breaking up. Despite its potential importance, there is very limited quantitative information concerning the hydrodynamic processes that control storage release during ice breakup in rivers, while the issue of climate change underscores the need for improved understanding of the relevant mechanisms. Quantitative analysis for assumed prismatic channels shows that ablation and sustained ice dislodgment and breaking can cause significant flow enhancement via storage release. The latter process is far more dynamic than ablation and, under certain conditions, may lead to formation of a self-sustaining wave (SSW). Analytical results are applied to natural stream conditions, using the Lower Peace River as a case study. Observed rates of ice recession typically indicate ice melt as the dominant process. A rare occurrence of rapid ice breaking over hundreds of kilometres (2014) indicated partial agreement with the SSW concept and revealed that discrepancies may arise from the characteristic irregularity of rivers. Impacts on storage release by climate-driven changes to river ice regimes are examined, along with their implications to ice-jam formation and associated flooding.
Breakup, celerity, climate, flow enhancement, ice cover, melt, self-sustaining wave, storage release
1
ACCEPTED MANUSCRIPT 1. INTRODUCTION The breakup of the ice cover in cold-region rivers is a brief but seminal period of their hydrologic regime with important ecological and socio-economic implications (Burrell 2008). Major ice jams that often form during breakup can cause extreme floods that imperil human life
T
and damage property and infrastructure (Ashton 1986). Navigation and hydro-power generation
IP
are also adversely affected by ice jamming. At the same time, ice-jam flooding can provide
CR
essential replenishment to the multitude of lakes and ponds characteristic of the northern Canadian deltas, which are havens for wildlife, especially waterfowl and aquatic animals. On the
US
other hand, the sharp waves produced by ice jam releases (also termed “javes”) can have
AN
negative ecological impacts, owing to amplified velocity and erosive capacity (Beltaos and
M
Burrell 2016, Beltaos 2016).
ED
The main driver of ice breakup processes is the flow discharge hydrograph. It is generated by runoff from snowmelt and rainfall but can be modified by highly dynamic processes resulting
PT
from the rapid release of water from storage, as the ice cover recedes by ablation and mechanical
CE
breaking up. Storage release has been postulated to play a role in rapid ice clearance that has occasionally been observed in rivers of northern Canada. An extreme occurrence took place on
AC
the Yukon River in 1983: Gerard et al (1984) reported apparent ice breaking over a 300 km stretch at an average celerity of 5.2 m/s, which reduced the intact ice cover to a 140 km long ice run. Such events could potentially be explained by formation and advance of non-attenuating waves that are sustained by storage release as the ice cover is being broken up. This effect was detected in exploratory numerical runs (Ferrick and Mulherin 1989) and quantified, also numerically, by Jasek et al (2005). The results of the latter study are used extensively in Section
2
ACCEPTED MANUSCRIPT 4 and Appendix B, where the analysis and theory of the self-sustaining wave (SSW for short) are developed.
The main objective of this paper is to systematically examine and, where possible, quantify
T
storage release mechanisms and their effects on flow magnitude and breakup progression. A
IP
secondary objective pertains to how changes in ice regime that may be caused by climatic
CR
warming can alter storage release processes.
US
Following a background section, storage release by melting of the ice cover is quantified for
AN
different water temperatures and channel characteristics. The much more dynamic process of storage release by ice breaking is considered next and the potential formation of self-sustaining
M
waves illustrated. Unique properties of the SSW are discussed and the links between ice-
ED
breaking criteria and wave celerity elucidated. Application to natural streams is illustrated next
PT
with a case study, followed by a discussion of potential impacts of changing ice regimes.
CE
2. BACKGROUND INFORMATION As has been described in previous literature (Gerard 1990, Prowse and Carter 2002, Beltaos
AC
2009), the formation of an ice cover causes river water levels to increase in order to accommodate the keel of the cover and the hydraulic resistance generated by the bottom boundary of the cover, in addition to that of the river bed. The volume of water between the prevailing and the initial (before ice cover formation) water levels is called “storage”, even though the “stored” water is not stationary but moves along the river with the rest of the flow,
3
ACCEPTED MANUSCRIPT albeit at reduced velocity. As water is being stored during the freezeup period, downstream flow is reduced and this effect is often captured in the records of various hydrometric gauges.
During the winter, the river flow typically decreases and so does the roughness of the bottom
T
surface of the ice cover. In turn, water levels decrease and some of the stored water is released
IP
back into the flow. By the time spring runoff begins, there is a relatively small amount of water
CR
in storage, since the prevailing flow and thence the “backwater” is small. The term backwater
US
denotes the difference between ice-covered and open-water flow depth at the same discharge.
AN
With the arrival of mild weather in the spring (or during a mid-winter thaw), the flow begins to increase while thermal effects cause the ice bottom surface to become rougher (Carey 1966,
M
1967). The backwater will increase as a result and so will the storage until the ice cover melts or
ED
is forcibly mobilized by the rising flow. In either case, there will be a release of storage, which will cause the downstream flow to increase over and above the “carrier” flow in the river. Carrier
PT
flow is the flow that is strictly attributable to runoff and is the flow that prevails in the river
CE
before dynamic effects, such as javes and storage-release, appear or after they cease. In the
AC
following, the term “excess” flow will be used to denote the increment caused by storage release.
When a segment of the ice cover is removed, flow resistance is significantly reduced. Local velocity and discharge increase, and thence “stored” water is released. To visualize this effect, a numerical experiment was carried out using the River1d model (2004 version; kindly supplied to the writer by Professor Faye Hicks). For simplicity, the open water version of the model was used while the backwater effect of the ice cover was simulated by doubling the bed Manning
4
ACCEPTED MANUSCRIPT coefficient, nb, from 0.03 to 0.06 in the “ice-covered” reach of a 32 km long prismatic rectangular channel (flow = 343 m3/s. width = 100 m, slope = 0.0005). Using Eq. A1 of Appendix A, it can be shown that open-water flow with nb = 0.06 is equivalent, in terms of depth and mean velocity, to ice-covered flow with nb = 0.03 and ni = Manning coefficient of the ice
IP
T
cover = 0.045).
CR
The “ice cover” extended from km 10 to 32. The resulting water level profile (Fig. 1) resembles what would be expected if an ice cover were present at and beyond km 10. At time 0, the
US
Manning coefficient is reduced to 0.03 (from 0.06) between kilometres 10 and 10.5, simulating
AN
removal of a 0.5 km segment of ice cover. This change generates a small wave (Fig. 2) that ultimately leads to an asymptotic configuration, which corresponds to a 0.5 km forward
M
translation of the initial water surface profile (Fig. 1). The computation indicated that it took 3.5
ED
h to attain a profile that is fairly close to the asymptotic one (Fig. 2). Moreover, it was found that the discharge increased along the wave-affected segment of the channel, always attaining a peak
PT
at km 10; the highest such peak (377 m3/s) occurred at ~0.1 h. The peak value gradually
AC
CE
decreased afterward, falling to 350 m3/s at 5.0 h.
5
ACCEPTED MANUSCRIPT 104 Channel bed
102
Water surface
98
T
96
IP
94
92 90 5
10
15
Channel distance (km)
20
US
0
CR
Elevation (m)
100
AN
Figure 1. Initial flow configuration used to numerically simulate the effects of storage release. Discharge = 343 m3/s, width = 100 m, slope = 0.0005, nb = 0.03 (0 to 10 km), 0.06 (10 to 32 km).
M
0.15
ED PT
0.05
CE
0.00 -0.05 -0.10
-0.15
0.1 h
0.4 h 1.5 h
AC
Wave disturbance (m)
0.10
3.5 h 6.0 h Asymptotic
-0.20 0
5
10
15
20
25
30
Channel distance (km)
Fig. 2. Spatial variation of change in water level at different times following simulated removal of a 0.5 km long “ice cover” segment from the channel of Fig. 1. 6
ACCEPTED MANUSCRIPT 3. ICE MELT: THERMAL BREAKUP Low spring runoff typically results in thermally driven ice clearance. The flow is not sufficiently high to cause mechanical breaking and mobilization of the ice cover. Instead the cover largely melts in place. The resulting process of storage release is discussed in this section by considering
T
the simplified situation depicted in Fig. 3. Steady flow in an open-water reach of a prismatic
IP
rectangular channel is subjected to thermal inputs and the water temperature increases with
CR
downstream distance until the flow encounters a long stretch of stationary ice cover. As the water dives under the ice cover, its temperature decreases gradually to zero because heat is
US
extracted by ice melt. Assuming that the incoming flow and water temperature do not change
AN
over time, the short transitional reach AB will simply be translated downstream at a rate (Hicks
ED
c pQTW i LFW
(1)
PT
Cmelt
M
et al 2008):
CE
in which TW = temperature of the water just upstream of the edge of the ice cover; Q = river discharge arriving at the ice edge; cp = specific heat of water (4,220 J/kgC at 0C); = thickness
AC
of the ice cover (m); = density of water (1000 kg/m3); i = density of freshwater ice (916 kg/m3); LF = latent heat of fusion (3.34105 J/kg); and W = width of the ice cover.
7
US
CR
IP
T
ACCEPTED MANUSCRIPT
AN
Figure 3. Schematic illustration of a wave generated by ice melt. The wave front is depicted as a straight line for simplicity. Melting results in zero thickness at the ice edge. The symbols C, Q, and Y represent celerity, discharge, and water depth, respectively.
M
The melting of the ice cover can be viewed as a series of successive removals of short ice
ED
segments at equal time intervals, such that each removal generates an elementary wave similar to that of Fig. 2. The resulting composite wave was numerically simulated by removing 0.5 km
PT
“ice” segments every 0.2 h (12 min) from the initial profile of Fig. 1. The speed of ice removal is
CE
thus 2.5 km/h or 0.69 m /s; this is a relatively large value that, depending on channel conditions, would require a water temperature of a few to several degrees C. The resulting wave forms at
AC
different times are illustrated in Fig. 4. The peak of the wave always coincides with the edge of the ice cover, that is, it moves at the melting rate. At the same time, the speed of advance of the leading edge of the wave is considerably larger (1.9 m/s indicated by the first few waves).
Owing to the relatively short length of the channel, the computed shape of the wave front was, after ~3 h, influenced by the downstream boundary condition (fixed water depth). Nevertheless, the results clearly indicate that the wave height tends towards an asymptotic value that is nearly 8
ACCEPTED MANUSCRIPT attained within hours from the start of melt. At the same time, the wave front elongates, owing to the faster advance of the leading edge relative to the ice edge. In turn, this implies that the water surface slope of the wave front will eventually approach the value of the bed slope (So) of the channel. It is only after the time when this configuration is established that the flow downstream
IP
T
of the ice edge can be considered approximately uniform with a flow equal to the peak value.
CR
5
US AN M
3
2
ED
Depth (m)
4
0
5
AC
0
CE
PT
1
Initial profile profiles at 0.6, 1.6, 2.4, 3.0, 3.8 h asymptote
10 15 20 River distance (km)
25
30
Figure 4. Temporal development of wave generated by sequential release of 0.5 km “ice segments” at 12 minute (0.2 h) intervals. The wave peak is always located at the “ice edge” and tends to an asymptotic value of 4.71 m. The “ice cover” is simulated by simply increasing the Manning coefficient from 0.03 to 0.06.
To determine the asymptotic flow properties, one may start with the differential equation of continuity:
9
ACCEPTED MANUSCRIPT
Y q 0 t x
(2)
T
in which t = time, x = downstream distance; Y = flow depth (measured from the water surface or
IP
from the bottom of the ice cover, as applicable); and q = discharge per unit channel width = Q
CR
(discharge)/W (width). When a waveform propagates at constant rate, C, a single coordinate can
US
be defined as = x-Ct and Eq. 2 transformed to an ordinary differential equation, i.e.
(3)
M
AN
d (CY q) 0 d
which implies that the quantity CY-q is constant. Because C changes downstream of the peak of
ED
the wave, it is only the portion upstream of the peak where Eq. 3 applies. The largest excess flow
PT
due to storage release occurs at the ice edge; it will remain approximately equal to this value downstream of the edge, once an adequate time has elapsed since the start of melt to ensure that
AC
CE
the frontal water surface slope is approximately equal to So.
Application of Eq. 3 between the ice edge and a section far upstream, where uniform open-water flow prevails, results in:
Q Qe Qu Cmelt (Ye Yu )W
(4)
10
ACCEPTED MANUSCRIPT in which Q is the excess flow (storage-release contribution to discharge) and the subscripts “e” and “u” denote edge and far upstream conditions. By defining Ye as water depth (water surface to channel bed), the volume of melted ice is included in Qe. One may note here that the depth Ye is greater than what prevailed in the ice-covered reach before the start of melt (Fig. 4) because it
T
is accommodating a greater discharge than the carrier discharge, Qu. In general, therefore, it is
IP
not possible to determine Qe from Eq. 4 alone, because both Ye and Qe are unknown. However,
CR
when the condition of near-uniformity of the flow downstream of the ice edge is established, an
US
additional relationship may be obtained via the Manning resistance equation:
(5)
M
AN
W (Ye si )5/ 3 So Qe 22 / 3 nc
ED
in which si = specific gravity of ice = 0.916; and nc is the composite-flow Manning coefficient,
CE
Sabaneev formula:
PT
which can be computed in terms of ni and nb using the simplified version of the Belokon-
(6)
AC
nc [(ni 3/ 2 nb3/ 2 ) / 2]2 / 3
Noting that Qu and Yu are related by a similar relationship to that of Eq. 5, it is possible to rearrange Eq. 4 into a dimensionless form, i.e:
3/ 5 Q Qe Cmelt 22 / 3 nc Qe si 1 1 Qu Qu U u nb Qu Yu
(7)
11
ACCEPTED MANUSCRIPT
in which Uu = mean velocity in the uniform-flow, open-water reach. In Eq. 1, the amount of heat arriving at the edge of the ice cover is expressed in terms of local discharge and local water temperature. In the present application, the excess discharge has zero temperature, since it
T
derives from under-ice storage. Consequently, the heat available for melt is taken as (cpQTW)u,
IP
i.e at point A rather than at point B where the discharge will be higher than Qu but the water
CR
temperature lower than the temperature at point A. [This approximation neglects any heat that
US
may be gained in the transitional river segment between points A and B in Fig. 3, which is considered small but can be taken into account using more elaborate equations]. Equation 7 can
AN
be solved by trial and error to determine Qe and thence the excess flow Q.
M
As an example, consider a 600 m wide channel with carrier flow of 3000 m3/s, slope of 0.00005,
ED
Manning coefficients nb = ni = nc= 0.03, and ice thickness = 0.7 m. This type of channel could be
PT
viewed as a prismatic approximation to the lower Peace River (Kellerhals et al 1972; Jasek and Pryse-Phillips 2015). Typically, the winter ice cover is smoother than the river bed but becomes
CE
rougher during the pre-breakup period, so that ni approaches nb (Carey 1966, 1967). This
AC
assumption is adopted here for illustration purposes; in real-life applications actual value of these coefficients should be determined based on in situ hydraulic data. By application of the flow resistance equation to the upstream open-water flow, the depth Yu is calculated as [0.033000/6000.00005]0.6 = 6.25 m; and thence Uu = 3000/600/6.25=0.80 m/s. The rate of melt, Cmelt, is then calculated from Eq. 1 as ~0.10 m/s (8.6 km/day) for each degree of the water temperature.
12
ACCEPTED MANUSCRIPT
40 bed slope = 0.00005 bed slope = 0.0005
T IP
20
CR
Excess flow (%)
30
US
10
0
1
AN
0
2
3
4
5
6
M
Water temperature (°C)
PT
ED
Figure 5. Addition to carrier flow caused by melt-induced storage release in two hypothetical channels, both with 0.7 m thick ice covers. Channel 1: So = 0.00005, W = 600 m, Qu = 3000 m3/s, Cmelt = 8.6 km/day per degree of water temperature; Channel 2: S o = 0.0005, W = 200 m, Qu = 800 m3/s, Cmelt = 6.8 km/day per degree of water temperature
CE
With the aid of Eq. 7, the excess flow can be calculated for different water temperatures and is plotted as a percentage in Fig. 5, where it is seen to exceed 10% for Tw > 2°C. Also shown in
AC
Fig. 5 are the results of a second calculation, this time for a channel that approximates average properties of the upper Saint John River. The higher slope greatly reduces the excess flow, which remains under 10% of the carrier flow, even if the water temperature were to increase to 5°C. In situ measurements that were carried out in upper Saint John River during the breakup in the 1990s have actually indicated water temperatures not exceeding 2°C. On the other hand, the large solar radiation flux associated with spring conditions at high latitudes and the extensive lengths of major rivers of Northern Canada suggest that water temperatures of several degrees C 13
ACCEPTED MANUSCRIPT may be realistic. Measurements in the lower Mackenzie River indicated values up to 10°C during the highly thermal breakup event of 1980 (Parkinson 1982). Considering that the excess flow can be substantial in the low-slope channel (Fig. 5), it is conceivable that the excess may enable
T
mechanical ice breaking farther downstream.
IP
To determine the time after which Eqs. 5 and 7 can be applied with some confidence beyond the
Ye Yu 1) Ye Yi
US
Cmelt (2
(8)
AN
CL
CR
ice edge, the leading-edge celerity CL is estimated from the approximate relationship:
M
in which Yi is the water depth in the ice-covered reach far downstream from the ice edge where
ED
the water level has not yet been affected by the wave; and the average celerity for the wave front is assumed to be simply equal to (Cmelt+CL)/2. To attain a wave slope that is within 5% of So, the
20
Ye Yi CL S o
(9)
AC
TL
CE
PT
leading edge of the wave would have to travel for a time TL:
As an example, the value of Yi for the low-slope channel is calculated as 8.89 m, while Ye = 9.37 m for TW = 2°C, which corresponds to Cmelt = 0.2 m/s. With Yu = 6.25 m, as noted earlier, CL works out to be 2.4 m/s and thence TL = 20 (9.37-8.89)/(2.40.000053600) ~ 22 h after the start of melt. During this time interval, the ice edge would be located 2236000.2/1000 ~ 16 km downstream from its initial location. 14
ACCEPTED MANUSCRIPT
Of course, the preceding calculations are predicated on constant water temperature at the ice edge. But as the open-water reach lengthens, the water temperature will increase and thence the melting rate would be expected to increase as well. Variable weather conditions, ice thickness,
T
and runoff intensity will further complicate the ablation-driven storage release process.
IP
Nevertheless, the present examples capture the main elements of the process and elucidate the
CR
magnitudes of the important variables such as the melt rate and the resulting excess discharge. The analysis also shows that considerable melting time must elapse before the flow underneath
AN
US
the ice cover becomes nearly uniform.
4. ICE DISLODGMENT AND BREAKING: THE SELF-SUSTAINING WAVE
M
Much more dynamic conditions will prevail in cases where the rising flow dislodges and
ED
mobilizes the ice cover, which then quickly breaks down into smaller slabs and blocks. At the same time, the remaining intact ice cover maintains a good portion of its mechanical strength
PT
(“mechanical” breakup) and can, under certain geomorphic conditions, arrest the ice run and
CE
initiate an ice jam. As stationary ice is being mobilized, new water volume is released from hydraulic storage and added to the flow. This effect counteracts the natural attenuation of the
AC
resulting wave and has been postulated via numerical experiments to render it self-sustaining (Ferrick and Mulherin 1989). Self-sustaining waves could perhaps explain reports of rapid breakup progression over very long reaches (Gerard et al. 1984).
Jasek et al. (2005) applied the CRISSP model to study wave propagation in ice-covered prismatic rectangular channels, assuming that ice dislodgment and breaking occur at a specified
15
ACCEPTED MANUSCRIPT discharge (Qbr), which exceeded the base flow (Qo) by varying amounts. Ice cover dislodgment and subsequent breaking into ice blocks is more complex than what is envisaged by singlevariable thresholds. It involves additional parameters, e.g. hydrodynamic forces applied on the ice cover, channel width and planform curvature as well as freezeup water level (Beltaos 2007,
T
2013b). Nevertheless, the single-variable criterion greatly simplifies the study of storage release
IP
by ice breaking and is rigorous for the assumed prismatic channel, for which the critical
CR
condition could equally well be expressed in terms of velocity, depth, or shear stress.
US
Dynamic conditions were generated by applying carrier flow hydrographs of different peaks and
AN
durations at the upstream boundary of the computational reach, which extended for hundreds of kilometres. If the incoming wave attained flows in excess of Qbr, the wave became self-
M
sustaining after an evolutionary period; it comprised a steep front of constant celerity, followed
ED
by a protracted crest of constant depth and discharge (herein termed sustained peak depth and sustained peak flow). The numerical runs indicated that the celerity of the front (Cbr) increased as
PT
Qbr decreased, becoming extremely large (~10 m/s) when Qbr was just a little higher than the base
AC
the channel width.
CE
flow, denoted by Qo. The results were also influenced by the slope of the channel bed but not by
The waves studied by Jasek et al (2005) were triggered by a flow hydrograph rather than by javes that result from ice jam releases, as is common in relatively flat rivers (Beltaos and Carter 2009; Beltaos 2013a). Some of these hydrographs were very steep while others were more gradual, approximating snowmelt runoff. Such differences did not seem to influence the properties of the SSWs, which often attained higher peak flows than the peaks of the respective
16
ACCEPTED MANUSCRIPT input hydrographs. To examine whether jave-induced ice breaking can also lead to formation of a SSW, the writer used the model River1d to study the release of a 9 km long jam in a 120 km prismatic channel of slope 0.5 m/km, width 100 m, base flow of 400 m3/s, and bed Manning coefficient, nb = 0.03. To simulate the backwater caused by the ice cover downstream of the jam,
T
nb was increased to 0.05 while the calculation was performed in open-water mode throughout the
IP
length of the channel (equivalent to ice-covered flow with nb = 0.03 and ni = 0.033; see also
CR
Appendix A). The shear stress in the “ice-covered” reach was calculated as one-half of gYSf (Sf
US
= friction slope; = density of water = 1000 kg/m3, g = gravitational acceleration = 9.81 m/s2). When the shear stress at both ends of the first 1 km long “ice segment” exceeded 20 Pa, the
AN
segment was removed, i.e. nb was reduced to 0.03 and the computation repeated using the latest water level profile as the initial condition. The 20 Pa threshold is higher than the pre-jave value
ED
M
but low enough to ensure that the jave could initiate ice breaking.
PT
Figure 6 shows that the jave quickly formed a front of constant shape and height, unlike the attenuating configurations of javes that do not dislodge the ice cover (Beltaos 2013b).
CE
Instantaneous depth profiles at different times from the release of the jam are shown in Figure 7.
AC
After 46 minutes, a short segment of practically constant depth has already developed upstream of the wave front; as time goes on, it is expected to get longer, much as was found by Jasek et al (2005). The sustained peak flow and depth were smaller than respective values associated with the initial jave (e.g. 5-minute graph). The magnitude of the celerity of the wave front (5.1 m/s) lies between those of the gravity and the kinematic waves (7.1 m/s and 1.6 m/s, respectively).
17
ACCEPTED MANUSCRIPT 205 initial WL Riverbed
200
T IP
190
CR
Elevation (m)
195
US
185
0
5
10
AN
180 15
20
25
30
35
40
45
M
River distance (km)
AC
CE
PT
ED
Figure 6. River1d simulation of ice jam release and formation of sell-sustaining jave (kms 0 to 15: open water; kms 15 to 24: ice jam; kms 24 to 45: intact “ice cover”, simulated using elevated value of nb. The coloured lines represent instantaneous water level profiles up to 46 minutes following release. Front celerity ~5.1 m/s; front height ~ 1.7 m above pre-jave water depth downstream of the jam.
18
ACCEPTED MANUSCRIPT 7 0 min.-jam 5 min.
9 min.
6
13 min. 26 min.
5
T
33 min.
IP
40 min. 46 min.
4
CR
Depth (m)
19 min.
US
3
0
5
10
AN
2 15
20
25
30
35
40
45
M
River distance (km)
CE
PT
ED
Figure 7. Instantaneous longitudinal variations of water depth at different times following release of the jam depicted in Figure 6 (kms 0 to 8: open water segment unaffected by jam backwater; kms 8 to 15: open water segment affected by jam backwater; jam extends from km 15 to km 24; intact “ice cover” beyond km 24).
The front and the crest of the wave remain essentially unchanged after ~30 minutes from the
AC
release of the jam, a much shorter development time than was found in the CRISSP runs (many hours). The constancy of the wave front and crest shape implies that all points along the front advance at the same celerity and there is no attenuation. This is an important difference from ordinary javes, which attenuate and spread over time, owing to increasing celerity along the waveform, from the trailing to the leading edge (Beltaos and Burrell 2005a). Unlike the case of ice melt, where the front of the wave increases in length and propagates under the intact ice
19
ACCEPTED MANUSCRIPT cover, the SSW wave has a sharp front and only a portion of this front is moving under intact ice. The flow depth at which the ice cover is dislodged (Ybr) is not much less than the sustained peak depth (Yp); ratio Ybr/Yp = 0.89.
T
A simplified sketch of the SSW is shown in Fig. 8. If the ice breaking “trigger” is a jave, there is
IP
no net inflow within the control volume (Q = Qu at both ends). The same applies to the CRISSP
CR
runs after the time when the input hydrograph returns to the base flow. In such instances,
US
continuity requires that
(10)
M
AN
Cbr (Yp Yi ) Ctail (Yp Yu )
in which Ctail = celerity of the tail end of the wave and Yi = depth of flow under the undisturbed
ED
ice cover downstream of the wave front, not to be confused with the water depth Yi (Fig. 3)
PT
which is measured from the water surface to the bed. Clearly, Ctail < Cbr. Therefore the uniformflow crest AB will lengthen as time goes on. The few wave shapes shown in Fig. 7 only hint at
AC
CE
this feature, which is more conspicuous in the extensive CRISSP runs (Jasek et al 2005).
20
US
CR
IP
T
ACCEPTED MANUSCRIPT
AN
Figure 8. Schematic illustration of wave generated by sustained ice dislodgment. The wave front is depicted as a straight line segment for simplicity.
M
The preceding results suggest further that the properties of the SSW do not depend on whether it is triggered by rising runoff or by jam release. Once ice breaking is initiated and the wave has
ED
travelled for a sufficiently long time, its front takes on a constant shape and moves at constant
PT
rate. These parameters are primarily dictated by the threshold criterion, be it a threshold
CE
discharge or a threshold shear stress.
AC
Figure 9, which derives from CRISSP results provided in Table 1 of Jasek et al (2005), shows how the celerity of the wave front (Cbr) decreases as the threshold discharge increases. Extremely rapid ice breaking is indicated when this threshold is just above the base flow while a strong negative influence of the bed slope is also evident. One may note further that Cbr is always greater than the celerity of the kinematic wave (~1.5Ui) and generally differs from the celerity of the gravity wave, Cg, which can be expressed in dimensionless form as:
21
ACCEPTED MANUSCRIPT Cg Ui
gYi Ui
1
(11)
30 CRISSP runs - bed slope = 0.0003 CRISSP runs - bed slope = 0.00005
T
25
CR
Cg/Ui - slope = 0.00005
15
US
Cbr / Ui
20
IP
Ice jam release - bed slope = 0.0005
AN
10
Cg/Ui - slope = 0.0003 0 1.4
1.6
ED
1.2
1.8
2
2.2
2.4
Qbr / Qo
PT
1
M
5
CE
Figure 9. Variation of ice breaking celerity with threshold discharge, both normalized with respective base quantities (the CRISSP runs are described in Jasek et al 2005; Ui was calculated by the writer via the Manning resistance equation)
AC
For So = 0.0003 and 0.00005, Cg/Ui amounts to 8.06 and 16.80, respectively. As shown in Appendix B, it is not physically correct for Cbr to exceed Cg and this is consistent with previous work on javes and open-water waves (Ponce and Simons, 1977; Ferrick and Goodman 1998; Beltaos 2013b). On the other hand, Fig. 9 suggests that the upper-most data point of each set exceeded the corresponding value of Cg/Ui. This could be the result of numerical-discretization, as discussed in Appendix B.
22
ACCEPTED MANUSCRIPT Figure 10 shows that the sustained peak flow, Qp, is primarily determined by the threshold discharge, as well as by the base flow and the channel slope. Since both Qp and Qbr exceed Qo, excess ratios are plotted in Fig. 10 for convenience. Within their respective ranges, the data points for each slope are well described by the corresponding linear fits. For the tested bed
T
slopes, the flow contributed by storage release is seen to range from ~70% to ~200% of the base
IP
flow. Moreover, it seems likely that higher percentages can be attained at higher bed slopes.
CR
These figures underscore the dynamic nature of the SSW relative to melt-induced waves, which
2.5 CRISSP runs- bed slope = 0.0003
2.0
AN
CRISSP runs - bed slope = 0.00005
Ice jam release - bed slope = 0.0005
M
1.5
ED
Q P / Qo - 1
equality line
CE
PT
1.0
AC
0.5
0.0
US
generate excess flows under ~40% of the base flow (Fig. 5).
0
0.2
0.4
0.6
0.8
1
1.2
1.4
Qbr / Qo - 1 Figure 10. Variation of sustained peak flow with threshold discharge; both normalized with the base discharge Qo and expressed as “excess” values by subtracting 1 from each ratio. The CRISSP runs are described in Jasek et al (2005).
23
ACCEPTED MANUSCRIPT Since the open-water flow in the sustained-peak reach is uniform, one can write:
Yp
2
nbQ p ncQo
3/ 5
(12)
T
Yi
2 / 5
IP
To derive Eq. 12, it was assumed that the ice blocks being carried by the flow in the sustained-
CR
peak reach do not impede water motion. Jasek et al (2005) considered possible impedance by simply increasing the local value of nb to 0.04 and 0.05 but the resulting wave profiles were
AN
US
deemed physically implausible. Consequently, the results of these runs are not considered herein.
ED
(13)
PT
Q p Qo Cbr (Y p Yi )W
M
Continuity (Eq. 3) requires that
CE
Combining Eqs 12 and 13 results in:
(14)
AC
(Q p / Qo ) 1 Cbr 2 / 5 U i 2 ( nbQ p / ncQo )3/ 5 1
which suggests that the normalized ice breaking celerity is uniquely defined by the normalized sustained peak flow and by the ratio of the bed and composite-flow Manning coefficients nb and nc.
24
ACCEPTED MANUSCRIPT The relationship of Eq. 14 is validated in Fig. 11 where the various data points collapse onto a single curve. If the line shown in Fig. 11 were extended beyond the right-hand margin, it would be seen to attain a minimum of Cbr/Ui = 4.03 and then begin to increase slowly. During the prebreakup phase of the ice cover, the velocity Ui is often in the range 0.7 to 1.0 m/s. Consequently,
T
SSWs are expected to advance at speeds of no less than ~2.8 m/s (for nb/nc=1). The latter value
IP
amounts to ~240 km/day, which is much greater than what can be expected from ice melt, and
US
30
CR
again illustrates the dynamic nature of the SSW.
CRISSP runs - bed slope = 0.0003
25
AN
CRISSP runs - bed slope = 0.00005 Ice jam release - bed slope = 0.0005
M ED
Predicted from continuity equation
15
PT
Cbr /Ui
20
CE
10
0
AC
5
1.5
2
2.5
3
3.5
Qp /Qo
Figure 11. Normalized ice breaking celerity versus normalized peak wave discharge for nb = nc (the CRISSP runs are described in Jasek et al 2005; Ui was calculated by the writer via the Manning resistance equation). The line is not physically meaningful where Cbr > Cg
25
ACCEPTED MANUSCRIPT Equation 14 indicates further that the ratio Qp/Qo must exceed 22/3(nc/nb)5/3. Below this limit, Cbr would be calculated as negative while at the limit it becomes infinite. For nb = nc, this limit takes on the value 22/3 or 1.587. This means that Qp/Qo – 1 has to exceed 0.587 which is seen to be the case for all data points in Fig. 10 and implies that SSW storage release contributes no less than
IP
T
~60% of the base flow.
CR
The ratio nb/nc has a strong effect on Cbr, as illustrated in Fig. 12. Increasing values of nb/nc result in lower values of Cbr/Uo for any given ratio Qp/Qo. The value of nb/nc would typically be
US
expected to be equal to 1 (e.g. nb = ni = nc = 0.03) or more than 1 (e.g. nb=0.03, ni = 0.02, nc =
AN
0.025). Cases in which nb/nc is less than 1 require that ni be greater than nb, which is not likely for sheet ice covers, but could occur where the ice cover has a very rough underside owing to
AC
CE
PT
ED
M
slush deposition during freezeup or possibly jamming of ice blocks during mid-winter breakup.
26
ACCEPTED MANUSCRIPT
20
nb/nc = 0.8
15
0.9
T IP
10 1.1
CR
Cbr /Ui
1.0
1.2
AN
US
5
0
2
3
4
5
Qp /Qo
ED
M
1
PT
Figure 12. Normalized ice breaking celerity versus normalized peak wave discharge for different values of the ratio nb /nc. The lines are not physically meaningful where Cbr > Cg
CE
As noted earlier, the jave-trigger simulation indicated that the ice cover was dislodged when the under-ice flow depth, Ybr, was not much less than the sustained peak depth Yp (ratio = 0.89).
AC
Jasek et al (2005) did not provide values of Yp and Ybr in the tabulation of their results but Yp was calculated by the writer assuming open-water flow, as indicated previously; while Ybr was determined from the equation of continuity (Eq. 3), which implies that
Ybr Yi (Qbr Qo ) / CbrW
(15)
27
ACCEPTED MANUSCRIPT Calculated values of Ybr were always slightly less than those of Yp while the ratio Ybr/Yp was approximately constant for a given slope, regardless of the value of Qbr (~0.97 and 0.90 for So = 0.00005 and 0.0003, respectively).
T
While the preceding quantitative results may be modified by different values of the Manning
IP
coefficients nb and ni, a reasonably well-defined picture of the SSW has emerged. The relatively
CR
sharp wave front enters the ice-covered region and lifts the ice cover, which does not break until the flow depth is a little lower than the sustained-peak depth. Upstream of the ice breaking
US
location, a transitional reach leads to the prolonged wave crest. Here, broken ice is carried along
AN
with the water but the Manning resistance equation is neither fully of the ice-covered nor of the open-water type (Appendix B). This effect likely dictates the increase in flow and depth over and
ED
M
above the respective ice-breaking values upstream of the “breaking” point.
The SSW has certain similarities with the well-known monoclinal wave (Ferrick and Goodman
PT
1998) but is unique in that its water surface profile combines open-water and ice-covered reaches
CE
which are joined by a hybrid reach. Appendix B contains a detailed discussion of this aspect and presents methodology for analytical determination of the profile of the wave front based on
AC
integration of the Saint-Venant equations of motion.
5. NATURAL STREAM CONDITIONS: PEACE RIVER CASE STUDY Since 1973, Alberta Environment (AENV) and BC Hydro have carried out detailed observations of ice formation and breakup throughout each ice season between the Bennett Dam and Vermilion Chutes (Fig. 13). The Chutes are located ~900 km below the Bennett Dam and some
28
ACCEPTED MANUSCRIPT 70 km downstream of Fort Vermilion. On rare occasions, observations have extended beyond Vermilion Chutes, all the way to the Mouth of Peace (MOP for short), which is also the head of the Slave River, and situated 1243 km from the dam. River slope generally decreases with distance from the dam. Jasek and Pryse-Phillips (2015) indicate slopes of 0.0003 between ~kms
T
300 and 500 (approx. Dunvegan to Sunny Valley), and 0.00005 between ~kms 700 and 900
AC
CE
PT
ED
M
AN
US
CR
IP
(approx. Tompkins Landing to Vermilion Chutes).
Figure 13. Map of Peace River main stem in BC and Alberta showing key locations and river kilometres from the W.A.C. Bennett Dam. From Jasek and Pryse-Phillips (2015), with changes.
The upstream advance of the ice front in the fall and winter and its recession in the spring are summarized graphically in a date-distance plot of the kind illustrated in Fig. 14. Ice front recession usually begins in late February or early March, slowly at first but accelerating as time
29
ACCEPTED MANUSCRIPT goes on. The rate of downstream advance of the front is typically less than ~ 50 km/day, while
ED
M
AN
US
CR
IP
T
there are 4 occasions where it exceeds 100 km/day over considerable distances.
PT
Figure 14. Observed locations of Peace River ice front downstream of Bennett Dam. The dashed green line has a slope of 50 km/day. From Alberta Environment archived material, with changes (accessed Dec. 19, 2016). http://www.environment.alberta.ca/forecasting/RiverIce/pubs/2013-
CE
2014_Peace_River_Ice_Obs_Rpr_No46.pdf
AC
With typical values of Q (3000 m3/s) and W (600 m) for the lower Peace River (slope ~0.00005), the rate of advance from Eq. 1 has been calculated as 8.6 km/day per degree C of water temperature. Observed rates of up to 50 km/day correspond to water temperatures of up to ~6°C and can generate storage-release flow enhancements of up to ~50% (Fig. 5). The increase in the speed of the ice front with the passage of time is consistent with the lengthening of the openwater reach upstream of the front. It thus appears that ice melt dominates the process of ice front recession. Note, however, that such rates of advance can also be experienced by intermittent 30
ACCEPTED MANUSCRIPT jamming and release; javes propagate much faster than tens of km per day but jamming may last long enough to significantly reduce the overall-average rate of advance. Not shown in Fig. 14 is an extremely rapid advance between Vermilion Chutes and MOP that occurred in 2014 and is
T
discussed in detail next.
IP
Unlike in the case of a prismatic channel, the resistance of the ice cover to dislodgement varies
CR
along a river owing to variable bathymetry, slope, and channel planform as well as variable ice thickness and strength (Beltaos 2013b). If a SSW encounters a highly resistant ice cover
US
segment, it may be unable to dislodge it and begin to attenuate as it moves under it. If this
AN
attenuation is small (e.g. because the resistant segment is short), the wave can resume its ice breaking action and again become self-sustaining. At the same time, the running ice that trails
M
the breaking front will be arrested at the edge of the resistant ice cover and a jam will be
ED
initiated. As backwater builds behind the jam, it is possible that the jam will release soon afterwards and produce a new jave that follows on the heels of the first one. Alternatively, the
PT
jam may remain in place for a long enough time for the first jave to be arrested somewhere
CE
downstream and form another jam. One would then encounter two jams at the same time, separated by a considerable stretch of open water. This configuration occurred in Peace River
AC
near Fort Vermilion in late April of 2014.
According to publically available AENV observation reports, ice breaking commenced near the Town of Peace River (TPR) on April 22. In the following days, the breakup progressed in a sequence of jams and releases, with major ice runs noticed on the 25th and 27th. In the morning of
31
ACCEPTED MANUSCRIPT April 28, ice conditions were as follows (archive accessed Dec. 19, 2016): (http://www.environment.alberta.ca/forecasting/RiverIce/PeaceRiverArchive1314.html)
Open water with brash ice above km 806; ice jam, km 806 to 830; open water, km 830 to 872;
T
ice jam, km 872 to 900 (approx. Vermilion Chutes); relatively intact ice cover below km 900.
IP
(The observation reports and associated diagrams indicate that Fort Vermilion is located at km
CR
832, which the writer assumes to be the site of the local WSC hydrometric gauge). By the afternoon of April 29, this configuration had changed to: ice jam, km 813 to 830; open water, km
US
830 to 878.5; ice jam, km 878.5 to 900; open water, km 900 to 915; brash ice, km 915 to 916;
AN
intact ice cover, km 916 to 935 and most likely beyond.
M
According to M. Jasek of BC Hydro (Pers. Comm. 2014) the jam upstream of Fort Vermilion
ED
released at 1815 h (Mountain Daylight Saving Time or MDT) on April 30. The record of the local gauge indicates that the water level began to decrease shortly after 1300 h (Mountain
PT
Standard Time or MST) and began to rise at 1716 h MST. [In what follows, time will be
CE
expressed as standard time and the MST designation will be implied, unless otherwise specified]. The decrease in stage (~0.6 m) was most probably caused by the negative wave generated by
AC
earlier release of the downstream jam. The travel time of this negative wave is estimated as 3.5 hours, based on the theoretical analysis of Henderson and Gerard (1981). Consequently, the release of the downstream jam would have occurred at about 1000 h. This estimate could be in error by an hour or so because relevant hydraulic properties of the reach between Fort Vermilion and Vermilion Chutes are unknown. On May 1, the river was observed by Mr. Jasek from the air and an ice run of variable surface concentration was seen all the way to the Mouth of Peace (M.
32
ACCEPTED MANUSCRIPT Jasek, Pers. Comm. 2014). This run was arrested by intact ice cover ~10 km into the Slave River; the resulting ice jam remained in place until May 9 and produced large-scale flooding of the Peace-Athabasca Delta. This kind of flooding is essential to the ecological maintenance of the higher-elevation lakes and ponds of the delta (perched basins) and can only be produced by ice
IP
T
jamming (Prowse and Conly 1998; Beltaos et al. 2006).
CR
The hydrometric gauge record at Peace Point for May 1 (Fig. 15) comprises data sets of different reliability, owing to ice effects on the primary and secondary orifice lines. The primary orifice
US
line was moved by ice on April 29, while the secondary line was likely pinched by ice in early
AN
morning of May 1. The primary record was corrected by WSC staff by applying a constant offset starting at 0825 h of May 1; earlier data were deemed unusable and thence not included in the
M
corrected record. The writer applied the same offset value as the one used by WSC to the earlier
ED
gauge readings and obtained the line designated to be of “low confidence”. The gauge heights indicated by this line are of unknown accuracy, but the pattern may be helpful in piecing together
PT
details of the breakup of the local ice cover. The flat line before ~0300 h May 1 indicates that
CE
the orifice had been tossed above the water level and was reading atmospheric pressure. At ~0320 h, it appears that a brief wave (wave 1) passed by the gauge site, re-submerging the
AC
orifice, which again began reading some water pressure, even though its actual elevation remained unknown.
33
ACCEPTED MANUSCRIPT 12 wave 2
wave 1
IP
T
10
orifice above water level
corrected primary water level; high confidence (WSC) corrected primary water level; low confidence secondary water level
US
8
CR
9
AN
Gauge height (m)
11
0:00
6:00
M
7
12:00
18:00
0:00
ED
May 1, 2014, MST
CE
PT
Figure 15. Hydrometric gauge record at Peace Point during the 2014 ice breakup.
Hourly WSC time-lapse images at Peace Point indicate that the sheet ice cover may have shifted
AC
during the night of Apr 30 to May 1, but remained stationary until at least 0701 h of May 1; at 0801 h rubble appears in the image and the ice sheet was likely in motion because its surface looks different from that of the previous photo (Fig. 16). The presence of intact sheet ice at 0801 h suggests that the cover had been dislodged recently; else it would have already been reduced to smaller slabs and blocks. Transport of ice blocks between the sheet ice cover and the near bank (lower photo, Fig. 16) might have caused the large water level fluctuations indicated by the lowconfidence line in Fig. 15 after 0700 h of May 1. 34
AC
CE
PT
ED
M
AN
US
CR
IP
T
ACCEPTED MANUSCRIPT
Fig. 16. Peace River ice conditions at the Peace Point gauge site, morning of May 1, 2014. Flow is from right to left. Images courtesy of Water Survey of Canada.
35
ACCEPTED MANUSCRIPT Full-width rubble appeared in the 0901 h image, while subsequent images suggested that the broken ice might have been running continuously at varying surface concentration, until 0400 h of May 2. Visual comparisons of ice surface elevations relative to the far river bank between successive images have indicated that the water level changed little during 0401 and 0701 h of
T
May 1, while there was a conspicuous rise between 0701 and 0801 h. The relative steadiness of
CR
IP
the water level prior to 0701 h in the photos is consistent with the low-confidence line in Fig. 15.
Since the photos indicate that ice was mobilized between 0701 and 0801 h, selection of 0730 h
US
would not be in error by more than ~30 min, even though it is likely that the actual time of
AN
mobilization was closer to 0800 h than to 0700. The estimated travel time between kms 900 and 1135 (Peace Point), would then be 21.5 h, indicating an average ice-breaking rate of 235/21.5/3.6
M
= 3.0 m/s. Available hydrometric data at Peace Point (Beltaos 2011) indicate a base velocity of
ED
0.75 m/s under the intact ice cover, rendering the ratio Cbr/Ui equal to 4.0. This value is compatible with the data points and line shown in Fig. 11; therefore it could have resulted from a
CE
PT
SSW.
Known open-water and ice-covered flow hydraulics in the vicinity of Peace Point (Beltaos 2011)
AC
were utilized next to calculate the base-flow as 3450 m3/s (Table 1), assuming a base gauge height of 9.0 m. The ice thickness on May 1 was estimated as 0.82 m with a keel of 0.75 m (based on a measured value of 0.92 m on April 8 and typical ablation rates through April). If the breakup at Peace Point was caused by a SSW, with peak stage of 11.70 m, the flow Qp can be estimated as follows: (a) subtract 0.5 m to account for rubble, assumed to be 1 m thick with a porosity of 0.5, to find “net” stage = 11.2 m; (b) calculate the difference from the stage at the
36
ACCEPTED MANUSCRIPT bottom of the ice cover as 11.2-(9-0.75) = 2.95 m; (c) multiply by Cbr (3.0 m/s) and by the local width of 680 m, to find the excess flow as 6020 m3/s (via continuity), for a total Qp of 6020+3450= 9470 m3/s (Qp/Qo = 2.75). At the beginning of May, the Manning coefficient of the ice cover is ~0.024, same as that of the river bed (Beltaos 2011); therefore nc ~ 0.024 and nb/nc ~
T
1, equal to the ratio used in the CRISSP runs. Therefore, the methodology outlined in Appendix
IP
B can be applied to determine the value of Qbr that would produce a SSW with peak flow of
CR
9470 m3/s. This operation resulted in Qbr = 8750 m3/s (Qbr/Qo = 2.54), while the predicted value
US
of Cbr is 3.4 m/s, in fair agreement with the observed value of 3.0 m/s.
AC
CE
PT
ED
M
AN
Table 1. SSW analysis of ice breaking at Peace Point Variable Value Water surface slope, So 0.000064 Manning coefficients, ni, nb, nc 0.024, 0.024, 0.024 Base flow, Qo 3450 m3/s 0.75 m/s, 6.76 m Base velocity, Ui, and flow depth, Yi Base gauge height 9.0 m Bottom of ice gauge height 8.25 m Peak and effective(1) gauge heights 11.7 m, 11.2 m Sustained peak flow, Qp 9470 m3/s Ice-breaking flow, Qbr 8790 m3/s Ice-breaking celerity, Cbr (observed) 3.0 m/s Ice-breaking celerity, Cbr (calculated) 3.4 m/s Expected duration of crest at Peace Point 9 hours (1) gauge height that would prevail if there were no moving ice
It is possible, however, that the ice run may have stalled at one or more locations between Vermilion Chutes and Peace Point, forming short-lived jams, which upon release would have generated large javes, propagating much faster than 3.0 m/s. Jave speeds of 8 to 9 m/s have been estimated on two separate occasions for Peace Point (Beltaos 2007), so that there would have been ample time for stoppage between the mornings of April 30 and May 1.
37
ACCEPTED MANUSCRIPT
To further explore this issue, the front profile of the SSW wave was calculated according to the methodology developed in Appendix B. Figure 17 shows that this front comprises a gradual rise that would have commenced very early in the morning of May 1; this is not consistent with the
T
available, albeit low-confidence, gauge data and with the time-lapse imagery (a sharp rise did not
IP
start until after 0700 h). It is more likely that the ice cover was dislodged by a large wave (wave
CR
2) of relatively high celerity, produced by the release of an upstream ice jam, possibly initiated at a sharp bend above Boyer Rapids (km ~1117), a known jamming site. If this is indeed what
US
occurred, wave 1 could have been the leading part of the original jave, which was “clipped” as
AN
the jam formed and water again went into storage. The breakup front arrived at the MOP at about 1310 h (according to time-lapse imagery, kindly provided to the writer by Dr. Daniel Peters of
M
Environment and Climate Change Canada). With the estimated time of ice mobilization at Peace
ED
Point (0730 h), the average ice-breaking celerity over the intervening 108 km works out to be 5.3 m/s, almost twice the calculated average value over the previous 235 km. For the entire travel
PT
reach (343 km), the average ice-breaking celerity is approximately equal to 343/(27.173.6) or
CE
3.5 m/s.
AC
The preceding interpretation is corroborated by the non-sustained nature of the water level peak (Fig. 17). Had the peak been produced by a SSW, it would have lasted much longer: using Eq. 10, the celerity of the tail of the SSW is calculated as 1.7 m/s, and thence the duration of the sustained crest is estimated to have been ~9 h, placing its end at about 1730 h on May 1, well beyond the right-hand border of Fig. 17. The observed water level variation, on the other hand, continues to decrease after the peak, attaining a value of 10.25 m by 1730 h.
38
ACCEPTED MANUSCRIPT
12 WSC corrected water level WSC secondary water level
IP
T
Calculated from SSW theory
CR
10
US
Gauge height (m)
11
8 3:00
6:00
9:00
12:00
ED
0:00
M
AN
9
May 1, 2014, MST
CE
PT
Figure 17. Calculated rising limb of SSW at Peace Point. Spatial wave configuration was transformed to temporal variation using the ice-breaking celerity Cbr; the starting time was selected so as to match the time of the peak. Predicted ice mobilization time and gauge height: 0720 h and 11.3 m, respectively.
AC
The preceding considerations suggest that the SSW theory does not match the details of how ice breakup was initiated at Peace Point and possibly elsewhere along the river. Nevertheless, the theory provides a helpful framework for examining breakup events in an “average” sense. Since the resistance to ice breaking varies along the river, the breaking front will speed up or slow down depending on local resistance characteristics; occasional brief jamming could still fit into this framework, so long as the average celerity of ice mobilization is not less than the minimum Cbr value indicated in Figs. 9 and 11. Consequently, the values of Qbr and Qp that are listed in 39
ACCEPTED MANUSCRIPT Table 1 should be viewed as magnitudes that would occur if breakup had been advancing at a constant celerity, equal to the average observed value. This may not reflect what actually happens at any particular site (such as Peace Point), but does provide a “feel” for average
T
conditions over a long travel reach.
IP
6. EFFECTS OF CHANGING THRESHOLD CONDITIONS ON STORAGE RELEASE
CR
It has been shown in the preceding sections that the ice-breaking threshold parameters have a controlling influence on storage release. Therefore, it is of interest to consider how climate-
US
induced changes to river ice regimes may affect storage release flows. For example, climate
AN
warming may bring about a significant reduction in the thickness of the ice cover of a particular river. This change would reduce the resistance to dislodgment and thence the value of the
ED
M
threshold flow Qbr.
Let QbrP and QbrF be the present and future values of Qbr, such that
QbrP > QbrF. Though it is
PT
likely that the spring flow hydrograph will also change, such change is ignored in order to focus
CE
on the threshold flow effects, assuming that other factors are equal. If the highest carrier flow of the spring runoff hydrograph (Qomax) is lower than both QbrP and QbrF, a thermal breakup event
AC
will result under both present and future conditions; the relatively small flow enhancement by storage release would be the same.
A second possibility is that of a mechanical event occurring under both present and future conditions (QbrF < QbrP < Qomax). In this case, Fig. 10 suggests that the flow enhancement due to storage release will be smaller under future than under current conditions. The near-linear
40
ACCEPTED MANUSCRIPT variations implied by the plotted data points indicate further that the difference will be equal to the slope of the applicable line, multiplied by (QbrP- QbrF). In the case of the flatter channel in Fig. 10, the slope of the line defined by the data points is ~0.55; if QbrP- QbrF = 1000 m3/s, the flow enhancement under future conditions will be 550 m3/s less than it is under present
T
conditions. This calculation applies so long as Qbr/Qo is within the range defined by the data
IP
points in Fig. 10. It is reasonable to expect, and has been shown in Appendix B, that the linear
CR
fits in Fig. 10 will eventually bend and become parallel to the equality line as Q br/Qo increases (Fig. B3). In the unlikely case that both QbrP and QbrF lie within the latter range, the difference in
AN
US
peak flows will be equal to QbrP- QbrF or 1000 m3/s in the present example.
A mixed situation arises when the carrier flow is intermediate between present and future
M
thresholds, i.e. QbrF < Qomax < QbrP. In that case, the breakup event would be thermal under
ED
present conditions and mechanical under future conditions. Since only mechanical breakup events can lead to major ice jams, the frequency of ice-jam flooding could increase in the future
PT
regime, despite the reduced ice thickness. Moreover, storage release under present conditions
CE
would be controlled by melt, rather than by ice breaking. Consequently, flow enhancement under future conditions (QF) would be considerably more than the “present” value (QP), as can be
AC
seen by comparing Figs. 5 and 10. Table 2 summarizes the effects of a reduced threshold flow on storage-release flow enhancement. Table 2. Effects of reduced ice-breaking threshold flow on breakup type and storage release Case Magnitude of carrier flow Breakup type Effect on flow number relative to ice-breaking enhancement due to Present Future thresholds storage release conditions conditions 1 Qomax < QbrF < QbrP thermal thermal no effect 2 QbrF < Qomax < QbrP thermal mechanical positive, large 3 QbrF < QbrP < Qomax mechanical mechanical negative, small ( QbrP- QbrF) 41
ACCEPTED MANUSCRIPT
Similar reasoning can be applied to situations where the threshold flow will increase, i.e. QbrP < QbrF. This could happen, for instance, if the fall season becomes wetter in the future and thence results in higher freezeup levels than under the present climate. This would increase
T
resistance to ice dislodgment and breaking in the spring and result in a higher value of Qbr. The
CR
IP
effects of this kind of change are summarized in Table 3.
M
AN
US
Table 3. Effects of increased ice-breaking threshold flow on breakup type and storage release Case Magnitude of carrier flow Breakup type Effect on flow number relative to ice-breaking enhancement due to Present Future thresholds storage release conditions conditions 1 Qomax < QbrP < QbrF thermal thermal no effect P max F 2 Qbr < Qo < Qbr mechanical thermal negative, large P F max 3 Qbr < Qbr < Qo mechanical mechanical positive, small ( QbrF- QbrP)
ED
Notable impacts of increased threshold flows are: (a) reduced frequency of ice jams and related
PT
flooding due to less frequent occurrence of mechanical breakup events; and (b) large reduction in flow enhancement for case 2 and small increase for case 3. This increase would typically be less
CE
than the difference in threshold flows, but could approach (QbrF - QbrU) if both QbrF and QbrU are
AC
larger than about 2Qo (Fig. B3).
An effect that is difficult to quantify at present pertains to case 3 in Table 3. Since QbrF > QbrP, ice breaking will be occurring at a higher water level and at lower speed under future than under present conditions (Fig. 9). Therefore, there will be increased attrition of the volume of brokenup ice that forms the ice run behind the ice-breaking front. Such attrition is primarily caused by stranding of ice blocks on islands and low banks, as well as by melt that may be occurring in 42
ACCEPTED MANUSCRIPT transit (Prowse 1986). Increased attrition will result in reduced ice volumes that are available to form ice jams when ice runs are eventually arrested. In turn, this would reduce the potential of the jam for large-scale flooding.
T
7. DISCUSSION
IP
Melting and sustained ice breaking have been identified as two processes that can enhance the
CR
flow during the breakup period. River slope is an important factor that has opposing effects, depending on the process. In the case of melting (Fig. 5), typical river slopes (e.g. 0.0005) seem
US
to produce just a slight flow enhancement, which becomes much more substantial in flat rivers
AN
(e.g. slope = 0.00005). In the case of ice breaking, higher slopes produce higher flow enhancement than do lower slopes (Fig. 10). Moreover, the flow enhancement caused by ice
ED
M
breaking is typically much larger than what is generated by ice melt.
A striking result of the CRISSP runs is that ice-breaking can produce SSWs with sustained peak
PT
flows well in excess of the peak of the incoming carrier flow hydrograph, up to ~60% higher
CE
according to Jasek et al (2005). It is not known whether peak amplification can occur when ice breaking is initiated by a jave, as is usually the case in flat rivers like the lower Peace or the
AC
lower Mackenzie (Beltaos and Carter 2009, Beltaos 2013a). The simulation with River1d indicated that the peak flow was reduced when the jave morphed into a SSW. Javes typically encounter intact ice cover as soon as they are generated because jams are held in place by intact ice. At that time the jave front is almost as high and steep as the toe region of the parent jam (Fig. 7). If the jave does break the ice cover, a SSW may eventually form but it seems unlikely
43
ACCEPTED MANUSCRIPT that its front can be as high and steep as that of the jam toe. This conjecture is, of course, subject to testing with further numerical modelling.
As noted in Section 4, Jasek et al (2005) found that enhancing the value of nb upstream of the
T
wave front in order to simulate possible impedance to flow by the moving broken ice led to
IP
implausible results. Consequently, this effect has not been considered herein. Since front celerity
CR
exceeds water velocity, as can be easily ascertained using Eq. 13, ice congestion is unlikely to
US
occur. Moreover, the broken ice is floating at a higher elevation than that of the intact ice cover and therefore is less confined by the river banks; any resistance that may develop near the river
AN
banks is likely to be minimal, relative to bed friction, which is applied throughout the channel width. On the other hand, the thickened ice pack from the released jam moves at a lower
M
elevation that than of the jam; it could conceivably generate significant side resistance before
ED
dispersion reduces its thickness to that of a single ice block. However, the moving pack trails
PT
well behind the wave front and does not affect front properties. These considerations are consistent with the findings of She and Hicks (2006), which showed that moving-ice impedance
AC
CE
is of secondary importance in jave propagation, especially near the wave front.
An important point that should be stressed here is that flow enhancement due to storage release is a highly transient phenomenon that ceases when the recession of the edge of the ice cover is discontinued. For instance, ice breaking will cease when a sufficiently resistant ice cover is encountered by the rising flow and a jam begins to form. Downstream of the jamming point (toe of jam), the wave will continue to advance but keep attenuating unless ice breaking resumes. Upstream of the toe, water levels will rise rapidly, creating new storage and thence reducing 44
ACCEPTED MANUSCRIPT downstream flows. When all of the moving ice comes to a stop, upstream and downstream flows will return to the value of the carrier river flow, which therefore governs the potential of the jam to cause flooding over a sustained period of, say, a day or more.
T
It follows that the main effect of storage release is in facilitating the occurrence of a mechanical
IP
event when the carrier flow in the river is lower than the threshold value. In the case of melt, this
CR
eventuality seems very unlikely in steeper channels but conceivable in flat ones, for which meltinduced flow enhancement can be sizeable (Fig. 5). This is particularly probable where the
US
threshold flow is not much greater than the carrier flow. In the case of ice breaking, this
AN
reasoning may appear circular for prismatic channels: if the carrier flow cannot initiate breakup at the beginning of the ice covered reach, how can it be enhanced at any location farther
M
downstream? In natural streams, however, one may envisage circumstances where the carrier
ED
flow can actually dislodge the ice cover in a low-resistance section of the river. The resulting enhanced flow could then exceed the threshold farther downstream, where the resistance to
CE
PT
dislodgment may be higher.
Climatic warming could reduce or increase the values of ice-breaking threshold flows, depending
AC
on basin-specific changes, such as thinner ice covers or higher freezeup levels. In turn, such changes can have significant effects on storage release. Using the analytical methodology developed in previous sections and in Appendix B, it was possible to assess the effects of changes in threshold flows on storage release contributions to the flow hydrograph. Such changes might also result from regulation; assessment of their impacts could also be made using the summaries presented in Tables 2 and 3.
45
ACCEPTED MANUSCRIPT
8. SUMMARY AND CONCLUSIONS Two mechanisms of storage release, ablation and ice breaking, have been examined in details and quantified. Melting of the ice cover can significantly enhance the flow, especially in
T
relatively flat channels. Excess flows may amount as much as 40% for water temperatures up to
IP
5°C. Storage release by ice breaking is much more dynamic and results in considerably greater
CR
excess flows than ice melt. Numerical modelling and analysis for prismatic channels suggest that ice breaking, once initiated, morphs into a SSW with constant front shape and celerity. The
US
properties of the SSW do not appear to depend on the “trigger” flow hydrograph or jave
AN
properties but are largely determined by the threshold discharge, as well as base flow, ice, and channel properties. Ice breaking celerities can be large if the threshold discharge is very close to
M
the base flow, but decrease to a minimum value of ~ 4 x (under-ice base velocity) as the
ED
threshold-to-base flow ratio increases. Application of the analytical and numerical results to the lower Peace River indicated that ice melt could explain typical observed rates of ice recession;
PT
the SSW theory was only partially applicable to the very rapid, and rare, ice clearance that was
CE
observed in the spring of 2014. Ice jamming at high ice resistance sites that are invariably encountered in natural streams can disrupt the orderly advance of the wave. Nevertheless, it was
AC
concluded that the SSW concept can be useful in terms of the overall average advance of highly dynamic breakup events. The effects of changing ice regimes, e.g. due to climate change or regulation on storage release, can be assessed in terms of changes in threshold ice-breaking flows by simple reasoning that was illustrated in two examples. It has been stressed that flow enhancement by storage release is highly transient and its main role is in facilitating ice breaking farther downstream.
46
ACCEPTED MANUSCRIPT
9. ACKNOWLEDGMENTS Support by Environment and Climate Change Canada to carry out this work is gratefully acknowledged. I would also like to thank: Angus Pippy and his colleagues at the WSC
T
Yellowknife office for supplying archived, but unpublished, hydrometric data and interpretation;
IP
Professor Faye Hicks (ret. University of Alberta) for making available the model River1d;
CR
Martin Jasek of BC Hydro for sharing observational information regarding the 2014 Peace River breakup; Nadia Kovachis Watson of Alberta Environment for sharing historical data on ice cover
US
recession rates in Peace River; and my colleague Dr. Daniel Peters for sharing 2014 time-lapse
AN
imagery at the mouth of Peace River .
ED
M
APPENDIX A. USING INCREASED BED MANNING COEFFICIENT TO SIMULATE ICECOVERED FLOW In the present applications of the River1d model, the open-water option has been used for
PT
simplicity, while the presence of a stationary ice cover of constant thickness and roughness has been simulated by locally enhancing the value of the bed Manning coefficient. Using the
CE
Belokon-Sabaneev formula (Eq. 6), it can be shown that this assumption is equivalent, in terms
AC
of flow depth and mean velocity, to having an actual ice cover of which the roughness coefficient (ni) is given by:
ni (nbE 3/ 2 nb3/ 2 )2 / 3
(A1)
with nb = actual roughness coefficient of the river bed; and nbE = enhanced bed roughness value used to simulate the “ice covered” reach. This simplification facilitates the simulation of 47
ACCEPTED MANUSCRIPT sequential ice breakup, which comprises successive runs of River1d using the output computed at the end of any one time step as input for the next time step but with a shortened “ice cover”. Rather than repeatedly modify and re-load the ice file, the open-water reach (lower value of nb)
T
is simply extended to the new location of the ice edge.
IP
APPENDIX B. FRONT PROFILE OF THE SELF-SUSTAINING WAVE
US
Y 1 U U U x g t g x
(B1)
AN
S f So
CR
The one-dimensional momentum equation reads (Henderson 1966):
M
in which U = flow velocity; and Sf = friction slope, defined as total boundary resistance force per
ED
unit channel length divided by flow area = average boundary shear stress divided by hydraulic
(B2)
CE
2
AC
nU S f 2/3 R
PT
radius. From the Manning resistance equation, Sf can be expressed in terms of U, i.e.
Here n = nc, R = Y/2 for ice-covered flow while n = nb, R = Y for open-water flow that exhibits the well-known logarithmic velocity distribution in the vertical direction.
With reference to the sketch of Fig. B1, three flow regions can be identified: (a) essentially openwater flow upstream of point A, with freely moving ice blocks on the surface; (b) ice-covered flow downstream of the breaking site (point B); and (c) transitional reach (AB) of “hybrid” flow. 48
ACCEPTED MANUSCRIPT Here, the surface velocity, which is zero just downstream of point B, will increase gradually in the upstream direction, from zero to that of the open-water reach at point A. Therefore the Manning resistance relationship, from which Eq. B2 derives, would not be accurate if the
ED
M
AN
US
CR
IP
T
hydraulic radius were set equal to Y or to Y/2.
CE
PT
Figure B1. Schematic illustration of vertical velocity distributions in different portions of the SSW front.
R
AC
To compute the transitional segment of the wave profile, the hydraulic radius is expressed as
Y 1 e x / L
(B3)
where x is distance measured from point B in the upstream direction. The hydraulic radius is equal to Y/2 at point B and asymptotically tends to Y as x increases. The length L is an empirical length scale, characterizing the extent of the transitional reach. The hydraulic radius 49
ACCEPTED MANUSCRIPT becomes practically equal to Y for x = 4L (denominator = 1.02), which can be considered to be the length of the segment AB in Fig. B1.
Introducing the transformations = x-Cbrt and = x-Cbrt (for the reaches BC and AB,
(reach BC);
dY d
S f So (C U ) 2 1 br gY
CR
So S f (C U ) 2 1 br gY
(reach AB)
(B4)
US
dY d
IP
T
respectively), and using also the equation of continuity (Eq. 3), Eq. B1 can be re-arranged to read
AN
Moreover, the flow velocity, U, can be obtained from the equation of continuity (Eq. 3), which is
M
first applied between any point along the front (Q = UYW) and the undisturbed flow downstream
PT
(B5)
CE
Y U C C 1 i Ui Ui Ui Y
ED
of the front (Qo = UiYiW), and then re-arranged to read:
AC
An obvious requirement for the profile of the wave front is that Y should decrease with increasing , i.e. dY/d < 0. For the present conditions, it can be shown that So-Sf < 0; therefore, the denominator of the first of equations B4 should be positive, which will be satisfied for
Cbr U i gYi Cg
(B6)
50
ACCEPTED MANUSCRIPT The quantity on the right-hand-side represents the celerity of a gravity wave (Cg) propagating in the ice-covered reach downstream of the SSW front.
The ordinary differential Eqs. B4 can be easily integrated in Excel, once So, Qo, W, Qbr, ni, nb,
T
are given and a trial value of Qp is selected. With the selected Qp, one can then compute Cbr from
IP
Eq. 14 and Ybr from Eq. 15; and subsequently compute Ubr as Qbr/YbrW and thence determine
CR
(Sf)br via Eq. B2 with x = 0. The computation of the profile is straightforward in the ice-covered
US
reach [with known (Sf)br determine dY/d and thence the value of Y at a small increment of downstream; calculate new values of U, Sf and dY/d; repeat with equal increments of , until Y
AN
no longer changes, becoming equal to Yi]. For the transitional reach, the same algorithm is
M
applied but in the upstream direction; here, the value of the hydraulic radius R is no longer equal to Y/2 but calculated from Eq. B3. The selected value of Qp is adjusted until the solution matches
ED
asymptotically the value of Yp that is obtained from the hydraulic resistance equation (Eq. 12).
PT
An example is shown in Fig. B2, where the wave front is seen to be slightly under 10 km long. This example corresponds to CRISSP Run 9 of Jasek et al. (2005), for which the development of
CE
the SSW was illustrated graphically; the length of the wave front was scaled off their
AC
corresponding figure and estimated as ~10 km.
51
ACCEPTED MANUSCRIPT 4.2 ice-covered reach
4
transitional reach
T IP
3.6
CR
4 3.4
3.2
0 -5000
US
2
0
3 -10000
5000
M
-5000
AN
Flow depth (m)
3.8
0
5000
ED
Channel distance (m)
PT
Figure B2. Wave profile for steeper channel (So = 0.0003), with W = 600 m, Qbase = 1600 m3/s, Qbr = 2400, ni=nb=0.03 and L ~ 1300 m. Ice breaking occurs at distance zero. The inset allows direct visual comparison of breaking and peak depths.
CE
The results of the CRISSP runs (Jasek et al 2005), were closely reproduced analytically (Fig. B3)
AC
using a constant value of L for each channel slope, i.e. ~3000 and ~1300 m for So = 0.00005 and 0.0003 respectively. For interpolation purposes, normalized values of L can be used, e.g. L/Yi = 529 (So = 0.00005) and 393 (So = 0.0003).
No solution to Eq. B4 could be found for the ice-covered portion of the profiles in the two instances (one for each slope) in which Cbr exceeded the value of Cg. These two profiles would be very short as they involve minute differences between the depths Ybr and Yi (less than 52
ACCEPTED MANUSCRIPT ~0.1%). It is possible that the discretization of the numerical model was too coarse to account for the equation of momentum along this short reach. The constraint imposed by Cg did not permit application of the analytical method to very small values of Qbr/Qo-1 (Fig. B3).
CR
IP
T
2.0
US
1.0
AN
QP / Qo - 1
1.5
M
CRISSP runs - bed slope = 0.0003
0.5
CRISSP runs - bed slope = 0.00005
ED
Analytical solution - bed slope = 0.0003
0.2
0.4
0.6
Equality line
0.8
1
1.2
1.4
Qbr / Qo - 1
CE
0
PT
0.0
Analytical solution - bed slope = 0.00005
AC
Figure B3. Results of present methodology versus CRISSP-generated relationships between sustained peak and threshold flows.
REFERENCES Ashton, G.D. 1986. River and Lake Ice Engineering. Water Resources Publications, Littleton, CO., USA.
53
ACCEPTED MANUSCRIPT Beltaos, S. 2007. The role of waves in ice-jam flooding of the Peace-Athabasca Delta. Hydrological Processes, 21(19), 2548-2559. Beltaos, S. 2009. River flow abstraction due to hydraulic storage at freezeup. Can. J. Civ. Eng. 36: 519-523.
T
Beltaos, S. 2011. Developing winter flow rating relationships using slope-area hydraulics. River
IP
Research and Applications. 27(9): 1076-1089.
CR
Beltaos, S. 2013a. Hydrodynamic and climatic drivers of ice breakup in the lower Mackenzie River. Cold Regions Science and Technology 95: 39-52.
US
Beltaos, S. 2013b. Hydrodynamic characteristics and effects of river waves caused by ice jam
AN
releases. Cold Regions Science and Technology 85: 42–55.
Beltaos, S. 2016. Extreme sediment pulses during ice breakup, Saint John River, Canada. Cold
M
Regions Science and Technology, 128: 38-46.
ED
Beltaos, S., and Burrell, B.C. 2005a. Field measurements of ice-jam-release surges. Canadian Journal of Civil Engineering, 32(4), 699-711.
PT
Beltaos, S., and Burrell, B.C. 2005b. Determining ice-jam surge characteristics from measured
CE
wave forms. Canadian Journal of Civil Engineering, 32(4), 687-698. Beltaos S. and Burrell, B.C. 2016. Transport of suspended sediment during the breakup of the ice
AC
cover, Saint John River, Canada. Cold Regions Science and Technology, 129: 1-13. Beltaos, S. and Carter, T. 2009. Field studies of ice breakup and jamming in lower Peace River, Canada. Cold Regions Science and Technology Journal. 56(2-3): 102-114.
Beltaos, S., Prowse, T.D. and Carter, T. 2006. Ice regime of the lower Peace River and ice-jam flooding of the Peace-Athabasca Delta.
Hydrological Processes, Volume 20, Issue 19,
pp 4009-4029.
54
ACCEPTED MANUSCRIPT Burrell, B.C. 2008. Chapter 1. Introduction. In: River Ice Breakup, Water Resources Publications, Highlands Ranch, Co., USA., pp 1-19. Carey, K.L. 1966. Observed Configuration and Computed Roughness of the Underside of River Ice, St. Croix River, Wisconsin, U.S. Geological Survey Professional Paper 550-B, B192-
T
B198.
IP
Carey, K.L. 1967. The Underside of River Ice, St. Croix River, Wisconsin, U.S. Geological
CR
Survey Professional Paper 575-C, C195-C199.
Ferrick, M.G., and N. D. Mulherin, 1989. Framework for control of dynamic ice breakup by river
US
regulation. Journal of Regulated Rivers: Research & Management, 3, pp. 79-92.
AN
Ferrick, M.G. and Goodman, N.J. 1998. Analysis of Linear and Monoclinal River Wave Solutions. USACE Cold Regions Research and Engineering Laboratory, CRREL Report 98-1,
M
Hanover, NH, USA, 24 p.
ED
Gerard, R. 1990. Hydrology of floating ice. In: Northern Hydrology, Canadian Perspectives, Prowse TD, Ommanney CSL (eds). NHRI Science Report No. 1, National Hydrology
PT
Research Institute, Environment Canada, Saskatoon, 103–134, plus references.
CE
Gerard, R., Kent, T.D., Janowicz, R., and Lyons, R.O. 1984. Ice regime reconnaissance, Yukon River, Yukon. Proceedings of the 3rd International Specialty Conference on Cold Regions
AC
Engineering, 4–6 April 1984, Edmonton, Alta. Compiled and edited by D.W. Smith. Canadian Society for Civil Engineering, Montréal, Quebec, pp. 1059–1073. Henderson, F.M. 1966. Open channel flow. The Macmillan Co., New York, Toronto. Henderson, F.M., and Gerard, R. 1981. Flood waves caused by ice jam formation and failure, Proceedings IAHR Symposium on Ice, Quebec, Canada, Vol. I, pp. 277-287.
55
ACCEPTED MANUSCRIPT Hicks, F., Cui, W, and Ashton, G. 2008. Chapter 4. Heat transfer and ice cover decay. In: River Ice Breakup, Water Resources Publications, Highlands Ranch, Co., USA., pp 67-123. Jasek, M. and Pryse-Phillips, A. 2015. Influence of the proposed Site C hydroelectric project on the ice regime of the Peace River. Canadian Journal of Civil Engineering, 2015, 42(9): 645-
T
655.
IP
Jasek, M., Ashton, G., Shen, H.T., and Chen, F. 2005. Numerical Modeling of Storage Release
CR
during Dynamic River Ice Break-up. Proceedings (CD-ROM) of 13th Workshop on the Hydraulics of Ice Covered Rivers, Hanover, NH, September 15-16, 2005, CGU-HS
US
Committee on River Ice Processes and the Environment, Edmonton, Canada, 421-439.
AN
Kellerhals, R., Neill, C.R. and Bray, D.I. 1972. Hydraulic and geomorphic characteristics of rivers in Alberta. River Engineering and Surface Hydrology Report 72-1, Research Council of Alberta,
M
Edmonton, Canada, 52 p.
ED
Parkinson, F.E. 1982. Water temperature observations during break-up on the Liard-Mackenzie River system. Workshop on the Hydraulics of Ice Covered Rivers, Edmonton, Canada, 261-295.
PT
Ponce, V.M., Simons, D.B., 1977. Shallow wave propagation in open channel flow. Journal of
CE
the Hydraulics Division, ASCE 103 (HY12), 1461–1476. Prowse, T.D. and Carter, T. 2002. Significance of ice-induced storage to spring runoff: a case
AC
study of the Mackenzie River. Hydrological Processes 16(4): 779–788. Prowse, T.D. 1986. Ice jam characteristics, Liard-Mackenzie rivers confluence. Canadian Journal of Civil Engineering, 13(6), 653-665. She, Y. and Hicks, F. 2006. Modeling ice jam release waves with consideration for ice effects. Cold Regions Science and Technology 45 (2006) 137–147.
56
ACCEPTED MANUSCRIPT HIGHLIGHTS
-Storage release during river ice breakup can generate large additions to flow -It can be caused by ablation or by ice breaking
T
-Hydrodynamic properties of each process studied analytically and numerically
IP
-Ice breaking is far more dynamic and may lead to formation of self-sustaining waves
AC
CE
PT
ED
M
AN
US
CR
-Results applied to a case study and to assess climate impacts on ice-jam flooding
57
Spyros Beltaos PII: DOI: Reference:
S0165-232X(16)30416-5 doi: 10.1016/j.coldregions.2017.04.009 COLTEC 2388
To appear in:
Cold Regions Science and Technology
Received date: Revised date: Accepted date:
19 December 2016 31 March 2017 29 April 2017
Please cite this article as: Spyros Beltaos , Hydrodynamics of storage release during river ice breakup, Cold Regions Science and Technology (2017), doi: 10.1016/ j.coldregions.2017.04.009
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
ACCEPTED MANUSCRIPT Hydrodynamics of storage release during river ice breakup Spyros Beltaos
IP
T
Watershed Hydrology and Ecology Research Division, Canada Centre for Inland Waters, Environment and Climate Change Canada, Burlington, Ontario, Canada [email protected]
CR
ABSTRACT
AC
KEY WORDS
CE
PT
ED
M
AN
US
The breakup of the ice cover in cold-region rivers is a brief but seminal period of their hydrologic regime, with important ecological and socio-economic implications. The main driver of ice breakup processes is the flow discharge hydrograph. It is generated by runoff from snowmelt and rainfall but can be modified by rapid release of water from storage as the ice cover recedes by ablation and mechanical breaking up. Despite its potential importance, there is very limited quantitative information concerning the hydrodynamic processes that control storage release during ice breakup in rivers, while the issue of climate change underscores the need for improved understanding of the relevant mechanisms. Quantitative analysis for assumed prismatic channels shows that ablation and sustained ice dislodgment and breaking can cause significant flow enhancement via storage release. The latter process is far more dynamic than ablation and, under certain conditions, may lead to formation of a self-sustaining wave (SSW). Analytical results are applied to natural stream conditions, using the Lower Peace River as a case study. Observed rates of ice recession typically indicate ice melt as the dominant process. A rare occurrence of rapid ice breaking over hundreds of kilometres (2014) indicated partial agreement with the SSW concept and revealed that discrepancies may arise from the characteristic irregularity of rivers. Impacts on storage release by climate-driven changes to river ice regimes are examined, along with their implications to ice-jam formation and associated flooding.
Breakup, celerity, climate, flow enhancement, ice cover, melt, self-sustaining wave, storage release
1
ACCEPTED MANUSCRIPT 1. INTRODUCTION The breakup of the ice cover in cold-region rivers is a brief but seminal period of their hydrologic regime with important ecological and socio-economic implications (Burrell 2008). Major ice jams that often form during breakup can cause extreme floods that imperil human life
T
and damage property and infrastructure (Ashton 1986). Navigation and hydro-power generation
IP
are also adversely affected by ice jamming. At the same time, ice-jam flooding can provide
CR
essential replenishment to the multitude of lakes and ponds characteristic of the northern Canadian deltas, which are havens for wildlife, especially waterfowl and aquatic animals. On the
US
other hand, the sharp waves produced by ice jam releases (also termed “javes”) can have
AN
negative ecological impacts, owing to amplified velocity and erosive capacity (Beltaos and
M
Burrell 2016, Beltaos 2016).
ED
The main driver of ice breakup processes is the flow discharge hydrograph. It is generated by runoff from snowmelt and rainfall but can be modified by highly dynamic processes resulting
PT
from the rapid release of water from storage, as the ice cover recedes by ablation and mechanical
CE
breaking up. Storage release has been postulated to play a role in rapid ice clearance that has occasionally been observed in rivers of northern Canada. An extreme occurrence took place on
AC
the Yukon River in 1983: Gerard et al (1984) reported apparent ice breaking over a 300 km stretch at an average celerity of 5.2 m/s, which reduced the intact ice cover to a 140 km long ice run. Such events could potentially be explained by formation and advance of non-attenuating waves that are sustained by storage release as the ice cover is being broken up. This effect was detected in exploratory numerical runs (Ferrick and Mulherin 1989) and quantified, also numerically, by Jasek et al (2005). The results of the latter study are used extensively in Section
2
ACCEPTED MANUSCRIPT 4 and Appendix B, where the analysis and theory of the self-sustaining wave (SSW for short) are developed.
The main objective of this paper is to systematically examine and, where possible, quantify
T
storage release mechanisms and their effects on flow magnitude and breakup progression. A
IP
secondary objective pertains to how changes in ice regime that may be caused by climatic
CR
warming can alter storage release processes.
US
Following a background section, storage release by melting of the ice cover is quantified for
AN
different water temperatures and channel characteristics. The much more dynamic process of storage release by ice breaking is considered next and the potential formation of self-sustaining
M
waves illustrated. Unique properties of the SSW are discussed and the links between ice-
ED
breaking criteria and wave celerity elucidated. Application to natural streams is illustrated next
PT
with a case study, followed by a discussion of potential impacts of changing ice regimes.
CE
2. BACKGROUND INFORMATION As has been described in previous literature (Gerard 1990, Prowse and Carter 2002, Beltaos
AC
2009), the formation of an ice cover causes river water levels to increase in order to accommodate the keel of the cover and the hydraulic resistance generated by the bottom boundary of the cover, in addition to that of the river bed. The volume of water between the prevailing and the initial (before ice cover formation) water levels is called “storage”, even though the “stored” water is not stationary but moves along the river with the rest of the flow,
3
ACCEPTED MANUSCRIPT albeit at reduced velocity. As water is being stored during the freezeup period, downstream flow is reduced and this effect is often captured in the records of various hydrometric gauges.
During the winter, the river flow typically decreases and so does the roughness of the bottom
T
surface of the ice cover. In turn, water levels decrease and some of the stored water is released
IP
back into the flow. By the time spring runoff begins, there is a relatively small amount of water
CR
in storage, since the prevailing flow and thence the “backwater” is small. The term backwater
US
denotes the difference between ice-covered and open-water flow depth at the same discharge.
AN
With the arrival of mild weather in the spring (or during a mid-winter thaw), the flow begins to increase while thermal effects cause the ice bottom surface to become rougher (Carey 1966,
M
1967). The backwater will increase as a result and so will the storage until the ice cover melts or
ED
is forcibly mobilized by the rising flow. In either case, there will be a release of storage, which will cause the downstream flow to increase over and above the “carrier” flow in the river. Carrier
PT
flow is the flow that is strictly attributable to runoff and is the flow that prevails in the river
CE
before dynamic effects, such as javes and storage-release, appear or after they cease. In the
AC
following, the term “excess” flow will be used to denote the increment caused by storage release.
When a segment of the ice cover is removed, flow resistance is significantly reduced. Local velocity and discharge increase, and thence “stored” water is released. To visualize this effect, a numerical experiment was carried out using the River1d model (2004 version; kindly supplied to the writer by Professor Faye Hicks). For simplicity, the open water version of the model was used while the backwater effect of the ice cover was simulated by doubling the bed Manning
4
ACCEPTED MANUSCRIPT coefficient, nb, from 0.03 to 0.06 in the “ice-covered” reach of a 32 km long prismatic rectangular channel (flow = 343 m3/s. width = 100 m, slope = 0.0005). Using Eq. A1 of Appendix A, it can be shown that open-water flow with nb = 0.06 is equivalent, in terms of depth and mean velocity, to ice-covered flow with nb = 0.03 and ni = Manning coefficient of the ice
IP
T
cover = 0.045).
CR
The “ice cover” extended from km 10 to 32. The resulting water level profile (Fig. 1) resembles what would be expected if an ice cover were present at and beyond km 10. At time 0, the
US
Manning coefficient is reduced to 0.03 (from 0.06) between kilometres 10 and 10.5, simulating
AN
removal of a 0.5 km segment of ice cover. This change generates a small wave (Fig. 2) that ultimately leads to an asymptotic configuration, which corresponds to a 0.5 km forward
M
translation of the initial water surface profile (Fig. 1). The computation indicated that it took 3.5
ED
h to attain a profile that is fairly close to the asymptotic one (Fig. 2). Moreover, it was found that the discharge increased along the wave-affected segment of the channel, always attaining a peak
PT
at km 10; the highest such peak (377 m3/s) occurred at ~0.1 h. The peak value gradually
AC
CE
decreased afterward, falling to 350 m3/s at 5.0 h.
5
ACCEPTED MANUSCRIPT 104 Channel bed
102
Water surface
98
T
96
IP
94
92 90 5
10
15
Channel distance (km)
20
US
0
CR
Elevation (m)
100
AN
Figure 1. Initial flow configuration used to numerically simulate the effects of storage release. Discharge = 343 m3/s, width = 100 m, slope = 0.0005, nb = 0.03 (0 to 10 km), 0.06 (10 to 32 km).
M
0.15
ED PT
0.05
CE
0.00 -0.05 -0.10
-0.15
0.1 h
0.4 h 1.5 h
AC
Wave disturbance (m)
0.10
3.5 h 6.0 h Asymptotic
-0.20 0
5
10
15
20
25
30
Channel distance (km)
Fig. 2. Spatial variation of change in water level at different times following simulated removal of a 0.5 km long “ice cover” segment from the channel of Fig. 1. 6
ACCEPTED MANUSCRIPT 3. ICE MELT: THERMAL BREAKUP Low spring runoff typically results in thermally driven ice clearance. The flow is not sufficiently high to cause mechanical breaking and mobilization of the ice cover. Instead the cover largely melts in place. The resulting process of storage release is discussed in this section by considering
T
the simplified situation depicted in Fig. 3. Steady flow in an open-water reach of a prismatic
IP
rectangular channel is subjected to thermal inputs and the water temperature increases with
CR
downstream distance until the flow encounters a long stretch of stationary ice cover. As the water dives under the ice cover, its temperature decreases gradually to zero because heat is
US
extracted by ice melt. Assuming that the incoming flow and water temperature do not change
AN
over time, the short transitional reach AB will simply be translated downstream at a rate (Hicks
ED
c pQTW i LFW
(1)
PT
Cmelt
M
et al 2008):
CE
in which TW = temperature of the water just upstream of the edge of the ice cover; Q = river discharge arriving at the ice edge; cp = specific heat of water (4,220 J/kgC at 0C); = thickness
AC
of the ice cover (m); = density of water (1000 kg/m3); i = density of freshwater ice (916 kg/m3); LF = latent heat of fusion (3.34105 J/kg); and W = width of the ice cover.
7
US
CR
IP
T
ACCEPTED MANUSCRIPT
AN
Figure 3. Schematic illustration of a wave generated by ice melt. The wave front is depicted as a straight line for simplicity. Melting results in zero thickness at the ice edge. The symbols C, Q, and Y represent celerity, discharge, and water depth, respectively.
M
The melting of the ice cover can be viewed as a series of successive removals of short ice
ED
segments at equal time intervals, such that each removal generates an elementary wave similar to that of Fig. 2. The resulting composite wave was numerically simulated by removing 0.5 km
PT
“ice” segments every 0.2 h (12 min) from the initial profile of Fig. 1. The speed of ice removal is
CE
thus 2.5 km/h or 0.69 m /s; this is a relatively large value that, depending on channel conditions, would require a water temperature of a few to several degrees C. The resulting wave forms at
AC
different times are illustrated in Fig. 4. The peak of the wave always coincides with the edge of the ice cover, that is, it moves at the melting rate. At the same time, the speed of advance of the leading edge of the wave is considerably larger (1.9 m/s indicated by the first few waves).
Owing to the relatively short length of the channel, the computed shape of the wave front was, after ~3 h, influenced by the downstream boundary condition (fixed water depth). Nevertheless, the results clearly indicate that the wave height tends towards an asymptotic value that is nearly 8
ACCEPTED MANUSCRIPT attained within hours from the start of melt. At the same time, the wave front elongates, owing to the faster advance of the leading edge relative to the ice edge. In turn, this implies that the water surface slope of the wave front will eventually approach the value of the bed slope (So) of the channel. It is only after the time when this configuration is established that the flow downstream
IP
T
of the ice edge can be considered approximately uniform with a flow equal to the peak value.
CR
5
US AN M
3
2
ED
Depth (m)
4
0
5
AC
0
CE
PT
1
Initial profile profiles at 0.6, 1.6, 2.4, 3.0, 3.8 h asymptote
10 15 20 River distance (km)
25
30
Figure 4. Temporal development of wave generated by sequential release of 0.5 km “ice segments” at 12 minute (0.2 h) intervals. The wave peak is always located at the “ice edge” and tends to an asymptotic value of 4.71 m. The “ice cover” is simulated by simply increasing the Manning coefficient from 0.03 to 0.06.
To determine the asymptotic flow properties, one may start with the differential equation of continuity:
9
ACCEPTED MANUSCRIPT
Y q 0 t x
(2)
T
in which t = time, x = downstream distance; Y = flow depth (measured from the water surface or
IP
from the bottom of the ice cover, as applicable); and q = discharge per unit channel width = Q
CR
(discharge)/W (width). When a waveform propagates at constant rate, C, a single coordinate can
US
be defined as = x-Ct and Eq. 2 transformed to an ordinary differential equation, i.e.
(3)
M
AN
d (CY q) 0 d
which implies that the quantity CY-q is constant. Because C changes downstream of the peak of
ED
the wave, it is only the portion upstream of the peak where Eq. 3 applies. The largest excess flow
PT
due to storage release occurs at the ice edge; it will remain approximately equal to this value downstream of the edge, once an adequate time has elapsed since the start of melt to ensure that
AC
CE
the frontal water surface slope is approximately equal to So.
Application of Eq. 3 between the ice edge and a section far upstream, where uniform open-water flow prevails, results in:
Q Qe Qu Cmelt (Ye Yu )W
(4)
10
ACCEPTED MANUSCRIPT in which Q is the excess flow (storage-release contribution to discharge) and the subscripts “e” and “u” denote edge and far upstream conditions. By defining Ye as water depth (water surface to channel bed), the volume of melted ice is included in Qe. One may note here that the depth Ye is greater than what prevailed in the ice-covered reach before the start of melt (Fig. 4) because it
T
is accommodating a greater discharge than the carrier discharge, Qu. In general, therefore, it is
IP
not possible to determine Qe from Eq. 4 alone, because both Ye and Qe are unknown. However,
CR
when the condition of near-uniformity of the flow downstream of the ice edge is established, an
US
additional relationship may be obtained via the Manning resistance equation:
(5)
M
AN
W (Ye si )5/ 3 So Qe 22 / 3 nc
ED
in which si = specific gravity of ice = 0.916; and nc is the composite-flow Manning coefficient,
CE
Sabaneev formula:
PT
which can be computed in terms of ni and nb using the simplified version of the Belokon-
(6)
AC
nc [(ni 3/ 2 nb3/ 2 ) / 2]2 / 3
Noting that Qu and Yu are related by a similar relationship to that of Eq. 5, it is possible to rearrange Eq. 4 into a dimensionless form, i.e:
3/ 5 Q Qe Cmelt 22 / 3 nc Qe si 1 1 Qu Qu U u nb Qu Yu
(7)
11
ACCEPTED MANUSCRIPT
in which Uu = mean velocity in the uniform-flow, open-water reach. In Eq. 1, the amount of heat arriving at the edge of the ice cover is expressed in terms of local discharge and local water temperature. In the present application, the excess discharge has zero temperature, since it
T
derives from under-ice storage. Consequently, the heat available for melt is taken as (cpQTW)u,
IP
i.e at point A rather than at point B where the discharge will be higher than Qu but the water
CR
temperature lower than the temperature at point A. [This approximation neglects any heat that
US
may be gained in the transitional river segment between points A and B in Fig. 3, which is considered small but can be taken into account using more elaborate equations]. Equation 7 can
AN
be solved by trial and error to determine Qe and thence the excess flow Q.
M
As an example, consider a 600 m wide channel with carrier flow of 3000 m3/s, slope of 0.00005,
ED
Manning coefficients nb = ni = nc= 0.03, and ice thickness = 0.7 m. This type of channel could be
PT
viewed as a prismatic approximation to the lower Peace River (Kellerhals et al 1972; Jasek and Pryse-Phillips 2015). Typically, the winter ice cover is smoother than the river bed but becomes
CE
rougher during the pre-breakup period, so that ni approaches nb (Carey 1966, 1967). This
AC
assumption is adopted here for illustration purposes; in real-life applications actual value of these coefficients should be determined based on in situ hydraulic data. By application of the flow resistance equation to the upstream open-water flow, the depth Yu is calculated as [0.033000/6000.00005]0.6 = 6.25 m; and thence Uu = 3000/600/6.25=0.80 m/s. The rate of melt, Cmelt, is then calculated from Eq. 1 as ~0.10 m/s (8.6 km/day) for each degree of the water temperature.
12
ACCEPTED MANUSCRIPT
40 bed slope = 0.00005 bed slope = 0.0005
T IP
20
CR
Excess flow (%)
30
US
10
0
1
AN
0
2
3
4
5
6
M
Water temperature (°C)
PT
ED
Figure 5. Addition to carrier flow caused by melt-induced storage release in two hypothetical channels, both with 0.7 m thick ice covers. Channel 1: So = 0.00005, W = 600 m, Qu = 3000 m3/s, Cmelt = 8.6 km/day per degree of water temperature; Channel 2: S o = 0.0005, W = 200 m, Qu = 800 m3/s, Cmelt = 6.8 km/day per degree of water temperature
CE
With the aid of Eq. 7, the excess flow can be calculated for different water temperatures and is plotted as a percentage in Fig. 5, where it is seen to exceed 10% for Tw > 2°C. Also shown in
AC
Fig. 5 are the results of a second calculation, this time for a channel that approximates average properties of the upper Saint John River. The higher slope greatly reduces the excess flow, which remains under 10% of the carrier flow, even if the water temperature were to increase to 5°C. In situ measurements that were carried out in upper Saint John River during the breakup in the 1990s have actually indicated water temperatures not exceeding 2°C. On the other hand, the large solar radiation flux associated with spring conditions at high latitudes and the extensive lengths of major rivers of Northern Canada suggest that water temperatures of several degrees C 13
ACCEPTED MANUSCRIPT may be realistic. Measurements in the lower Mackenzie River indicated values up to 10°C during the highly thermal breakup event of 1980 (Parkinson 1982). Considering that the excess flow can be substantial in the low-slope channel (Fig. 5), it is conceivable that the excess may enable
T
mechanical ice breaking farther downstream.
IP
To determine the time after which Eqs. 5 and 7 can be applied with some confidence beyond the
Ye Yu 1) Ye Yi
US
Cmelt (2
(8)
AN
CL
CR
ice edge, the leading-edge celerity CL is estimated from the approximate relationship:
M
in which Yi is the water depth in the ice-covered reach far downstream from the ice edge where
ED
the water level has not yet been affected by the wave; and the average celerity for the wave front is assumed to be simply equal to (Cmelt+CL)/2. To attain a wave slope that is within 5% of So, the
20
Ye Yi CL S o
(9)
AC
TL
CE
PT
leading edge of the wave would have to travel for a time TL:
As an example, the value of Yi for the low-slope channel is calculated as 8.89 m, while Ye = 9.37 m for TW = 2°C, which corresponds to Cmelt = 0.2 m/s. With Yu = 6.25 m, as noted earlier, CL works out to be 2.4 m/s and thence TL = 20 (9.37-8.89)/(2.40.000053600) ~ 22 h after the start of melt. During this time interval, the ice edge would be located 2236000.2/1000 ~ 16 km downstream from its initial location. 14
ACCEPTED MANUSCRIPT
Of course, the preceding calculations are predicated on constant water temperature at the ice edge. But as the open-water reach lengthens, the water temperature will increase and thence the melting rate would be expected to increase as well. Variable weather conditions, ice thickness,
T
and runoff intensity will further complicate the ablation-driven storage release process.
IP
Nevertheless, the present examples capture the main elements of the process and elucidate the
CR
magnitudes of the important variables such as the melt rate and the resulting excess discharge. The analysis also shows that considerable melting time must elapse before the flow underneath
AN
US
the ice cover becomes nearly uniform.
4. ICE DISLODGMENT AND BREAKING: THE SELF-SUSTAINING WAVE
M
Much more dynamic conditions will prevail in cases where the rising flow dislodges and
ED
mobilizes the ice cover, which then quickly breaks down into smaller slabs and blocks. At the same time, the remaining intact ice cover maintains a good portion of its mechanical strength
PT
(“mechanical” breakup) and can, under certain geomorphic conditions, arrest the ice run and
CE
initiate an ice jam. As stationary ice is being mobilized, new water volume is released from hydraulic storage and added to the flow. This effect counteracts the natural attenuation of the
AC
resulting wave and has been postulated via numerical experiments to render it self-sustaining (Ferrick and Mulherin 1989). Self-sustaining waves could perhaps explain reports of rapid breakup progression over very long reaches (Gerard et al. 1984).
Jasek et al. (2005) applied the CRISSP model to study wave propagation in ice-covered prismatic rectangular channels, assuming that ice dislodgment and breaking occur at a specified
15
ACCEPTED MANUSCRIPT discharge (Qbr), which exceeded the base flow (Qo) by varying amounts. Ice cover dislodgment and subsequent breaking into ice blocks is more complex than what is envisaged by singlevariable thresholds. It involves additional parameters, e.g. hydrodynamic forces applied on the ice cover, channel width and planform curvature as well as freezeup water level (Beltaos 2007,
T
2013b). Nevertheless, the single-variable criterion greatly simplifies the study of storage release
IP
by ice breaking and is rigorous for the assumed prismatic channel, for which the critical
CR
condition could equally well be expressed in terms of velocity, depth, or shear stress.
US
Dynamic conditions were generated by applying carrier flow hydrographs of different peaks and
AN
durations at the upstream boundary of the computational reach, which extended for hundreds of kilometres. If the incoming wave attained flows in excess of Qbr, the wave became self-
M
sustaining after an evolutionary period; it comprised a steep front of constant celerity, followed
ED
by a protracted crest of constant depth and discharge (herein termed sustained peak depth and sustained peak flow). The numerical runs indicated that the celerity of the front (Cbr) increased as
PT
Qbr decreased, becoming extremely large (~10 m/s) when Qbr was just a little higher than the base
AC
the channel width.
CE
flow, denoted by Qo. The results were also influenced by the slope of the channel bed but not by
The waves studied by Jasek et al (2005) were triggered by a flow hydrograph rather than by javes that result from ice jam releases, as is common in relatively flat rivers (Beltaos and Carter 2009; Beltaos 2013a). Some of these hydrographs were very steep while others were more gradual, approximating snowmelt runoff. Such differences did not seem to influence the properties of the SSWs, which often attained higher peak flows than the peaks of the respective
16
ACCEPTED MANUSCRIPT input hydrographs. To examine whether jave-induced ice breaking can also lead to formation of a SSW, the writer used the model River1d to study the release of a 9 km long jam in a 120 km prismatic channel of slope 0.5 m/km, width 100 m, base flow of 400 m3/s, and bed Manning coefficient, nb = 0.03. To simulate the backwater caused by the ice cover downstream of the jam,
T
nb was increased to 0.05 while the calculation was performed in open-water mode throughout the
IP
length of the channel (equivalent to ice-covered flow with nb = 0.03 and ni = 0.033; see also
CR
Appendix A). The shear stress in the “ice-covered” reach was calculated as one-half of gYSf (Sf
US
= friction slope; = density of water = 1000 kg/m3, g = gravitational acceleration = 9.81 m/s2). When the shear stress at both ends of the first 1 km long “ice segment” exceeded 20 Pa, the
AN
segment was removed, i.e. nb was reduced to 0.03 and the computation repeated using the latest water level profile as the initial condition. The 20 Pa threshold is higher than the pre-jave value
ED
M
but low enough to ensure that the jave could initiate ice breaking.
PT
Figure 6 shows that the jave quickly formed a front of constant shape and height, unlike the attenuating configurations of javes that do not dislodge the ice cover (Beltaos 2013b).
CE
Instantaneous depth profiles at different times from the release of the jam are shown in Figure 7.
AC
After 46 minutes, a short segment of practically constant depth has already developed upstream of the wave front; as time goes on, it is expected to get longer, much as was found by Jasek et al (2005). The sustained peak flow and depth were smaller than respective values associated with the initial jave (e.g. 5-minute graph). The magnitude of the celerity of the wave front (5.1 m/s) lies between those of the gravity and the kinematic waves (7.1 m/s and 1.6 m/s, respectively).
17
ACCEPTED MANUSCRIPT 205 initial WL Riverbed
200
T IP
190
CR
Elevation (m)
195
US
185
0
5
10
AN
180 15
20
25
30
35
40
45
M
River distance (km)
AC
CE
PT
ED
Figure 6. River1d simulation of ice jam release and formation of sell-sustaining jave (kms 0 to 15: open water; kms 15 to 24: ice jam; kms 24 to 45: intact “ice cover”, simulated using elevated value of nb. The coloured lines represent instantaneous water level profiles up to 46 minutes following release. Front celerity ~5.1 m/s; front height ~ 1.7 m above pre-jave water depth downstream of the jam.
18
ACCEPTED MANUSCRIPT 7 0 min.-jam 5 min.
9 min.
6
13 min. 26 min.
5
T
33 min.
IP
40 min. 46 min.
4
CR
Depth (m)
19 min.
US
3
0
5
10
AN
2 15
20
25
30
35
40
45
M
River distance (km)
CE
PT
ED
Figure 7. Instantaneous longitudinal variations of water depth at different times following release of the jam depicted in Figure 6 (kms 0 to 8: open water segment unaffected by jam backwater; kms 8 to 15: open water segment affected by jam backwater; jam extends from km 15 to km 24; intact “ice cover” beyond km 24).
The front and the crest of the wave remain essentially unchanged after ~30 minutes from the
AC
release of the jam, a much shorter development time than was found in the CRISSP runs (many hours). The constancy of the wave front and crest shape implies that all points along the front advance at the same celerity and there is no attenuation. This is an important difference from ordinary javes, which attenuate and spread over time, owing to increasing celerity along the waveform, from the trailing to the leading edge (Beltaos and Burrell 2005a). Unlike the case of ice melt, where the front of the wave increases in length and propagates under the intact ice
19
ACCEPTED MANUSCRIPT cover, the SSW wave has a sharp front and only a portion of this front is moving under intact ice. The flow depth at which the ice cover is dislodged (Ybr) is not much less than the sustained peak depth (Yp); ratio Ybr/Yp = 0.89.
T
A simplified sketch of the SSW is shown in Fig. 8. If the ice breaking “trigger” is a jave, there is
IP
no net inflow within the control volume (Q = Qu at both ends). The same applies to the CRISSP
CR
runs after the time when the input hydrograph returns to the base flow. In such instances,
US
continuity requires that
(10)
M
AN
Cbr (Yp Yi ) Ctail (Yp Yu )
in which Ctail = celerity of the tail end of the wave and Yi = depth of flow under the undisturbed
ED
ice cover downstream of the wave front, not to be confused with the water depth Yi (Fig. 3)
PT
which is measured from the water surface to the bed. Clearly, Ctail < Cbr. Therefore the uniformflow crest AB will lengthen as time goes on. The few wave shapes shown in Fig. 7 only hint at
AC
CE
this feature, which is more conspicuous in the extensive CRISSP runs (Jasek et al 2005).
20
US
CR
IP
T
ACCEPTED MANUSCRIPT
AN
Figure 8. Schematic illustration of wave generated by sustained ice dislodgment. The wave front is depicted as a straight line segment for simplicity.
M
The preceding results suggest further that the properties of the SSW do not depend on whether it is triggered by rising runoff or by jam release. Once ice breaking is initiated and the wave has
ED
travelled for a sufficiently long time, its front takes on a constant shape and moves at constant
PT
rate. These parameters are primarily dictated by the threshold criterion, be it a threshold
CE
discharge or a threshold shear stress.
AC
Figure 9, which derives from CRISSP results provided in Table 1 of Jasek et al (2005), shows how the celerity of the wave front (Cbr) decreases as the threshold discharge increases. Extremely rapid ice breaking is indicated when this threshold is just above the base flow while a strong negative influence of the bed slope is also evident. One may note further that Cbr is always greater than the celerity of the kinematic wave (~1.5Ui) and generally differs from the celerity of the gravity wave, Cg, which can be expressed in dimensionless form as:
21
ACCEPTED MANUSCRIPT Cg Ui
gYi Ui
1
(11)
30 CRISSP runs - bed slope = 0.0003 CRISSP runs - bed slope = 0.00005
T
25
CR
Cg/Ui - slope = 0.00005
15
US
Cbr / Ui
20
IP
Ice jam release - bed slope = 0.0005
AN
10
Cg/Ui - slope = 0.0003 0 1.4
1.6
ED
1.2
1.8
2
2.2
2.4
Qbr / Qo
PT
1
M
5
CE
Figure 9. Variation of ice breaking celerity with threshold discharge, both normalized with respective base quantities (the CRISSP runs are described in Jasek et al 2005; Ui was calculated by the writer via the Manning resistance equation)
AC
For So = 0.0003 and 0.00005, Cg/Ui amounts to 8.06 and 16.80, respectively. As shown in Appendix B, it is not physically correct for Cbr to exceed Cg and this is consistent with previous work on javes and open-water waves (Ponce and Simons, 1977; Ferrick and Goodman 1998; Beltaos 2013b). On the other hand, Fig. 9 suggests that the upper-most data point of each set exceeded the corresponding value of Cg/Ui. This could be the result of numerical-discretization, as discussed in Appendix B.
22
ACCEPTED MANUSCRIPT Figure 10 shows that the sustained peak flow, Qp, is primarily determined by the threshold discharge, as well as by the base flow and the channel slope. Since both Qp and Qbr exceed Qo, excess ratios are plotted in Fig. 10 for convenience. Within their respective ranges, the data points for each slope are well described by the corresponding linear fits. For the tested bed
T
slopes, the flow contributed by storage release is seen to range from ~70% to ~200% of the base
IP
flow. Moreover, it seems likely that higher percentages can be attained at higher bed slopes.
CR
These figures underscore the dynamic nature of the SSW relative to melt-induced waves, which
2.5 CRISSP runs- bed slope = 0.0003
2.0
AN
CRISSP runs - bed slope = 0.00005
Ice jam release - bed slope = 0.0005
M
1.5
ED
Q P / Qo - 1
equality line
CE
PT
1.0
AC
0.5
0.0
US
generate excess flows under ~40% of the base flow (Fig. 5).
0
0.2
0.4
0.6
0.8
1
1.2
1.4
Qbr / Qo - 1 Figure 10. Variation of sustained peak flow with threshold discharge; both normalized with the base discharge Qo and expressed as “excess” values by subtracting 1 from each ratio. The CRISSP runs are described in Jasek et al (2005).
23
ACCEPTED MANUSCRIPT Since the open-water flow in the sustained-peak reach is uniform, one can write:
Yp
2
nbQ p ncQo
3/ 5
(12)
T
Yi
2 / 5
IP
To derive Eq. 12, it was assumed that the ice blocks being carried by the flow in the sustained-
CR
peak reach do not impede water motion. Jasek et al (2005) considered possible impedance by simply increasing the local value of nb to 0.04 and 0.05 but the resulting wave profiles were
AN
US
deemed physically implausible. Consequently, the results of these runs are not considered herein.
ED
(13)
PT
Q p Qo Cbr (Y p Yi )W
M
Continuity (Eq. 3) requires that
CE
Combining Eqs 12 and 13 results in:
(14)
AC
(Q p / Qo ) 1 Cbr 2 / 5 U i 2 ( nbQ p / ncQo )3/ 5 1
which suggests that the normalized ice breaking celerity is uniquely defined by the normalized sustained peak flow and by the ratio of the bed and composite-flow Manning coefficients nb and nc.
24
ACCEPTED MANUSCRIPT The relationship of Eq. 14 is validated in Fig. 11 where the various data points collapse onto a single curve. If the line shown in Fig. 11 were extended beyond the right-hand margin, it would be seen to attain a minimum of Cbr/Ui = 4.03 and then begin to increase slowly. During the prebreakup phase of the ice cover, the velocity Ui is often in the range 0.7 to 1.0 m/s. Consequently,
T
SSWs are expected to advance at speeds of no less than ~2.8 m/s (for nb/nc=1). The latter value
IP
amounts to ~240 km/day, which is much greater than what can be expected from ice melt, and
US
30
CR
again illustrates the dynamic nature of the SSW.
CRISSP runs - bed slope = 0.0003
25
AN
CRISSP runs - bed slope = 0.00005 Ice jam release - bed slope = 0.0005
M ED
Predicted from continuity equation
15
PT
Cbr /Ui
20
CE
10
0
AC
5
1.5
2
2.5
3
3.5
Qp /Qo
Figure 11. Normalized ice breaking celerity versus normalized peak wave discharge for nb = nc (the CRISSP runs are described in Jasek et al 2005; Ui was calculated by the writer via the Manning resistance equation). The line is not physically meaningful where Cbr > Cg
25
ACCEPTED MANUSCRIPT Equation 14 indicates further that the ratio Qp/Qo must exceed 22/3(nc/nb)5/3. Below this limit, Cbr would be calculated as negative while at the limit it becomes infinite. For nb = nc, this limit takes on the value 22/3 or 1.587. This means that Qp/Qo – 1 has to exceed 0.587 which is seen to be the case for all data points in Fig. 10 and implies that SSW storage release contributes no less than
IP
T
~60% of the base flow.
CR
The ratio nb/nc has a strong effect on Cbr, as illustrated in Fig. 12. Increasing values of nb/nc result in lower values of Cbr/Uo for any given ratio Qp/Qo. The value of nb/nc would typically be
US
expected to be equal to 1 (e.g. nb = ni = nc = 0.03) or more than 1 (e.g. nb=0.03, ni = 0.02, nc =
AN
0.025). Cases in which nb/nc is less than 1 require that ni be greater than nb, which is not likely for sheet ice covers, but could occur where the ice cover has a very rough underside owing to
AC
CE
PT
ED
M
slush deposition during freezeup or possibly jamming of ice blocks during mid-winter breakup.
26
ACCEPTED MANUSCRIPT
20
nb/nc = 0.8
15
0.9
T IP
10 1.1
CR
Cbr /Ui
1.0
1.2
AN
US
5
0
2
3
4
5
Qp /Qo
ED
M
1
PT
Figure 12. Normalized ice breaking celerity versus normalized peak wave discharge for different values of the ratio nb /nc. The lines are not physically meaningful where Cbr > Cg
CE
As noted earlier, the jave-trigger simulation indicated that the ice cover was dislodged when the under-ice flow depth, Ybr, was not much less than the sustained peak depth Yp (ratio = 0.89).
AC
Jasek et al (2005) did not provide values of Yp and Ybr in the tabulation of their results but Yp was calculated by the writer assuming open-water flow, as indicated previously; while Ybr was determined from the equation of continuity (Eq. 3), which implies that
Ybr Yi (Qbr Qo ) / CbrW
(15)
27
ACCEPTED MANUSCRIPT Calculated values of Ybr were always slightly less than those of Yp while the ratio Ybr/Yp was approximately constant for a given slope, regardless of the value of Qbr (~0.97 and 0.90 for So = 0.00005 and 0.0003, respectively).
T
While the preceding quantitative results may be modified by different values of the Manning
IP
coefficients nb and ni, a reasonably well-defined picture of the SSW has emerged. The relatively
CR
sharp wave front enters the ice-covered region and lifts the ice cover, which does not break until the flow depth is a little lower than the sustained-peak depth. Upstream of the ice breaking
US
location, a transitional reach leads to the prolonged wave crest. Here, broken ice is carried along
AN
with the water but the Manning resistance equation is neither fully of the ice-covered nor of the open-water type (Appendix B). This effect likely dictates the increase in flow and depth over and
ED
M
above the respective ice-breaking values upstream of the “breaking” point.
The SSW has certain similarities with the well-known monoclinal wave (Ferrick and Goodman
PT
1998) but is unique in that its water surface profile combines open-water and ice-covered reaches
CE
which are joined by a hybrid reach. Appendix B contains a detailed discussion of this aspect and presents methodology for analytical determination of the profile of the wave front based on
AC
integration of the Saint-Venant equations of motion.
5. NATURAL STREAM CONDITIONS: PEACE RIVER CASE STUDY Since 1973, Alberta Environment (AENV) and BC Hydro have carried out detailed observations of ice formation and breakup throughout each ice season between the Bennett Dam and Vermilion Chutes (Fig. 13). The Chutes are located ~900 km below the Bennett Dam and some
28
ACCEPTED MANUSCRIPT 70 km downstream of Fort Vermilion. On rare occasions, observations have extended beyond Vermilion Chutes, all the way to the Mouth of Peace (MOP for short), which is also the head of the Slave River, and situated 1243 km from the dam. River slope generally decreases with distance from the dam. Jasek and Pryse-Phillips (2015) indicate slopes of 0.0003 between ~kms
T
300 and 500 (approx. Dunvegan to Sunny Valley), and 0.00005 between ~kms 700 and 900
AC
CE
PT
ED
M
AN
US
CR
IP
(approx. Tompkins Landing to Vermilion Chutes).
Figure 13. Map of Peace River main stem in BC and Alberta showing key locations and river kilometres from the W.A.C. Bennett Dam. From Jasek and Pryse-Phillips (2015), with changes.
The upstream advance of the ice front in the fall and winter and its recession in the spring are summarized graphically in a date-distance plot of the kind illustrated in Fig. 14. Ice front recession usually begins in late February or early March, slowly at first but accelerating as time
29
ACCEPTED MANUSCRIPT goes on. The rate of downstream advance of the front is typically less than ~ 50 km/day, while
ED
M
AN
US
CR
IP
T
there are 4 occasions where it exceeds 100 km/day over considerable distances.
PT
Figure 14. Observed locations of Peace River ice front downstream of Bennett Dam. The dashed green line has a slope of 50 km/day. From Alberta Environment archived material, with changes (accessed Dec. 19, 2016). http://www.environment.alberta.ca/forecasting/RiverIce/pubs/2013-
CE
2014_Peace_River_Ice_Obs_Rpr_No46.pdf
AC
With typical values of Q (3000 m3/s) and W (600 m) for the lower Peace River (slope ~0.00005), the rate of advance from Eq. 1 has been calculated as 8.6 km/day per degree C of water temperature. Observed rates of up to 50 km/day correspond to water temperatures of up to ~6°C and can generate storage-release flow enhancements of up to ~50% (Fig. 5). The increase in the speed of the ice front with the passage of time is consistent with the lengthening of the openwater reach upstream of the front. It thus appears that ice melt dominates the process of ice front recession. Note, however, that such rates of advance can also be experienced by intermittent 30
ACCEPTED MANUSCRIPT jamming and release; javes propagate much faster than tens of km per day but jamming may last long enough to significantly reduce the overall-average rate of advance. Not shown in Fig. 14 is an extremely rapid advance between Vermilion Chutes and MOP that occurred in 2014 and is
T
discussed in detail next.
IP
Unlike in the case of a prismatic channel, the resistance of the ice cover to dislodgement varies
CR
along a river owing to variable bathymetry, slope, and channel planform as well as variable ice thickness and strength (Beltaos 2013b). If a SSW encounters a highly resistant ice cover
US
segment, it may be unable to dislodge it and begin to attenuate as it moves under it. If this
AN
attenuation is small (e.g. because the resistant segment is short), the wave can resume its ice breaking action and again become self-sustaining. At the same time, the running ice that trails
M
the breaking front will be arrested at the edge of the resistant ice cover and a jam will be
ED
initiated. As backwater builds behind the jam, it is possible that the jam will release soon afterwards and produce a new jave that follows on the heels of the first one. Alternatively, the
PT
jam may remain in place for a long enough time for the first jave to be arrested somewhere
CE
downstream and form another jam. One would then encounter two jams at the same time, separated by a considerable stretch of open water. This configuration occurred in Peace River
AC
near Fort Vermilion in late April of 2014.
According to publically available AENV observation reports, ice breaking commenced near the Town of Peace River (TPR) on April 22. In the following days, the breakup progressed in a sequence of jams and releases, with major ice runs noticed on the 25th and 27th. In the morning of
31
ACCEPTED MANUSCRIPT April 28, ice conditions were as follows (archive accessed Dec. 19, 2016): (http://www.environment.alberta.ca/forecasting/RiverIce/PeaceRiverArchive1314.html)
Open water with brash ice above km 806; ice jam, km 806 to 830; open water, km 830 to 872;
T
ice jam, km 872 to 900 (approx. Vermilion Chutes); relatively intact ice cover below km 900.
IP
(The observation reports and associated diagrams indicate that Fort Vermilion is located at km
CR
832, which the writer assumes to be the site of the local WSC hydrometric gauge). By the afternoon of April 29, this configuration had changed to: ice jam, km 813 to 830; open water, km
US
830 to 878.5; ice jam, km 878.5 to 900; open water, km 900 to 915; brash ice, km 915 to 916;
AN
intact ice cover, km 916 to 935 and most likely beyond.
M
According to M. Jasek of BC Hydro (Pers. Comm. 2014) the jam upstream of Fort Vermilion
ED
released at 1815 h (Mountain Daylight Saving Time or MDT) on April 30. The record of the local gauge indicates that the water level began to decrease shortly after 1300 h (Mountain
PT
Standard Time or MST) and began to rise at 1716 h MST. [In what follows, time will be
CE
expressed as standard time and the MST designation will be implied, unless otherwise specified]. The decrease in stage (~0.6 m) was most probably caused by the negative wave generated by
AC
earlier release of the downstream jam. The travel time of this negative wave is estimated as 3.5 hours, based on the theoretical analysis of Henderson and Gerard (1981). Consequently, the release of the downstream jam would have occurred at about 1000 h. This estimate could be in error by an hour or so because relevant hydraulic properties of the reach between Fort Vermilion and Vermilion Chutes are unknown. On May 1, the river was observed by Mr. Jasek from the air and an ice run of variable surface concentration was seen all the way to the Mouth of Peace (M.
32
ACCEPTED MANUSCRIPT Jasek, Pers. Comm. 2014). This run was arrested by intact ice cover ~10 km into the Slave River; the resulting ice jam remained in place until May 9 and produced large-scale flooding of the Peace-Athabasca Delta. This kind of flooding is essential to the ecological maintenance of the higher-elevation lakes and ponds of the delta (perched basins) and can only be produced by ice
IP
T
jamming (Prowse and Conly 1998; Beltaos et al. 2006).
CR
The hydrometric gauge record at Peace Point for May 1 (Fig. 15) comprises data sets of different reliability, owing to ice effects on the primary and secondary orifice lines. The primary orifice
US
line was moved by ice on April 29, while the secondary line was likely pinched by ice in early
AN
morning of May 1. The primary record was corrected by WSC staff by applying a constant offset starting at 0825 h of May 1; earlier data were deemed unusable and thence not included in the
M
corrected record. The writer applied the same offset value as the one used by WSC to the earlier
ED
gauge readings and obtained the line designated to be of “low confidence”. The gauge heights indicated by this line are of unknown accuracy, but the pattern may be helpful in piecing together
PT
details of the breakup of the local ice cover. The flat line before ~0300 h May 1 indicates that
CE
the orifice had been tossed above the water level and was reading atmospheric pressure. At ~0320 h, it appears that a brief wave (wave 1) passed by the gauge site, re-submerging the
AC
orifice, which again began reading some water pressure, even though its actual elevation remained unknown.
33
ACCEPTED MANUSCRIPT 12 wave 2
wave 1
IP
T
10
orifice above water level
corrected primary water level; high confidence (WSC) corrected primary water level; low confidence secondary water level
US
8
CR
9
AN
Gauge height (m)
11
0:00
6:00
M
7
12:00
18:00
0:00
ED
May 1, 2014, MST
CE
PT
Figure 15. Hydrometric gauge record at Peace Point during the 2014 ice breakup.
Hourly WSC time-lapse images at Peace Point indicate that the sheet ice cover may have shifted
AC
during the night of Apr 30 to May 1, but remained stationary until at least 0701 h of May 1; at 0801 h rubble appears in the image and the ice sheet was likely in motion because its surface looks different from that of the previous photo (Fig. 16). The presence of intact sheet ice at 0801 h suggests that the cover had been dislodged recently; else it would have already been reduced to smaller slabs and blocks. Transport of ice blocks between the sheet ice cover and the near bank (lower photo, Fig. 16) might have caused the large water level fluctuations indicated by the lowconfidence line in Fig. 15 after 0700 h of May 1. 34
AC
CE
PT
ED
M
AN
US
CR
IP
T
ACCEPTED MANUSCRIPT
Fig. 16. Peace River ice conditions at the Peace Point gauge site, morning of May 1, 2014. Flow is from right to left. Images courtesy of Water Survey of Canada.
35
ACCEPTED MANUSCRIPT Full-width rubble appeared in the 0901 h image, while subsequent images suggested that the broken ice might have been running continuously at varying surface concentration, until 0400 h of May 2. Visual comparisons of ice surface elevations relative to the far river bank between successive images have indicated that the water level changed little during 0401 and 0701 h of
T
May 1, while there was a conspicuous rise between 0701 and 0801 h. The relative steadiness of
CR
IP
the water level prior to 0701 h in the photos is consistent with the low-confidence line in Fig. 15.
Since the photos indicate that ice was mobilized between 0701 and 0801 h, selection of 0730 h
US
would not be in error by more than ~30 min, even though it is likely that the actual time of
AN
mobilization was closer to 0800 h than to 0700. The estimated travel time between kms 900 and 1135 (Peace Point), would then be 21.5 h, indicating an average ice-breaking rate of 235/21.5/3.6
M
= 3.0 m/s. Available hydrometric data at Peace Point (Beltaos 2011) indicate a base velocity of
ED
0.75 m/s under the intact ice cover, rendering the ratio Cbr/Ui equal to 4.0. This value is compatible with the data points and line shown in Fig. 11; therefore it could have resulted from a
CE
PT
SSW.
Known open-water and ice-covered flow hydraulics in the vicinity of Peace Point (Beltaos 2011)
AC
were utilized next to calculate the base-flow as 3450 m3/s (Table 1), assuming a base gauge height of 9.0 m. The ice thickness on May 1 was estimated as 0.82 m with a keel of 0.75 m (based on a measured value of 0.92 m on April 8 and typical ablation rates through April). If the breakup at Peace Point was caused by a SSW, with peak stage of 11.70 m, the flow Qp can be estimated as follows: (a) subtract 0.5 m to account for rubble, assumed to be 1 m thick with a porosity of 0.5, to find “net” stage = 11.2 m; (b) calculate the difference from the stage at the
36
ACCEPTED MANUSCRIPT bottom of the ice cover as 11.2-(9-0.75) = 2.95 m; (c) multiply by Cbr (3.0 m/s) and by the local width of 680 m, to find the excess flow as 6020 m3/s (via continuity), for a total Qp of 6020+3450= 9470 m3/s (Qp/Qo = 2.75). At the beginning of May, the Manning coefficient of the ice cover is ~0.024, same as that of the river bed (Beltaos 2011); therefore nc ~ 0.024 and nb/nc ~
T
1, equal to the ratio used in the CRISSP runs. Therefore, the methodology outlined in Appendix
IP
B can be applied to determine the value of Qbr that would produce a SSW with peak flow of
CR
9470 m3/s. This operation resulted in Qbr = 8750 m3/s (Qbr/Qo = 2.54), while the predicted value
US
of Cbr is 3.4 m/s, in fair agreement with the observed value of 3.0 m/s.
AC
CE
PT
ED
M
AN
Table 1. SSW analysis of ice breaking at Peace Point Variable Value Water surface slope, So 0.000064 Manning coefficients, ni, nb, nc 0.024, 0.024, 0.024 Base flow, Qo 3450 m3/s 0.75 m/s, 6.76 m Base velocity, Ui, and flow depth, Yi Base gauge height 9.0 m Bottom of ice gauge height 8.25 m Peak and effective(1) gauge heights 11.7 m, 11.2 m Sustained peak flow, Qp 9470 m3/s Ice-breaking flow, Qbr 8790 m3/s Ice-breaking celerity, Cbr (observed) 3.0 m/s Ice-breaking celerity, Cbr (calculated) 3.4 m/s Expected duration of crest at Peace Point 9 hours (1) gauge height that would prevail if there were no moving ice
It is possible, however, that the ice run may have stalled at one or more locations between Vermilion Chutes and Peace Point, forming short-lived jams, which upon release would have generated large javes, propagating much faster than 3.0 m/s. Jave speeds of 8 to 9 m/s have been estimated on two separate occasions for Peace Point (Beltaos 2007), so that there would have been ample time for stoppage between the mornings of April 30 and May 1.
37
ACCEPTED MANUSCRIPT
To further explore this issue, the front profile of the SSW wave was calculated according to the methodology developed in Appendix B. Figure 17 shows that this front comprises a gradual rise that would have commenced very early in the morning of May 1; this is not consistent with the
T
available, albeit low-confidence, gauge data and with the time-lapse imagery (a sharp rise did not
IP
start until after 0700 h). It is more likely that the ice cover was dislodged by a large wave (wave
CR
2) of relatively high celerity, produced by the release of an upstream ice jam, possibly initiated at a sharp bend above Boyer Rapids (km ~1117), a known jamming site. If this is indeed what
US
occurred, wave 1 could have been the leading part of the original jave, which was “clipped” as
AN
the jam formed and water again went into storage. The breakup front arrived at the MOP at about 1310 h (according to time-lapse imagery, kindly provided to the writer by Dr. Daniel Peters of
M
Environment and Climate Change Canada). With the estimated time of ice mobilization at Peace
ED
Point (0730 h), the average ice-breaking celerity over the intervening 108 km works out to be 5.3 m/s, almost twice the calculated average value over the previous 235 km. For the entire travel
PT
reach (343 km), the average ice-breaking celerity is approximately equal to 343/(27.173.6) or
CE
3.5 m/s.
AC
The preceding interpretation is corroborated by the non-sustained nature of the water level peak (Fig. 17). Had the peak been produced by a SSW, it would have lasted much longer: using Eq. 10, the celerity of the tail of the SSW is calculated as 1.7 m/s, and thence the duration of the sustained crest is estimated to have been ~9 h, placing its end at about 1730 h on May 1, well beyond the right-hand border of Fig. 17. The observed water level variation, on the other hand, continues to decrease after the peak, attaining a value of 10.25 m by 1730 h.
38
ACCEPTED MANUSCRIPT
12 WSC corrected water level WSC secondary water level
IP
T
Calculated from SSW theory
CR
10
US
Gauge height (m)
11
8 3:00
6:00
9:00
12:00
ED
0:00
M
AN
9
May 1, 2014, MST
CE
PT
Figure 17. Calculated rising limb of SSW at Peace Point. Spatial wave configuration was transformed to temporal variation using the ice-breaking celerity Cbr; the starting time was selected so as to match the time of the peak. Predicted ice mobilization time and gauge height: 0720 h and 11.3 m, respectively.
AC
The preceding considerations suggest that the SSW theory does not match the details of how ice breakup was initiated at Peace Point and possibly elsewhere along the river. Nevertheless, the theory provides a helpful framework for examining breakup events in an “average” sense. Since the resistance to ice breaking varies along the river, the breaking front will speed up or slow down depending on local resistance characteristics; occasional brief jamming could still fit into this framework, so long as the average celerity of ice mobilization is not less than the minimum Cbr value indicated in Figs. 9 and 11. Consequently, the values of Qbr and Qp that are listed in 39
ACCEPTED MANUSCRIPT Table 1 should be viewed as magnitudes that would occur if breakup had been advancing at a constant celerity, equal to the average observed value. This may not reflect what actually happens at any particular site (such as Peace Point), but does provide a “feel” for average
T
conditions over a long travel reach.
IP
6. EFFECTS OF CHANGING THRESHOLD CONDITIONS ON STORAGE RELEASE
CR
It has been shown in the preceding sections that the ice-breaking threshold parameters have a controlling influence on storage release. Therefore, it is of interest to consider how climate-
US
induced changes to river ice regimes may affect storage release flows. For example, climate
AN
warming may bring about a significant reduction in the thickness of the ice cover of a particular river. This change would reduce the resistance to dislodgment and thence the value of the
ED
M
threshold flow Qbr.
Let QbrP and QbrF be the present and future values of Qbr, such that
QbrP > QbrF. Though it is
PT
likely that the spring flow hydrograph will also change, such change is ignored in order to focus
CE
on the threshold flow effects, assuming that other factors are equal. If the highest carrier flow of the spring runoff hydrograph (Qomax) is lower than both QbrP and QbrF, a thermal breakup event
AC
will result under both present and future conditions; the relatively small flow enhancement by storage release would be the same.
A second possibility is that of a mechanical event occurring under both present and future conditions (QbrF < QbrP < Qomax). In this case, Fig. 10 suggests that the flow enhancement due to storage release will be smaller under future than under current conditions. The near-linear
40
ACCEPTED MANUSCRIPT variations implied by the plotted data points indicate further that the difference will be equal to the slope of the applicable line, multiplied by (QbrP- QbrF). In the case of the flatter channel in Fig. 10, the slope of the line defined by the data points is ~0.55; if QbrP- QbrF = 1000 m3/s, the flow enhancement under future conditions will be 550 m3/s less than it is under present
T
conditions. This calculation applies so long as Qbr/Qo is within the range defined by the data
IP
points in Fig. 10. It is reasonable to expect, and has been shown in Appendix B, that the linear
CR
fits in Fig. 10 will eventually bend and become parallel to the equality line as Q br/Qo increases (Fig. B3). In the unlikely case that both QbrP and QbrF lie within the latter range, the difference in
AN
US
peak flows will be equal to QbrP- QbrF or 1000 m3/s in the present example.
A mixed situation arises when the carrier flow is intermediate between present and future
M
thresholds, i.e. QbrF < Qomax < QbrP. In that case, the breakup event would be thermal under
ED
present conditions and mechanical under future conditions. Since only mechanical breakup events can lead to major ice jams, the frequency of ice-jam flooding could increase in the future
PT
regime, despite the reduced ice thickness. Moreover, storage release under present conditions
CE
would be controlled by melt, rather than by ice breaking. Consequently, flow enhancement under future conditions (QF) would be considerably more than the “present” value (QP), as can be
AC
seen by comparing Figs. 5 and 10. Table 2 summarizes the effects of a reduced threshold flow on storage-release flow enhancement. Table 2. Effects of reduced ice-breaking threshold flow on breakup type and storage release Case Magnitude of carrier flow Breakup type Effect on flow number relative to ice-breaking enhancement due to Present Future thresholds storage release conditions conditions 1 Qomax < QbrF < QbrP thermal thermal no effect 2 QbrF < Qomax < QbrP thermal mechanical positive, large 3 QbrF < QbrP < Qomax mechanical mechanical negative, small ( QbrP- QbrF) 41
ACCEPTED MANUSCRIPT
Similar reasoning can be applied to situations where the threshold flow will increase, i.e. QbrP < QbrF. This could happen, for instance, if the fall season becomes wetter in the future and thence results in higher freezeup levels than under the present climate. This would increase
T
resistance to ice dislodgment and breaking in the spring and result in a higher value of Qbr. The
CR
IP
effects of this kind of change are summarized in Table 3.
M
AN
US
Table 3. Effects of increased ice-breaking threshold flow on breakup type and storage release Case Magnitude of carrier flow Breakup type Effect on flow number relative to ice-breaking enhancement due to Present Future thresholds storage release conditions conditions 1 Qomax < QbrP < QbrF thermal thermal no effect P max F 2 Qbr < Qo < Qbr mechanical thermal negative, large P F max 3 Qbr < Qbr < Qo mechanical mechanical positive, small ( QbrF- QbrP)
ED
Notable impacts of increased threshold flows are: (a) reduced frequency of ice jams and related
PT
flooding due to less frequent occurrence of mechanical breakup events; and (b) large reduction in flow enhancement for case 2 and small increase for case 3. This increase would typically be less
CE
than the difference in threshold flows, but could approach (QbrF - QbrU) if both QbrF and QbrU are
AC
larger than about 2Qo (Fig. B3).
An effect that is difficult to quantify at present pertains to case 3 in Table 3. Since QbrF > QbrP, ice breaking will be occurring at a higher water level and at lower speed under future than under present conditions (Fig. 9). Therefore, there will be increased attrition of the volume of brokenup ice that forms the ice run behind the ice-breaking front. Such attrition is primarily caused by stranding of ice blocks on islands and low banks, as well as by melt that may be occurring in 42
ACCEPTED MANUSCRIPT transit (Prowse 1986). Increased attrition will result in reduced ice volumes that are available to form ice jams when ice runs are eventually arrested. In turn, this would reduce the potential of the jam for large-scale flooding.
T
7. DISCUSSION
IP
Melting and sustained ice breaking have been identified as two processes that can enhance the
CR
flow during the breakup period. River slope is an important factor that has opposing effects, depending on the process. In the case of melting (Fig. 5), typical river slopes (e.g. 0.0005) seem
US
to produce just a slight flow enhancement, which becomes much more substantial in flat rivers
AN
(e.g. slope = 0.00005). In the case of ice breaking, higher slopes produce higher flow enhancement than do lower slopes (Fig. 10). Moreover, the flow enhancement caused by ice
ED
M
breaking is typically much larger than what is generated by ice melt.
A striking result of the CRISSP runs is that ice-breaking can produce SSWs with sustained peak
PT
flows well in excess of the peak of the incoming carrier flow hydrograph, up to ~60% higher
CE
according to Jasek et al (2005). It is not known whether peak amplification can occur when ice breaking is initiated by a jave, as is usually the case in flat rivers like the lower Peace or the
AC
lower Mackenzie (Beltaos and Carter 2009, Beltaos 2013a). The simulation with River1d indicated that the peak flow was reduced when the jave morphed into a SSW. Javes typically encounter intact ice cover as soon as they are generated because jams are held in place by intact ice. At that time the jave front is almost as high and steep as the toe region of the parent jam (Fig. 7). If the jave does break the ice cover, a SSW may eventually form but it seems unlikely
43
ACCEPTED MANUSCRIPT that its front can be as high and steep as that of the jam toe. This conjecture is, of course, subject to testing with further numerical modelling.
As noted in Section 4, Jasek et al (2005) found that enhancing the value of nb upstream of the
T
wave front in order to simulate possible impedance to flow by the moving broken ice led to
IP
implausible results. Consequently, this effect has not been considered herein. Since front celerity
CR
exceeds water velocity, as can be easily ascertained using Eq. 13, ice congestion is unlikely to
US
occur. Moreover, the broken ice is floating at a higher elevation than that of the intact ice cover and therefore is less confined by the river banks; any resistance that may develop near the river
AN
banks is likely to be minimal, relative to bed friction, which is applied throughout the channel width. On the other hand, the thickened ice pack from the released jam moves at a lower
M
elevation that than of the jam; it could conceivably generate significant side resistance before
ED
dispersion reduces its thickness to that of a single ice block. However, the moving pack trails
PT
well behind the wave front and does not affect front properties. These considerations are consistent with the findings of She and Hicks (2006), which showed that moving-ice impedance
AC
CE
is of secondary importance in jave propagation, especially near the wave front.
An important point that should be stressed here is that flow enhancement due to storage release is a highly transient phenomenon that ceases when the recession of the edge of the ice cover is discontinued. For instance, ice breaking will cease when a sufficiently resistant ice cover is encountered by the rising flow and a jam begins to form. Downstream of the jamming point (toe of jam), the wave will continue to advance but keep attenuating unless ice breaking resumes. Upstream of the toe, water levels will rise rapidly, creating new storage and thence reducing 44
ACCEPTED MANUSCRIPT downstream flows. When all of the moving ice comes to a stop, upstream and downstream flows will return to the value of the carrier river flow, which therefore governs the potential of the jam to cause flooding over a sustained period of, say, a day or more.
T
It follows that the main effect of storage release is in facilitating the occurrence of a mechanical
IP
event when the carrier flow in the river is lower than the threshold value. In the case of melt, this
CR
eventuality seems very unlikely in steeper channels but conceivable in flat ones, for which meltinduced flow enhancement can be sizeable (Fig. 5). This is particularly probable where the
US
threshold flow is not much greater than the carrier flow. In the case of ice breaking, this
AN
reasoning may appear circular for prismatic channels: if the carrier flow cannot initiate breakup at the beginning of the ice covered reach, how can it be enhanced at any location farther
M
downstream? In natural streams, however, one may envisage circumstances where the carrier
ED
flow can actually dislodge the ice cover in a low-resistance section of the river. The resulting enhanced flow could then exceed the threshold farther downstream, where the resistance to
CE
PT
dislodgment may be higher.
Climatic warming could reduce or increase the values of ice-breaking threshold flows, depending
AC
on basin-specific changes, such as thinner ice covers or higher freezeup levels. In turn, such changes can have significant effects on storage release. Using the analytical methodology developed in previous sections and in Appendix B, it was possible to assess the effects of changes in threshold flows on storage release contributions to the flow hydrograph. Such changes might also result from regulation; assessment of their impacts could also be made using the summaries presented in Tables 2 and 3.
45
ACCEPTED MANUSCRIPT
8. SUMMARY AND CONCLUSIONS Two mechanisms of storage release, ablation and ice breaking, have been examined in details and quantified. Melting of the ice cover can significantly enhance the flow, especially in
T
relatively flat channels. Excess flows may amount as much as 40% for water temperatures up to
IP
5°C. Storage release by ice breaking is much more dynamic and results in considerably greater
CR
excess flows than ice melt. Numerical modelling and analysis for prismatic channels suggest that ice breaking, once initiated, morphs into a SSW with constant front shape and celerity. The
US
properties of the SSW do not appear to depend on the “trigger” flow hydrograph or jave
AN
properties but are largely determined by the threshold discharge, as well as base flow, ice, and channel properties. Ice breaking celerities can be large if the threshold discharge is very close to
M
the base flow, but decrease to a minimum value of ~ 4 x (under-ice base velocity) as the
ED
threshold-to-base flow ratio increases. Application of the analytical and numerical results to the lower Peace River indicated that ice melt could explain typical observed rates of ice recession;
PT
the SSW theory was only partially applicable to the very rapid, and rare, ice clearance that was
CE
observed in the spring of 2014. Ice jamming at high ice resistance sites that are invariably encountered in natural streams can disrupt the orderly advance of the wave. Nevertheless, it was
AC
concluded that the SSW concept can be useful in terms of the overall average advance of highly dynamic breakup events. The effects of changing ice regimes, e.g. due to climate change or regulation on storage release, can be assessed in terms of changes in threshold ice-breaking flows by simple reasoning that was illustrated in two examples. It has been stressed that flow enhancement by storage release is highly transient and its main role is in facilitating ice breaking farther downstream.
46
ACCEPTED MANUSCRIPT
9. ACKNOWLEDGMENTS Support by Environment and Climate Change Canada to carry out this work is gratefully acknowledged. I would also like to thank: Angus Pippy and his colleagues at the WSC
T
Yellowknife office for supplying archived, but unpublished, hydrometric data and interpretation;
IP
Professor Faye Hicks (ret. University of Alberta) for making available the model River1d;
CR
Martin Jasek of BC Hydro for sharing observational information regarding the 2014 Peace River breakup; Nadia Kovachis Watson of Alberta Environment for sharing historical data on ice cover
US
recession rates in Peace River; and my colleague Dr. Daniel Peters for sharing 2014 time-lapse
AN
imagery at the mouth of Peace River .
ED
M
APPENDIX A. USING INCREASED BED MANNING COEFFICIENT TO SIMULATE ICECOVERED FLOW In the present applications of the River1d model, the open-water option has been used for
PT
simplicity, while the presence of a stationary ice cover of constant thickness and roughness has been simulated by locally enhancing the value of the bed Manning coefficient. Using the
CE
Belokon-Sabaneev formula (Eq. 6), it can be shown that this assumption is equivalent, in terms
AC
of flow depth and mean velocity, to having an actual ice cover of which the roughness coefficient (ni) is given by:
ni (nbE 3/ 2 nb3/ 2 )2 / 3
(A1)
with nb = actual roughness coefficient of the river bed; and nbE = enhanced bed roughness value used to simulate the “ice covered” reach. This simplification facilitates the simulation of 47
ACCEPTED MANUSCRIPT sequential ice breakup, which comprises successive runs of River1d using the output computed at the end of any one time step as input for the next time step but with a shortened “ice cover”. Rather than repeatedly modify and re-load the ice file, the open-water reach (lower value of nb)
T
is simply extended to the new location of the ice edge.
IP
APPENDIX B. FRONT PROFILE OF THE SELF-SUSTAINING WAVE
US
Y 1 U U U x g t g x
(B1)
AN
S f So
CR
The one-dimensional momentum equation reads (Henderson 1966):
M
in which U = flow velocity; and Sf = friction slope, defined as total boundary resistance force per
ED
unit channel length divided by flow area = average boundary shear stress divided by hydraulic
(B2)
CE
2
AC
nU S f 2/3 R
PT
radius. From the Manning resistance equation, Sf can be expressed in terms of U, i.e.
Here n = nc, R = Y/2 for ice-covered flow while n = nb, R = Y for open-water flow that exhibits the well-known logarithmic velocity distribution in the vertical direction.
With reference to the sketch of Fig. B1, three flow regions can be identified: (a) essentially openwater flow upstream of point A, with freely moving ice blocks on the surface; (b) ice-covered flow downstream of the breaking site (point B); and (c) transitional reach (AB) of “hybrid” flow. 48
ACCEPTED MANUSCRIPT Here, the surface velocity, which is zero just downstream of point B, will increase gradually in the upstream direction, from zero to that of the open-water reach at point A. Therefore the Manning resistance relationship, from which Eq. B2 derives, would not be accurate if the
ED
M
AN
US
CR
IP
T
hydraulic radius were set equal to Y or to Y/2.
CE
PT
Figure B1. Schematic illustration of vertical velocity distributions in different portions of the SSW front.
R
AC
To compute the transitional segment of the wave profile, the hydraulic radius is expressed as
Y 1 e x / L
(B3)
where x is distance measured from point B in the upstream direction. The hydraulic radius is equal to Y/2 at point B and asymptotically tends to Y as x increases. The length L is an empirical length scale, characterizing the extent of the transitional reach. The hydraulic radius 49
ACCEPTED MANUSCRIPT becomes practically equal to Y for x = 4L (denominator = 1.02), which can be considered to be the length of the segment AB in Fig. B1.
Introducing the transformations = x-Cbrt and = x-Cbrt (for the reaches BC and AB,
(reach BC);
dY d
S f So (C U ) 2 1 br gY
CR
So S f (C U ) 2 1 br gY
(reach AB)
(B4)
US
dY d
IP
T
respectively), and using also the equation of continuity (Eq. 3), Eq. B1 can be re-arranged to read
AN
Moreover, the flow velocity, U, can be obtained from the equation of continuity (Eq. 3), which is
M
first applied between any point along the front (Q = UYW) and the undisturbed flow downstream
PT
(B5)
CE
Y U C C 1 i Ui Ui Ui Y
ED
of the front (Qo = UiYiW), and then re-arranged to read:
AC
An obvious requirement for the profile of the wave front is that Y should decrease with increasing , i.e. dY/d < 0. For the present conditions, it can be shown that So-Sf < 0; therefore, the denominator of the first of equations B4 should be positive, which will be satisfied for
Cbr U i gYi Cg
(B6)
50
ACCEPTED MANUSCRIPT The quantity on the right-hand-side represents the celerity of a gravity wave (Cg) propagating in the ice-covered reach downstream of the SSW front.
The ordinary differential Eqs. B4 can be easily integrated in Excel, once So, Qo, W, Qbr, ni, nb,
T
are given and a trial value of Qp is selected. With the selected Qp, one can then compute Cbr from
IP
Eq. 14 and Ybr from Eq. 15; and subsequently compute Ubr as Qbr/YbrW and thence determine
CR
(Sf)br via Eq. B2 with x = 0. The computation of the profile is straightforward in the ice-covered
US
reach [with known (Sf)br determine dY/d and thence the value of Y at a small increment of downstream; calculate new values of U, Sf and dY/d; repeat with equal increments of , until Y
AN
no longer changes, becoming equal to Yi]. For the transitional reach, the same algorithm is
M
applied but in the upstream direction; here, the value of the hydraulic radius R is no longer equal to Y/2 but calculated from Eq. B3. The selected value of Qp is adjusted until the solution matches
ED
asymptotically the value of Yp that is obtained from the hydraulic resistance equation (Eq. 12).
PT
An example is shown in Fig. B2, where the wave front is seen to be slightly under 10 km long. This example corresponds to CRISSP Run 9 of Jasek et al. (2005), for which the development of
CE
the SSW was illustrated graphically; the length of the wave front was scaled off their
AC
corresponding figure and estimated as ~10 km.
51
ACCEPTED MANUSCRIPT 4.2 ice-covered reach
4
transitional reach
T IP
3.6
CR
4 3.4
3.2
0 -5000
US
2
0
3 -10000
5000
M
-5000
AN
Flow depth (m)
3.8
0
5000
ED
Channel distance (m)
PT
Figure B2. Wave profile for steeper channel (So = 0.0003), with W = 600 m, Qbase = 1600 m3/s, Qbr = 2400, ni=nb=0.03 and L ~ 1300 m. Ice breaking occurs at distance zero. The inset allows direct visual comparison of breaking and peak depths.
CE
The results of the CRISSP runs (Jasek et al 2005), were closely reproduced analytically (Fig. B3)
AC
using a constant value of L for each channel slope, i.e. ~3000 and ~1300 m for So = 0.00005 and 0.0003 respectively. For interpolation purposes, normalized values of L can be used, e.g. L/Yi = 529 (So = 0.00005) and 393 (So = 0.0003).
No solution to Eq. B4 could be found for the ice-covered portion of the profiles in the two instances (one for each slope) in which Cbr exceeded the value of Cg. These two profiles would be very short as they involve minute differences between the depths Ybr and Yi (less than 52
ACCEPTED MANUSCRIPT ~0.1%). It is possible that the discretization of the numerical model was too coarse to account for the equation of momentum along this short reach. The constraint imposed by Cg did not permit application of the analytical method to very small values of Qbr/Qo-1 (Fig. B3).
CR
IP
T
2.0
US
1.0
AN
QP / Qo - 1
1.5
M
CRISSP runs - bed slope = 0.0003
0.5
CRISSP runs - bed slope = 0.00005
ED
Analytical solution - bed slope = 0.0003
0.2
0.4
0.6
Equality line
0.8
1
1.2
1.4
Qbr / Qo - 1
CE
0
PT
0.0
Analytical solution - bed slope = 0.00005
AC
Figure B3. Results of present methodology versus CRISSP-generated relationships between sustained peak and threshold flows.
REFERENCES Ashton, G.D. 1986. River and Lake Ice Engineering. Water Resources Publications, Littleton, CO., USA.
53
ACCEPTED MANUSCRIPT Beltaos, S. 2007. The role of waves in ice-jam flooding of the Peace-Athabasca Delta. Hydrological Processes, 21(19), 2548-2559. Beltaos, S. 2009. River flow abstraction due to hydraulic storage at freezeup. Can. J. Civ. Eng. 36: 519-523.
T
Beltaos, S. 2011. Developing winter flow rating relationships using slope-area hydraulics. River
IP
Research and Applications. 27(9): 1076-1089.
CR
Beltaos, S. 2013a. Hydrodynamic and climatic drivers of ice breakup in the lower Mackenzie River. Cold Regions Science and Technology 95: 39-52.
US
Beltaos, S. 2013b. Hydrodynamic characteristics and effects of river waves caused by ice jam
AN
releases. Cold Regions Science and Technology 85: 42–55.
Beltaos, S. 2016. Extreme sediment pulses during ice breakup, Saint John River, Canada. Cold
M
Regions Science and Technology, 128: 38-46.
ED
Beltaos, S., and Burrell, B.C. 2005a. Field measurements of ice-jam-release surges. Canadian Journal of Civil Engineering, 32(4), 699-711.
PT
Beltaos, S., and Burrell, B.C. 2005b. Determining ice-jam surge characteristics from measured
CE
wave forms. Canadian Journal of Civil Engineering, 32(4), 687-698. Beltaos S. and Burrell, B.C. 2016. Transport of suspended sediment during the breakup of the ice
AC
cover, Saint John River, Canada. Cold Regions Science and Technology, 129: 1-13. Beltaos, S. and Carter, T. 2009. Field studies of ice breakup and jamming in lower Peace River, Canada. Cold Regions Science and Technology Journal. 56(2-3): 102-114.
Beltaos, S., Prowse, T.D. and Carter, T. 2006. Ice regime of the lower Peace River and ice-jam flooding of the Peace-Athabasca Delta.
Hydrological Processes, Volume 20, Issue 19,
pp 4009-4029.
54
ACCEPTED MANUSCRIPT Burrell, B.C. 2008. Chapter 1. Introduction. In: River Ice Breakup, Water Resources Publications, Highlands Ranch, Co., USA., pp 1-19. Carey, K.L. 1966. Observed Configuration and Computed Roughness of the Underside of River Ice, St. Croix River, Wisconsin, U.S. Geological Survey Professional Paper 550-B, B192-
T
B198.
IP
Carey, K.L. 1967. The Underside of River Ice, St. Croix River, Wisconsin, U.S. Geological
CR
Survey Professional Paper 575-C, C195-C199.
Ferrick, M.G., and N. D. Mulherin, 1989. Framework for control of dynamic ice breakup by river
US
regulation. Journal of Regulated Rivers: Research & Management, 3, pp. 79-92.
AN
Ferrick, M.G. and Goodman, N.J. 1998. Analysis of Linear and Monoclinal River Wave Solutions. USACE Cold Regions Research and Engineering Laboratory, CRREL Report 98-1,
M
Hanover, NH, USA, 24 p.
ED
Gerard, R. 1990. Hydrology of floating ice. In: Northern Hydrology, Canadian Perspectives, Prowse TD, Ommanney CSL (eds). NHRI Science Report No. 1, National Hydrology
PT
Research Institute, Environment Canada, Saskatoon, 103–134, plus references.
CE
Gerard, R., Kent, T.D., Janowicz, R., and Lyons, R.O. 1984. Ice regime reconnaissance, Yukon River, Yukon. Proceedings of the 3rd International Specialty Conference on Cold Regions
AC
Engineering, 4–6 April 1984, Edmonton, Alta. Compiled and edited by D.W. Smith. Canadian Society for Civil Engineering, Montréal, Quebec, pp. 1059–1073. Henderson, F.M. 1966. Open channel flow. The Macmillan Co., New York, Toronto. Henderson, F.M., and Gerard, R. 1981. Flood waves caused by ice jam formation and failure, Proceedings IAHR Symposium on Ice, Quebec, Canada, Vol. I, pp. 277-287.
55
ACCEPTED MANUSCRIPT Hicks, F., Cui, W, and Ashton, G. 2008. Chapter 4. Heat transfer and ice cover decay. In: River Ice Breakup, Water Resources Publications, Highlands Ranch, Co., USA., pp 67-123. Jasek, M. and Pryse-Phillips, A. 2015. Influence of the proposed Site C hydroelectric project on the ice regime of the Peace River. Canadian Journal of Civil Engineering, 2015, 42(9): 645-
T
655.
IP
Jasek, M., Ashton, G., Shen, H.T., and Chen, F. 2005. Numerical Modeling of Storage Release
CR
during Dynamic River Ice Break-up. Proceedings (CD-ROM) of 13th Workshop on the Hydraulics of Ice Covered Rivers, Hanover, NH, September 15-16, 2005, CGU-HS
US
Committee on River Ice Processes and the Environment, Edmonton, Canada, 421-439.
AN
Kellerhals, R., Neill, C.R. and Bray, D.I. 1972. Hydraulic and geomorphic characteristics of rivers in Alberta. River Engineering and Surface Hydrology Report 72-1, Research Council of Alberta,
M
Edmonton, Canada, 52 p.
ED
Parkinson, F.E. 1982. Water temperature observations during break-up on the Liard-Mackenzie River system. Workshop on the Hydraulics of Ice Covered Rivers, Edmonton, Canada, 261-295.
PT
Ponce, V.M., Simons, D.B., 1977. Shallow wave propagation in open channel flow. Journal of
CE
the Hydraulics Division, ASCE 103 (HY12), 1461–1476. Prowse, T.D. and Carter, T. 2002. Significance of ice-induced storage to spring runoff: a case
AC
study of the Mackenzie River. Hydrological Processes 16(4): 779–788. Prowse, T.D. 1986. Ice jam characteristics, Liard-Mackenzie rivers confluence. Canadian Journal of Civil Engineering, 13(6), 653-665. She, Y. and Hicks, F. 2006. Modeling ice jam release waves with consideration for ice effects. Cold Regions Science and Technology 45 (2006) 137–147.
56
ACCEPTED MANUSCRIPT HIGHLIGHTS
-Storage release during river ice breakup can generate large additions to flow -It can be caused by ablation or by ice breaking
T
-Hydrodynamic properties of each process studied analytically and numerically
IP
-Ice breaking is far more dynamic and may lead to formation of self-sustaining waves
AC
CE
PT
ED
M
AN
US
CR
-Results applied to a case study and to assess climate impacts on ice-jam flooding
57