Water Quality and Flow Management Scenarios in the Qu’Appelle River–Reservoir System Using Loosely Coupled WASP and CE-QUAL-W2 Models

Abstract

The water bodies of the Saskatchewan Prairies suffer multiple stressors, and demand for water is expected to increase. Water quality models can help evaluate water management strategies and risks such as climate change. This study assesses the impact of interbasin water transfers on the water quality of a strategic, eutrophic prairie reservoir that receives poor-quality watershed run-off. A one-dimensional WASP model was used to estimate nutrient transformations in the transfers along a 97 km river channel. The WASP model was then loosely coupled to a two-dimensional CE-QUAL-W2 model of the downstream receiving reservoir. Output from the WASP model was manually transformed into boundary conditions for the CE-QUAL-W2 reservoir model. This method improves on an earlier attempt to estimate nutrient transformations in the transfers using linear regression. Results from the loosely coupled models suggest Buffalo Pound Lake would respond well to the interbasin transfers. The number of threshold exceedances decreased for all modeled water quality variables. Nutrient concentrations were most influenced in the open water season following spring freshet. Any additional reduction of threshold exceedances during winter was minimal in comparison. These results are interesting from a management perspective as increased transfers under winter operations risk ice damage to the river channel.

1. Introduction

The river–reservoir networks of the Saskatchewan Prairies are the powerhouses of water supply in the province. This is especially true across the heavily farmed semi-arid grasslands of central and south Saskatchewan. This fertile region has long relied on irrigation to counteract hot, dry summers and low precipitation rates [1,2]. The construction of the Gardiner and Qu’Appelle River dams, and the creation of Lake Diefenbaker in 1967, added nearly 220,000 irrigated acres to the province alone [2]. Gardiner Dam provides water for Coteau Creek Hydroelectric Station and passes water north along the South Saskatchewan River towards Saskatoon, the largest city in the province. The Qu’Appelle River Dam passes water southeast along the Qu’Appelle River towards the drinking water reservoir Buffalo Pound Lake (BPL). An impounded natural lake, BPL is relatively small compared to many Canadian water bodies yet is more important than the size would suggest, supplying approximately 25% of the Saskatchewan population with their water requirements [3].

BPL has been the focus of several studies over the years due to problems with water quality. The contributing watershed is extensively farmed [4,5], and the reservoir is plagued by algal blooms and high concentrations of nutrients and organic carbon [6,7,8]. The annual operating cost is high for the on-site Buffalo Pound Water Treatment Plant (BPWTP). Like many waterbodies on the prairies (e.g., [9,10,11,12]), BPL has multiple stressors that are creating challenges for the supply of safe drinking water. Water demand is expected to increase in the future, and questions have been raised about how the reservoir will respond to changes in water management and climate change.

In recent years, water quality models have been an emerging tool in the prairie river basins. Most recently, a two-dimensional (2D) CE-QUAL-W2 (W2) water quality model of BPL ran scenarios of increased water transfers from the supply reservoir Lake Diefenbaker [7]. The authors of this study describe how several wet years of large amounts of watershed run-off deteriorated the water quality in BPL. Lake Diefenbaker, in contrast, has a different watershed area and is fed by water originating in the Rocky Mountains. The water of Lake Diefenbaker is of better quality than that of the downstream BPL. The flat landscape of the Qu’Appelle system means that any increase in water transfers must be carefully planned to avoid flooding the surrounding area. Management agencies were interested in modeling how quickly BPL might improve if inflows were augmented. The modeling study indicated BPL would respond rapidly to water transfers and that the flow rate and seasonal timing were crucial factors in the reservoir response. However, there was a question remaining about the water quality of the simulated water transfers. When adding the boundary conditions to the BPL model, it was uncertain how to estimate any nutrient transformations occurring between Lake Diefenbaker and the inflow to the reservoir. The study attempted a transfer function based on a linear regression between observations at an upstream and downstream location between the two reservoirs (a stretch known as the Upper Qu’Appelle River, or UQR). The observation dataset for the regression was infrequent, and the transfer function was limited in its application. The study concluded that an improved method of estimating the nutrient transformations would be required before proceeding with further work.

Models are effective tools to understand and predict future surface water quality. When data are unavailable or infrequent, water quality models can examine scenarios and be used for assessment when conditions are not directly measurable. An earlier research project by Hosseini et al. [13] applied a one-dimensional (1D) model of the UQR-BPL system in WASP. The project assessed linking a river and reservoir in a single model structure, and a 1D model was deemed sufficient for this purpose. The UQR section of the model returned good calibration results for the modeled variables. The BPL section did not calibrate so well, with mean absolute errors approximately doubled for the nutrients. The authors concluded a 2D modeling approach would have been ideal for the river–reservoir chain.

This study advanced the BPL water management scenarios by coupling the 2D W2 reservoir model to the updated 1D WASP model of the UQR. The objective was to improve the estimation of nutrient transformations in the water transferred along the UQR and into BPL. As WASP and W2 have different structures, this study tested a loose coupling setup where water transfer scenarios were first simulated in the WASP model. The resulting water quality and flow outputs from the river section of WASP were then passed manually as inputs to the W2 reservoir model. The impacts of the water transfers on BPL were then determined by comparing the outputs of the W2 reservoir model against predetermined threshold limits for water quality.

2. Methods

2.1. Study Site

The Qu’Appelle River flows 430 km east from Lake Diefenbaker to the Assiniboine River in Manitoba. The upper stretch (UQR) flows from Lake Diefenbaker through BPL to the town of Craven, Saskatchewan, where it enters the Lower Qu’Appelle River Watershed.

The UQR from Lake Diefenbaker to BPL is a 97 km channel that has both improved and natural meandering sections. The first 35 km are channelized with a design capacity of 14 m³/s [14]. The channel has experienced loss of flow capacity in the past due to erosion, siltation, and macrophyte growth. Mean daily discharges have ranged from 1.4 m³/s to 5.9 m³/s over a 20 year period (min 0.013 m³/s/max 14.2 m³/s; 1999–2019, Water Security Agency hydrometric data) at flow gauge 05JG006 (Elbow Diversion) situated 3.3 km below Lake Diefenbaker (Figure 1). Winter flow is reduced even further to ensure safe ice cover formation [14]. In recent years, the Water Security Agency (WSA) has worked to improve the channel flow capacity by deepening sections and through erosion control. The current channel capacity is considered sufficient to meet anticipated water demand [15].

Gauged and ungauged tributaries and ephemeral streams feed into the UQR along the 97 km stretch. The largest of these are the two tributaries, Iskwao Creek and Ridge Creek (Figure 1). The part of the Upper Qu’Appelle watershed that drains into the UQR is shown in Figure 1.

BPL was a natural shallow lake formation on the UQR that was impounded during the construction of the Buffalo Pound Dam in 1939. The reservoir was created to supply the water needs of the provincial capital Regina, the city of Moose Jaw, surrounding communities, industry, and potash mines. The Buffalo Pound Water Treatment Plant (BPWTP) is located at the southern end of the reservoir, close to the dam. BPL has a total length of 30 km, an average width of 890 m, and a mean and maximum depth of 3.8 m and 6.0 m, respectively. A bathymetry map is available in Terry et al. [7]. The reservoir is ice-covered, on average, from November to April and is wind-mixed readily when ice-free.

The reservoir is essentially split into two sections by a provincial highway (Hwy #2), as shown in Figure 1. The highway dams the reservoir down to the bed save for a 45 m overpass that the water squeezes under (three × 15 m sections between in-channel piles). The old highway now acts as a wave break to protect the new highway and has a 53 m gap that helps slow inflowing water before it reaches the overpass. The top section of the reservoir above these highways is comparatively shallow.

The prairie region endures large temperature changes over the course of a year, with temperatures at Moose Jaw A climate station ranging from an average daily maximum of 26.2 °C to an average daily minimum of −17.7 °C [16]. Approximately 30% of annual mean precipitation falls as snowfall. The region is semi-arid, although both floods and droughts are frequent occurrences [9,17]. Uncontrolled run-off is common in very wet years; however, the precipitation-to run-off ratio is complicated by antecedent soil moisture and the fullness of wetlands within the watershed [7]. In dry years, most of the inflows into BPL are controlled releases from Lake Diefenbaker other than at spring freshet.

Discharge within the UQR can be rapid, with a calculated travel time of under two days in high flows [7]. BPL residence time is extremely variable at 6 to 36 months [18] due to the variability of watershed run-off and the need to carefully control river flows. On infrequent occasions, flood waters of poorer water quality from Moose Jaw Creek downstream can backflow over the gates of the dam and into the reservoir.

2.2. Model Description

2.2.1. WASP

The WASP water quality model was first released into the public domain in 1981 and has been in continuous development since then by the U.S. Environmental Protection Agency [19]. The model was most recently restructured with the release of WASP8 in 2016. Changes are listed in Wool et al. [19] and include a new programming environment, changes to the coding and software, a new GUI, and a new kinetic structure. The original UQR model was built and calibrated in WASP7 [13], which is a well-tested and stable version of the WASP model framework [20,21]. This study continued with the WASP7 version for consistency.

A full description of the structure of the WASP model, as setup for 1D modeling of the UQR and BPL, is available in Hosseini et al. [13]. To summarize here, the original model grid consisted of 187 longitudinal segments. The UQR section contained 165 longitudinal segments of approximately 600 to 800 m. The river segmentation was originally determined from a HEC-RAS model provided by the WSA, developed from 770 survey locations along the 97 km river stretch. Additional inflows from Iskwao Creek, Ridge Creek, and seven ungauged tributaries were added to the main river flow network. The reservoir section of the model in WASP consisted of 22 longitudinal segments of around 1400 m. The cross-sectional profiles were determined from a digital elevation model (DEM) created in ArcGIS 10.2.2 for an earlier W2 modeling project of BPL (see [22]). One-dimensional kinematic wave routing was selected to model flows through the river–reservoir system. This study used the same hydrological setup as the previous authors but incorporated updated boundary data.

2.2.2. CE-QUAL-W2

Built on the simple reservoir model (LARM), W2 was first released in the complex CE-QUAL-W2 (W2) format (version 1.0) in 1986 by the U.S. Army Corps of Engineers [23]. Since then, the model has been applied to numerous waterbodies worldwide over the decades [21,24] and is in continuous development by Portland State University. The 2D model is laterally averaged and suitable for a long, narrow waterbody such as BPL. The latest version (V4.5) was released in 2021. The release included various new features and upgrades and was not file-compatible with earlier versions [23]. This study used CE-QUAL-W2 V4.2.2, which was the latest well-tested version available at the time of developing the existing water quality model of BPL.

A full description of the structure of the W2 model, as setup for 2D modeling of BPL, is available in Terry et al. [7]. Condensed here, the reservoir model consisted of 100 longitudinal segments of around 300 m and up to 28 vertical layers of 0.25 m depth. Constrictions were added to the model grid to reflect the narrow water passages at the location of the highways (Figure 1). A schematic of the BPL model grid and bathymetry is available in Terry et al. [7]. Inflows included the main UQR channel, a distributed tributary to account for ungauged flows from the watershed, and a tributary into the most downstream segment to simulate any backflows from Moose Jaw Creek. An outflow dam structure, one intake pipe for withdrawals to the water treatment plant, and one intake pipe supplying the surrounding industries were all included in the model.

2.3. Data and Calibration

2.3.1. WASP (River)

A full description of the calibration and testing of the WASP model for the period from January 2012 to December 2015 is described in Hosseini et al. [13]. For this current study, the first step was to update the full WASP model with newer data up to December 2019. This brought the UQR model in line with the simulation period of the W2 reservoir model. The original simulation period of the WASP model coincided with a three-year sampling initiative in the Upper Qu’Appelle watershed by the WSA that ended in March 2016 [25]. During this earlier period, additional water quality monitoring data had been available for both Iskwao Creek and Ridge Creek, the two main tributaries of UQR. After the WSA stopped monitoring the sites, the only water quality stations available were located in the main UQR channel. As this meant fewer data were available for the new time period, the additional four years were treated as a validation period for the existing model parameter set.

The water quality stations used in this model are shown in Figure 1. Data were provided by both the WSA and the BPWTP for this study. Inflowing concentrations to the UQR stretch were taken from the WSA station Highway 19 (Hwy #19). The data frequency at Hwy #19 averaged from daily to monthly in summer and intermittent in winter due to the logistics of taking water samples under ice cover. As no water quality stations were available for the UQR tributaries after March 2016, water quality concentrations had to be estimated for the new time period. Monthly concentrations between March 2016—December 2019 were assumed to be the same as the monthly mean values over the years of the WSA sampling initiative, with the same values repeated each year. Concentrations for the seven ungauged tributaries added to the WASP flow network were assumed to be the same as Ridge Creek. Observations from WSA stations Tugaske and Marquis were used for checking model performance in the river section of the model.

Observations from the WSA station Hwy #2 and the BPWTP for their intake location were used for checking model performance in the reservoir. Of the datasets, the BPWTP intake observations had the most winter data. Weekly sampling was performed at the pumping house before water was passed for treatment. Water withdrawals were through an intake pipe situated one meter off the bottom of the reservoir bed. Water traveled 300 m before reaching the pumping house, and concentrations underwent small changes between the sample location and reservoir due to transportation.

Daily inflow data for the main UQR were from the WSA gauge Elbow Diversion (05JG006). For the gauged tributaries, daily flows were available for WSA gauge Ridge Creek (05JG013) for March–October annually. Flows are not recorded by the WSA between November–February, and the tributaries are usually frozen over in the winter months. The seven ungauged tributaries were assumed to have the same flows as Ridge Creek. The Iskwao Creek gauge (05JG014) was discontinued in 2011, but estimated daily flows were provided by the WSA for the purpose of this study. Multiple time series can be added to WASP for water temperature, with the user able to select which time series is allocated to each segment of the model grid. Three water temperature time series were included in the WASP model based on observations at stations Tugaske (upstream river segments), Marquis (downstream river segments), and the BPWTP (reservoir segments).

For the meteorological data, air temperature and wind speed were downloaded from the Environmental Canada and Climate Change (ECCC) online database for stations Moose Jaw A and Moose Jaw CS. ECCC were in the process of phasing out the Moose Jaw A station in preference of Moose Jaw CS over the time period of the study. Solar radiation was downloaded from NOAA gridded climate datasets (NCEP-NCAR Reanalysis 1 data, NOAA PSL, Boulder, CO, USA, https://psl.noaa.gov (accessed on 28 December 2021)). The fraction of daylight was estimated from online sunrise and sunset times [26].

2.3.2. CE-QUAL-W2 (Reservoir)

A full description of the calibration and testing of the W2 reservoir model for the period April 2013–December 2019 is described in Terry et al. [7]. Water quality data was provided by the WSA for the Marquis station, 13 in-reservoir stations, and one station at the dam outflow, as marked in Figure 1. The BPWTP again provided their weekly water sampling data for the intake location.

Inflows were taken from daily estimates provided by the WSA for the period April 2013–May 2015, and then Keeler Bridge (05JG004), upstream of BPL, after the gauge was reinstated in June 2015 (had ceased operation in 1995). Ungauged inflow estimates, backflow estimates, dam outflows, and in-reservoir water level observations were also provided by the WSA. Water temperatures are allocated to inflow sources in the W2 model. Temperatures from the Marquis station were used for the main inflows and tributaries and from measurements just below the Buffalo Pound Dam for backflow estimates.

Hourly meteorological data (air and dewpoint temperature, wind speed and direction, and cloud cover) were downloaded for the same Moose Jaw A and Moose Jaw CS climate stations used for the WASP model. Precipitation data were downloaded for the ECCC station at Buffalo Pound Lake, located at the site of the BPWTP. Precipitation temperatures were set at dewpoint temperature (or zero if dewpoint was negative) per W2 recommendations, and precipitation water quality constituent concentrations were not modeled.

2.3.3. WASP and CE-QUAL-W2 (River–Reservoir Base Model)

The WASP and W2 models remained separate structures but were applied together in a loose coupling, where output water quality data from WASP became input boundary data for W2. After the WASP model was validated with the new data, the reservoir segments were no longer required for loose coupling, as only the UQR section was needed for modeling the nutrient transformations in the water transfers. The WASP model ran as normal with simulation results output at segment 165, the most downstream river segment in the model grid. This output was then manually converted to the W2 input file format.

The data requirements for the two models were mostly compatible, although some estimations were required. Algae were modeled as one group in WASP and three groups in W2. Algal biomass per WASP was split as 45% diatoms, 45% green algae, and 10% cyanobacteria for W2. Total organic carbon from WASP was used to estimate the dissolved and particulate organic matter compartments required by W2 by assuming a constant organic carbon content of 45% in organic matter. The organic nitrogen (N) and phosphorus (P) content was directly available from WASP and split as 20% labile and 80% refractory for the organic subcompartments in W2. Note that when these subcompartments are switched on, W2 no longer converts the main organic pool into N and P to avoid double-counting.

One difference between the two model setups was the approach to water balance. WASP is primarily a nutrient transport and cycling model that can be linked to hydrodynamic and sediment transport models for complex flows, velocities, and fluxes. For the WASP setup, flow data were entered at the upstream boundary and at model segments where flows increased, for example, from the joining of a tributary. WASP calculated the cumulative flow appropriate to the 1D kinematic wave routing. The outflows essentially are the cumulative flow leaving the final segment in the model grid. Withdrawals for the BPWTP and precipitation were not incorporated into the WASP setup, and Hosseini et al. [13] did not investigate the water balance of the UQR-BPL system as part of their study. Flow estimates along the UQR are subject to a great degree of uncertainty due to the complex hydrology of the watershed. Prior to the reinstallation of the Keeler Bridge flow gauge (05JG004) in June 2015, the WSA’s regression models only explained a percentage of the variance of ungauged flow estimates provided by the WSA (see [25]). During the wet year of 2014, the flow estimates were further complicated by flooding in the Qu’Appelle Valley.

The W2 model is a purpose-built coupled hydrological–ecological model. In the general setup the user provides all the inflow and outflow data (including tributaries, withdrawals, structures, etc.) in advance of the simulation. Having the model grid discretized into layers as well as segments means the user must ensure the water balance agrees before they can progress to the water quality stage. When too many layers lose water at once, the simulation is terminated with an error code. The water balance can often be an iterative process where inputs and outputs are added to the model until the model runs from beginning to end. The model developers recommend using the included water balance tool that comes with the W2 package for this purpose. The tool compares simulated water levels against observations and generates a file of adjustment flows required to balance the water volume that the user should investigate. The water balance procedure for the 2013–2019 W2 reservoir model is discussed in Terry et al. [7].

For the purposes of the scenario work, it was elected to output WASP flows at segment 165 and incorporate them directly as inflows into the W2 model. For the base scenario, the W2 tool was then used to create the additional inflows and outflows required for the W2 water balance. The adjustment flows for the base model simulation are presented in Figure 2. The flows were added as a distributed tributary around the perimeter of the reservoir. Water temperatures were assumed to be the same as the main inflow file, and no constituents were modeled. This method allowed the new WASP flow and concentration data to be output together for each scenario.

Some of the calibration coefficients used in both models are included in

Table 1. A summary of coefficients used in both the WASP river section and CE-QUAL-W2 (W2) reservoir models.

Description	WASP River (Section)	W2 Reservoir
Sediment oxygen demand (SOD)	0.41 (g O2 m2/day)	0.1–1.2 (g O2 m2/day)
Sediment release rate of phosphorus	0 (mg/m2/day)	0.015 (fraction of SOD)
Sediment release rate of ammonium	0 (mg/m2/day)	0.2 (fraction of SOD)
Nitrification rate	0.01 (1/day)	0.12 (1/day)
Denitrification rate	0.09 (1/day)	0.1 (1/day)
Phytoplankton	1 group	3 groups
Maximum growth rate	3 (1/day)	0.5–2 (1/day)
Death rate	0.006 (1/day)	0.08–0.1 (1/day)
Half-saturation constant for phosphorus	0.0048 (mg P/L)	0.003 (mg P/L)
Half-saturation constant for nitrogen	0.26 (mg N/L)	0–0.014 (mg N/L)
Ratio of algal biomass to carbon Ratio of algal biomass to chlorophyll-a		0.45 (fraction) 0.04–0.1 (mg algae/µg chl-a)
Carbon to chlorophyll-a ratio	20 (mg C/mg chl-a)
Fraction of phytoplankton recycled at death	0.95	0.8 (all groups)

for comparative interest.

2.4. Water Management Scenarios

Five water management scenarios were selected to simulate the release of additional water from Lake Diefenbaker over the 2013–2019 period (

Table 2. Description of the five scenarios for augmented water transfers from Lake Diefenbaker to Buffalo Pound Lake.

Scenario	Description
SC1	In WASP: Add 4 m3/s daily to the main UQR inflows for the full simulation period. In W2: Add 4 m3/s daily to outflows to help with the water balance.
SC2	As above, with 8 m3/s.
SC3	As above, with 12 m3/s.
SC4	The maximum design capacity in the improved upper section of the UQR channel is 14 m3/s in open water periods and 6 m3/s in winter [27]. Assuming capacity is the same in the lower channel, and that open water is May–October and winter is November–April: In WASP: Round-up daily flows at both Elbow Diversion (upstream) and Keeler Bridge (downstream) to the nearest integer and determine the larger of the two values; Subtract that value from the max capacity (14 or 6 m3/s); If the answer is positive, then there is spare capacity in the channel, so add the additional water transfers to the base model inflows; If the number is negative, then flow in the channel already exceeds capacity, so leave the inflows as per the base model. In W2: Add the equivalent flows to outflows to help with the water balance.
SC5	As above, but where water transfers are simulated in open water periods only. No water transfers are added between the winter months of November–April.

). The first three scenarios augmented flows by a selected fixed amount year-round. These flows were of theoretical interest to the WSA (Davies, J.-M., WSA, personal communication) and used to determine whether the model was behaving as expected. The fourth and fifth scenarios were chosen based on the maximum design capacity of the improved upper section of the UQR channel. The design capacity is seasonal, 14 m³/s in open water, but with a suggested maximum of 6 m³/s in ice-covered conditions [27]. These last two scenarios were more realistic in that they considered the current capacity of the channel and infrastructure in the UQR.

As per the base model, the model scenarios ran in two stages to work with the loose coupling setup. The first stage ran the five scenarios in the WASP model. For each scenario, the new outflows and water quality concentrations were then output at segment 165 and converted into boundary data for the W2 model. The second stage then ran the W2 model for each of the five scenarios with the relevant new input files. The results of the W2 model simulations were then compared against predetermined concentration thresholds to measure the impact of each management option.

Concentration Thresholds

Both TDS and DOC are problematic for drinking water quality and are variables of interest for the BPWTP. Health Canada has set an objective aesthetic guideline of ≤500 mg/L for TDS in drinking water, stating concentrations above this result in excessive scaling and unpalatability [28].

Organic matter can also result in unpleasant taste and odor and affect the efficiency of drinking water treatment processes and disinfection [29]. Although no current guideline is in place, the Saskatchewan authorities have previously suggested DOC concentrations greater than 5 mg/L will complicate drinking water treatment [30].

Total phosphorus (TP), total nitrogen (TN), and chlorophyll-a (CHLA) are indicators of interest for the health of aquatic ecosystems. The Canadian Council of Ministers of the Environment (CCME) recommends using TP as a meaningful measurement of phosphorus in water [31]. No federal quality guideline exists for TP as it is not a toxic substance to aquatic organisms [31]. The CCME recommends trigger ranges that define a desirable concentration range for individual waterbodies based on their reference conditions. When the upper limit of a range is exceeded, it ‘triggers’ the need for further assessment. TP concentrations at the intake location averaged 0.08 mg/L as P (BPWTP data 1980–2019, 1216 observations). This fits BPL into the eutrophic category with a trigger range of 0.035–0.1 [31]. The CCME suggested the second criterion would be when concentrations fell within the trigger range yet were 50% above baseline. The baseline concentration for TP was 0.08 mg/L at the BPWTP intake location. The 50% above the baseline value of 0.12 mg/L was already higher than the trigger range, so this test was not required.

No Canadian federal or provincial water quality guideline or standard exists for total nitrogen (TN), and the United States Environmental Protection Agency is still working with individual states to approve criteria for TN concentrations. Using similar criteria to TP, the same 50% above baseline method was chosen to determine a concentration threshold for TN. Concentrations at the intake location averaged 0.91 mg/L with a 50% above the baseline threshold value of 1.36 mg/L.

Similarly, no standardized water quality guidelines exist for chlorophyll-a concentrations, and the 50% above baseline was chosen as a concentration threshold. Concentrations at the intake location averaged 20.8 µg/L with a 50% above the baseline threshold value of 31.2 µg/L.

3. Results and Discussion

3.1. Base Model Calibration

Figure 3 plots the WASP results in their entirety, with both the river and reservoir sections of the 1D model. The red lines indicate the date at which the observed water quality data ended for the two main tributaries to the UQR (Iskwao Creek and Ridge Creek). Model performance statistics are provided in

Table 3. Mean absolute error and R² adjusted for the WASP river–reservoir model and the CE-QUAL-W2 reservoir model for 2013–2019. Tugaske and Marquis are river water quality stations. Hwy2 and Intake are reservoir water quality stations. Concentrations mg/L (µg/L CHLA). #obs = number of observations.

WASP		TDS	DOC	DO	TP	TN	CHLA
Tugaske	MAE	139.91	1.85	1.39	0.03	0.20	-
	Adjusted R2	0.68	0.41	0.36	0.81	0.75	-
	#obs	133	133	129	133	133	-
Marquis	MAE	195.16	1.75	1.17	0.12	0.27	-
	Adjusted R2	0.50	0.81	0.54	0.22	0.69	-
	#obs	143	142	147	143	143	-
Hwy#2	MAE	127.53	3.09	1.85	0.07	0.31	13.25
	Adjusted R2	0.52	0.55	0.21	0.26	0.30	0.02
	#obs	19	19	26	26	26	26
Intake	MAE	174.55	2.40	1.45	0.06	0.31	20.43
	Adjusted R2	0.01	0.67	0.31	−0.01	0.06	0.04
	#obs	119	222	157	121	116	142
CE-QUAL-W2		TDS	DOC	DO	TP	TN	CHLA
Hwy#2	MAE	133.32	3.43	1.61	0.06	0.23	14.88
	Adjusted R2	0.48	0.34	0.49	0.51	0.51	−0.02
	#obs	19	19	26	26	26	26
Intake	MAE	140.94	1.39	1.39	0.04	0.23	17.22
	Adjusted R2	0.38	0.58	0.22	0.04	0.18	0.31
	#obs	119	222	157	121	116	142

. As was found by the original Hosseini et al. [13] study, the WASP model performed better in the river section than in the reservoir. This was also the case for both TDS and DOC, which were new variables included in this study. This is largely expected for the UQR-BPL system. UQR responds quickly to watershed conditions, with water quality concentrations largely connected to loading from Lake Diefenbaker and the watershed. BPL, in contrast, has a residence time of up to three years, and the key environmental drivers of change are not known [8]. BPL has not shown typical patterns of reservoir ontogeny over the years (Baulch, H. M., GIWS, personal communication). An earlier W2 model of BPL performed poorly for the period 1986–1993, and the challenges of modeling shallow lakes and modeling BPL were thoroughly discussed by the authors of the study (see [32]).

Per the plots, it is likely the loading for TDS concentrations was overestimated during the freshet of 2018. For this period, the tributary loading concentrations were estimated based on earlier years, with no further adjustments made to the input values.

Figure 4 plots the results for the same locations for the loosely coupled base model, where the reservoir section now shows the output of the W2 model. Model performance statistics for the W2 reservoir model are presented in Table 3. The observed water quality data are the same in both Figure 3 and Figure 4. While the W2 model performed similarly to the WASP model at the highway location, the number of observations for that station is small. The more indicative station is the BPWTP intake, which had a greater number of observed data available for calibration purposes. Here, the 2D W2 model outperforms the 1D WASP model for all variables based on mean absolute error. The location of the intake is the area of interest from a management perspective, as this is where the drinking water withdrawals occur.

The W2 results also indicate that TDS concentrations may have been overestimated in the inflow boundary data during the 2018 freshet.

3.2. Scenarios

Table 4. Number of threshold exceedances at the treatment plant intake location over the simulation period 2013–2019.

Variable	Threshold	Base	SC1	SC2	SC3	SC4	SC5
TDS	>500 mg/L	1838	1001	811	671	829	875
DOC	>5 mg/L	1877	1116	739	566	627	637
TP	>0.1 mg/L as P	1176	526	378	325	346	360
TN	>50% above baseline (1.36 mg/L)	51	8	33	32	31	31
CHLA	>50% above baseline (31.2 µg/L)	878	773	739	700	726	724

presents the number of times that concentrations exceeded threshold values per scenario at the location of the BPWTP intake pipe.

Scenarios 1–3 were the test simulations. Lake Diefenbaker had better water quality than BPL, and the expected response of the model was a decrease in threshold exceedances as the simulated volume of additional water increased. This occurred for all variables except TN, where the lowest amount of threshold exceedances was SC1 at 4 m³/s. Figure 5 plots the results at the BPWTP intake location for TP and TN. By contrasting the two nutrients, it can be seen that the threshold selection method is more insightful for TP (trigger range) than TN (50% above baseline). In addition, as can be seen for scenarios 1–3, while the overall nutrient concentrations decrease with the increases in flow volume, the initial pulse of nutrients arriving at the intake has the opposite behavior (i.e., the maximum pulse increases with the flow). Sediments will be mobilized quickly as flow rates increase, and the inflowing sediment plume will be carried further and faster into BPL, arriving at the intake earlier and in a more concentrated pulse with each step-up in flow rate. This and the elevated threshold level explain why TN showed better results at only 4 m³/s, as only one nutrient pulse breached the threshold in the scenarios.

Scenarios 4 and 5 were the more practical options in that they assumed that water transfers would be maximized flows but kept within the design carrying capacity of the existing UQR channel. The relative closeness of the results was interesting, considering scenario 5 only simulated additional releases from Lake Diefenbaker in the open water periods. Augmenting the Lake Diefenbaker flows in winter had little impact on the threshold exceedances. Most loading of poor water quality occurs during freshet and after rainfall in the open water season. It would coincide that there are more nutrient inputs in the open water seasons that can be offset by the addition of cleaner water. The previous modeling attempt of Terry et al. [7] also suggested the optimum timing would be to commence water transfers once the spring freshet had subsided and terminate transfers before the subsequent spring melt. These current results indicate transfers could be terminated even earlier before the onset of winter.

This finding is valuable from a management perspective. Augmenting UQR flows in winter increases the risk of ice scouring and increased erosion due to the shallow water depth [33]. Flows should generally be ramped down over winter as ice thickens [33]. From the modeling results, we suggest that increasing Lake Diefenbaker transfers in open water periods only would see sufficient improvement in BPL without risking damage to the armoring and bed liner in the improved upper section in winter.

In the previous attempt to model increased water transfers into BPL (see [7]), changes in nutrient concentrations in the 97 km stretch between Lake Diefenbaker and BPL were estimated for a given flow rate using linear regression and observed data. Several challenges were discussed with regard to the method in Terry et al. [7], and the amount of variance explained by the regression models was low (Adj. R² at Marquis station: TDS (0.123), DOC (0.133), TP (0.540), TN (0.032), and DO was not modeled). Of these, TP concentrations were better explained by the regression model based on 40 observations; however, one concern was that the linear regression method returned a positive relationship between TP and flow. The result was TP concentrations increased in the water transfers in line with the flow rate. The resulting reservoir TP concentrations increased with the Lake Diefenbaker transfers, which was against expectations (Lake Diefenbaker concentrations being lower than those of BPL), and in contrast to the other modeled variables. Soluble P is strongly sorbed to suspended solids and sediments [34,35], so concentrations would likely increase with sediment mobilization along the UQR as flow increased. In a real aquatic environment, however, the stores of P would eventually be depleted and increasing water volume would dilute concentrations.

With the aim of improving this method, an older 1D WASP model of the UQR was updated to the same time period to provide boundary flows and nutrient concentrations at the inflow to BPL. Marquis station is close to the BPL inflow and to segment 165, the output of the WASP model. Per Table 3, the WASP model showed an improved fit of data at Marquis station to the linear regression model for most variables (Adj. R² at Marquis station: TDS (0.5), DOC (0.81), TP (0.22), TN (0.69), and DO (0.54)). TP was the exception, although WASP showed a negative relationship between TP and flow (see Table 4), as was expected. Although the adjusted R-squares at Marquis are not particularly strong, the UQR is a difficult system to model due to the complex hydrology of the watershed and the need to estimate data for the tributaries joining the stretch. This is especially true after March 2016, when the WSA stopped sampling the two main tributaries. Compared to the regression model based on available observations, the WASP model provided more information for estimating the nutrient transformations, such as the changing relationships over time.

4. Conclusions

The modeled results show a one-dimensional WASP river model improved the estimation of nutrient transformations over a linear regression model based on observations of water quality and flow. The WASP modelled flows and concentrations then created more realistic boundary conditions for the downstream two-dimensional CE-QUAL-W2 reservoir model. Running the WASP and CE-QUAL-W2 models together in a loose coupling for water management scenarios, results suggest maximizing water transfers in the open water season would benefit Buffalo Pound Lake in terms of water quality improvement. Most nutrients arrive at Buffalo Pound via the watershed during the spring freshet. The interbasin transfers of better-quality water from Lake Diefenbaker show more impact on threshold exceedances during the subsequent open-water season. In contrast, continuing the transfers at reduced flow capacity over winter showed only a small reduction in the number of threshold exceedances. This is useful knowledge for winter operations as the risk of ice damage to the channel grows with augmented winter flows.