Extension and refinement of a stochastic modelling approach to assess ice-jam flood hazard

Data from recent extreme events and the additional data available since the last assessment was carried out will lead to adjustments in the frequency distributions of flows and water levels since such distributions are sensitive to the length of their time series and the severity of the extreme events within the series (Hu et al. 2020). This not only affects flows and water-level frequency distributions for open-water floods but also affects flow and stage frequency distributions of ice-induced flooding of northern rivers. One such example is the town of Peace River along the Peace River in western Canada which responded to the 1992 and 1997 ice-jam events, the two largest events on record, by raising dike crest elevations and inserting extra flood defense walls in keeping with the revised design flood threshold levels (although it should be noted that the dike elevations were designed to a 1:100 AEP (annual exceedance probability) open-water flood level, not to ice-jam flood levels) (Lindenschmidt et al. 2016). In the town of Peace River, ice-jams lead to more severe floods than those of open-water alone; hence, the probability of overtopping is higher from an ice-induced flood than an open-water flood. Another example where new flood hazard assessments are periodically carried out is the town of Fort McMurray along the Athabasca River, also in western Canada, where the most extreme flood levels are also induced by ice-jams.

To capture the chaotic nature of ice-jams, the author has developed a stochastic modelling approach to mimic the probabilistic nature of ice-jam flooding. In this approach, a deterministic (repeatable output for a single set of inputs) river ice hydraulics model is embedded in a Monte-Carlo framework which allows many simulations of the model to be executed automatically. Each simulation has a different set of input values for the parameters and boundary conditions, which are extracted randomly from frequency distributions. The Monte-Carlo runs yield an ensemble of backwater level profiles, from which profiles of floodwater level exceedance probabilities can be extracted. The methodology was first developed for the assessment and mapping of flood hazards and risks for the town of Peace River (Lindenschmidt et al. 2015, 2016) and was then extended to forecast ice-jam flooding along the Athabasca River (Lindenschmidt et al. 2019), lower Red River (Williams et al. 2021), and Saint John River in eastern Canada (Das et al. submitted). The method has been applied operationally to forecast ice-jam flooding along the lower Churchill River in Labrador (Lindenschmidt et al. 2020, 2021) and the Exploits River in Newfoundland (Warren et al. 2017). Applying the method for flood hazard and risk assessments proved to be useful in evaluating the best options to mitigate ice-jam floods in Fort McMurray (Das & Lindenschmidt 2021a, 2021b, 2021c).

The method required refinement so that the water-level exceedance probabilities from the modelled ensembles could attain the new and higher elevations of the return periods calculated from the gauged and surveyed data. Refinements included the following:

• Relating the volume of ice constituting ice-jams to the river discharge; there can be much scatter in this relationship (Beltaos 2021); a confidence band of the volume rate of ice vs. discharge helps attain more extreme ensemble backscatter level profiles when the flooding effect of high flows is compounded by higher volumes of ice.

• In the previous development, due to lack of data, the locations of the ice-jam lodgements were assumed to be uniformly distributed along the river chainage; more data has allowed lodgement locations to be concentrated in certain reaches where more ice-jams have occurred in the past.

• Setting the water level at the downstream boundary condition of the model was made to be a function of the upstream discharge.

The three important objectives are as follows:

• include historical data in the stage frequency analysis to capture the severity of these events in the analysis;

• introduce refinements to the stochastic modelling approach to help attain the level of flood severity reflected in the frequency analysis of the observed data; and

• demonstrate the refined modelling framework's applicability by simulating scenarios of artificial breakage for ice-jam flood mitigation.

Ice-jamming is a regular occurrence along the Athabasca River during the spring ice-cover breakup with some of the most severe flooding having occurred in the springs of 1979, 1997, and 2020, as labelled in the water-level elevation vs. discharge plot in Figure 2. The most recent extreme event of 2020 caused insured damages amounting to over a half billion dollars (IBC 2020).

Rubble ice describes the ice floes arriving at the ice-jam front from upriver to supplement the volume of ice in the jam and extend the length of the ice-jam cover. This ice can be thicker than the intact (thermal) ice-cover thickness h downstream of the ice-jam since it may stem from consolidated ice covers, which is often the case for ice formed along the Athabasca River upstream of Fort McMurray. Hence, the front FT of the ice-jam cover can be quite thick. The volume of ice is expressed as the rate of running ice (with voids characterised by porosity PS) that arrives at the ice-jam from upriver. This volume does become the volume of ice within the ice-jam cover. The simulation stops when the ice-jam front reaches a user-specified location, generally just downstream of the upstream boundary.

Measured and estimated values of ice-jam porosities are reported by White (1999, Table 4 and Figure 12) and can vary widely, between approximately 0.2 and 0.9. There is a tendency for higher porosities to be found in frazil-deposited ice-jams compared to breakup ice-jams and the values can vary greatly along the length and thickness of the ice-jam (White 1999). Values of porosity are generally smaller (between 0.35 and 4.5) for shoved sections of the ice-jam but can have a wider range for hanging dams (between 0.24 and 0.89) (White 1999, p. 16–17). Ice porosities for the 1977 and 2020 ice-jam events were calibrated somewhat higher when air temperatures dropped below freezing potentially allowing frazil ice to have been formed and accumulated in the ice-jams.

Once a set of parameter and boundary condition values has been extracted from each of their corresponding frequency distributions or functions, a simulation is carried out. The process was repeated 1680 times (Monte-Carlo analysis labelled [vi] in Figure 7) to yield an ensemble of backwater level profiles ([vii] in Figure 7). The high number of runs is deemed sufficient to provide statistics on the backwater profiles that can be summarised using 0.5, 1 and 2 percentile profiles, which correspond to the 1:50, 1:100, and 1:200 AEP profiles ([viii] in Figure 7). The elevations of the percentile profiles are compared to the return periods calculated from the instantaneous water-level maxima observed at the Clearwater River mouth and the Athabasca River (WSC) gauge ([ix] in Figure 7). If they do not coincide, the confidence level of the V_ice = f(Q) relationship ([ii] in Figure 7) is adjusted and the Monte-Carlo process is repeated until the ensemble percentiles ([viii] in Figure 7) and the observed return periods ([ix] in Figure 7) coincide. To gain confidence in the results, the process is repeated several times, with each time outputting a slightly different set of percentile profiles, which are then averaged to yield the final profile. At this point, the framework is deemed calibrated.

For the 2015 ice-jam event, the ice-jam toe was reported to be 12 km downstream from the bridge (chainage = 36,600 m). However, the highwater mark at the WSC gauge is more than 3 m lower than the highwater marks at the Clearwater River mouth and the bridge. With such a large difference in these highwater marks, the model could only be calibrated assuming that there was a second ice-jam toe further upstream between the mouth and the gauge. Having two lodgements within an ice-jam is not uncommon for the Athabasca River as the 16 April 1977 ice-jam event demonstrates (Doyle 1977). Much more ice accumulated at the upstream lodgement than at the downstream lodgement.

For the 2020 ice-jam event, the 2,700 m³/s flow at the upstream boundary of the model is 400 m³/s more than that used for a simulation of the ice-jam by H&G (2022) using HEC-RAS, elaborated in Nafziger et al. (2021). The modelling carried out by H&G (2022) used highwater marks to calibrate their mode; two data sets, highwater marks and the extent of the ice-jam cover were used to calibrate the modelling in my study. Because the location of the ice-jam head and toe are known, the ice-jam volume, the ice roughness, discharge, and porosity of rubble ice, and the ice-jam become important calibration variables for model calibration. The volume of ice within the ice-jam was simulated to be 31 million m³, which is in line with H&G (2022)’s estimate of the ice volume being greater than 28 million m³.

The data from the Clearwater River mouth include historical data, the largest event being from 1875, which has been estimated to have reached an elevation of 252 m a.s.l. (Hatch 2017). Figure 10 shows the observations for both recent and historical data. The programme HQ-EX 3.0 from the consulting firm DHI WASY GmbH in Berlin, Germany was used to construct the stage frequency distributions with historical data, which uses an algorithm of repeated blocks of the recent time series to fill in the gaps between the recent time series and historical events. Distributions available in the software include Gumbel, generalised extreme value, Rossi, log-normal, Pearson type 3, log-Pearson type 3 and Weibull fitted to plotting positions using the method of moments, maximum likelihood, and mixing distribution. For the analysis with historical data, gaps within the recent time series are not allowed; hence, the analysis was carried out with the recent time series of 38 consecutively years falling within two 38-year time frames, an earlier one of 1962–1999 and a later one of 1983–2000. The parameters of the distributions – location, scale, and shape – of the two analyses were averaged to obtain the resulting stage frequency distribution.

A total of 1,680 model runs were carried out using the Monte-Carlo framework, to yield an ensemble of 1,680 backwater level profiles. The 0.5, 1, and 2 percentile profiles, respectively corresponding to the 1:200, 1:100, and 1:50 AEP, are plotted in Figure 12 for a short reach of the model domain extending from the Clearwater River mouth to the Athabasca River gauge. The profiles overlay the return periods from the observed data, indicating a well-calibrated stochastic modelling framework and represent the flood hazard status quo for Fort McMurray. For reference, the 2020 ice-jam profile is also included and lies somewhat higher than the 1:50 AEP level at the Clearwater River mouth, which coincides with what is reported in Nafziger et al. (2021), that ‘the 2020 ice-jam produced water levels between the 1:50 and 1:100 AEP ice-jam water levels at Fort McMurray’. An important variable in the modelling framework that was fine-tuned to improve the calibration was the confidence level of the confidence band of the volume rate of ice vs. flow data (V_ice vs. Q), indicated by [ii] in Figure 7, which was set from 95% initially, to a final setting of 85%. The narrower confidence band was required to constrain the scatter between V_ice and Q. Also, concentrating more ice-jam toe locations in the stretch immediately downstream of the Clearwater River mouth helped to obtain more extreme backwater level profiles to raise the percentile profiles to the corresponding observed return period elevations.

An important assumption for these artificial breakage scenarios is that there is no possibility of ice-jamming occurring within broken-up stretches. This is not always the case, however, and experience on the Red River has shown that ice-jams do occur even in stretches of the river which were previously broken up, as described in Lindenschmidt (2020). Another assumption is that the volume rate of ice distribution remains the same as for the status quo. This may not be the case, however, for it is not known (i) how much of the extra water exposed to the freezing air temperatures, after the ice is broken up, freezes to create additional ice that can accumulate in ice-jams, or (ii) how much melting of the ice rubble occurs when more ice surfaces from the broken-up pieces come in contact with the water during artificial breakage. Hence, no change in V_ice is assumed for the scenario simulations.

For the 8-km breakage scenario, additional drops in backwater levels were simulated at the Clearwater River mouth (see bottom panel of Figure 16). The additional drops are not as substantial as the initial drop incurred by the 4-km-only breakage scenario (see Figure 17). Again, there are elevation increases at the gauge but not as high for the 8-km breakage scenario as they were for the 4-km breakage scenario. This again reveals the shifting of more extreme flood hazard further downstream as artificial breakage is carried out further downstream. It appears that longer stretches of breakage downstream of the Clearwater River mouth move the toe location and the extreme ice-jam backwater levels further downstream.

An amphibex machine has difficulty breaking ice thicker than ∼0.6 m, as is the case near Fort McMurray; moreover, the Athabasca River is much wider than the Red River. Depending on how many machines are deployed, it might require several weeks to break 4 or 8 km of ice. During this time, below-zero air temperatures may cause ice fragments to freeze together again. Altogether, it is not clear that ice breaking is a feasible mitigation measure near Fort McMurray, an important point from this study is that the updated MOCA model is capable of quantifying the benefits of ice breaking in situations where ice breaking would be feasible.

Refinements to a stochastic modelling framework were successfully applied to capture the severity of the 2020 ice-jam flood in Fort McMurray. The refinements included establishing a dependency between the volume of ice accumulating in ice-jams and upstream discharge. Concentrating the ice-jam toe locations in areas where more ice-jams have occurred in the past helped to attain the extreme exceedance flood levels reflected in the return periods of the observed water-level elevations. Including historical data in the frequency analyses provided a better description of the return periods. The framework allowed scenarios to be run for mitigating ice-jam flooding by artificial breakage of the ice cover downstream of the Clearwater River mouth. There is a substantial drop in the flood hazard at the Clearwater River mouth when 4 km of the reach's ice cover broke up. Further breakage of the cover seems to have only marginal additional impact in reducing the flood hazard in Fort McMurray. At the gauge, however, the flood hazard increases, with the increase becoming less as more ice is broken up further downstream. This shows the shift of flood hazard further downstream along the Athabasca River as more of its ice is artificially broken downstream.

The stochastic framework is a means to evaluate ice-jam flood risk reductions due to the artificial breakage mitigation option. However, there is no guarantee that ice-jams will not lodge within the broken-up areas and a probability of lodgements occurring within those areas could be integrated within the modelling framework, a topic of future work. The cost should also include the likelihood of ice-jams still occurring in the broken-up areas of the ice cover.

I would like to thank Jennifer Nafziger from Alberta Environment and Parks for providing the new cross-sections of the Athabasca River and corrected gauge data.

All relevant data are included in the paper or its Supplementary Information.

The authors declare there is no conflict.