Representation of Water Management in Hydrological and Land Surface Models

Reservoirs significantly affect flow regimes in watershed systems by changing the magnitude and timing of streamflows. Failure to represent these effects limits the performance of hydrological and land surface models (H-LSMs) in the many highly regulated basins across the globe and limits the applicability of such models to investigate the futures of watershed systems through scenario analysis (e.g., scenarios of climate, land use, or reservoir regulation changes). An 15 adequate representation of reservoirs and their operation in an H-LSM is therefore essential for a realistic representation of the downstream flow regime. In this paper, we present a general parametric reservoir operation model based on piecewise linear relationships between reservoir storage, inflow, and release, to approximate actual reservoir operations. For the identification of the model parameters, we propose two strategies: (a) a “generalized” parameterization that requires a relatively limited amount of data; and (b) direct calibration via multi-objective optimization when more data on historical 20 storage and release are available. We use data from 37 reservoir case studies located in several regions across the globe for developing and testing the model. We further build this reservoir operation model into the MESH modelling system, which is a large-scale H-LSM. Our results across the case studies show that the proposed reservoir model with both of the parameter identification strategies leads to improved simulation accuracy compared with the other widely used approaches for reservoir operation simulation. We further show the significance of enabling MESH with this reservoir model and discuss 25 the interdependent effects of the simulation accuracy of natural processes and that of reservoir operation on the overall model performance. The reservoir operation model is generic and can be integrated into any H-LSM. Hydrol. Earth Syst. Sci. Discuss., https://doi.org/10.5194/hess-2019-7 Manuscript under review for journal Hydrol. Earth Syst. Sci. Discussion started: 16 January 2019 c © Author(s) 2019. CC BY 4.0 License.


Background and Motivation
Human interventions in natural hydrologic systems, through damming and storing water, diversion, surface and groundwater abstraction, irrigation, and land use change, have significantly altered the natural river flow regimes and the terrestrial water cycle of many river basins (Vörösmarty et al., 1997(Vörösmarty et al., , 2003;;Oki and Kanae, 2006;Wisser et al., 2010;Haddeland et al., 2014).These interventions are to fulfil different types of demands such as domestic, industrial, irrigation, and hydropower demands, and to meet other needs such as flood control and conservation of aquatic habitats.With a total storage volume of more than 8000 km3 (ICOLD, 2003;Hanasaki et al., 2006;Vörösmarty et al., 2003), more than 50,000 dams have been constructed globally to regulate more than half of the world's large river systems (Nilsson et al., 2005).The aggregate storage volume of these dams is greater than 20% of the global mean annual runoff (Vörösmarty et al., 1997) and is three times the annual average water storage in world's river channels (Hanasaki et al., 2006).
Despite the benefits in terms of enhancing water availability in support of food security, power supply, etc., dams result in several negative environmental and social consequences.Adverse environmental effects include changes in natural river dynamics in terms of water temperature, sediment and nutrient transport, etc. and the fragmentation and loss of biodiversity (Vörösmarty et al., 2010).Reservoirs can also intensify evaporation, by increasing the surface area of water exposed to direct sunlight and air, and through water supply for irrigation (de Rosnay et al., 2003;Pokhrel et al., 2012).Other environmental impacts of dams include the alteration of landscape due to dam construction and changes to land-atmosphere interaction that can have a profound impact on local/regional climate (Hossain et al., 2012).Adverse social effects include the displacement of people living near the dam site, changes to fishing patterns, and downstream erosion (Strobl and Strobl, 2011, p. 449).
There are research gaps remaining in evaluating both positive and negative social impacts of dams (Kirchherr et al., 2016).
Such gaps have been the subject of many studies in both academia and industry for years, and recently, have led to the formalization of this study area of "socio-hydrology" (Sivapalan et al., 2012;Sivakumar, 2012).
Dams and reservoirs change the natural flow regimes in rivers, both in terms of magnitude and timing of flows.As a result, for rivers that contain large or small dams and reservoirs, flow regimes are a combination of natural and managed flows.
Various modeling communities manage this mix of natural and managed flows differently.Archfield et al. (2015) compare three communities of models that can be used at continental scales: catchment models (CM), global water security models (GWSM), and land-surface models (LSM).CMs generally ignore water management and focus on unmanaged headwater catchments.GWSMs have been utilized in global-scale streamflow simulations and generally consider large-scale water management, which is made difficult by a lack of data on large-scale water management operational decisions.LSMs have traditionally focused on providing lower boundary conditions for atmospheric models, but are increasingly being used for hydrological applications in which they are referred to as Hydrologic Land Surface Models (H-LSMs).LSMs generally ignore water management (Clark et al., 2015;Davison et al., 2016), with a few exceptions (e.g.Voisin et al. 2013aVoisin et al. , 2013b)).
A fourth community of water models, that is relevant to the work presented here, are water management models (WMM) Hydrol.Earth Syst.Sci.Discuss., https://doi.org/10.5194/hess-2019-7Manuscript under review for journal Hydrol.Earth Syst.Sci. Discussion started: 16 January 2019 c Author(s) 2019.CC BY 4.0 License.(Yates et al., 2005;Labadie, 1995).Water modellers who know how the water is managed within their basins of interest generally use WMMs (Lund and Guzman, 1999;Labadie, 2004;Kasprzyk et al., 2013).These models contain very detailed representations of water management decisions, but often consider natural flow processes in a much more rudimentary fashion than CMs.
Modeling the many managed basins around the world using CMs or LSMs can result in models with limited fidelity and question the credibility of their predictions of future water resources in basins with dams and reservoirs.Therefore, there is a pressing need for better characterization and integration of the operation of dams and reservoirs into hydrological modelling frameworks using CMs and LSMs (Nazemi andWheater, 2015a and2015b;Wada et al., 2017;Pokhrel et al., 2016).This need motivated the objectives of this study, described in Section 1.2, and some previous research, outlined in Section 2. The integration of reservoir regulation into hydrological modeling frameworks will improve our ability to simulate highly regulated basins around the globe, leading to better understanding of historical conditions of water resource systems and improved assessment and prediction of their future vulnerability to climate and environmental change.

Objectives
Building upon previous research, this study aims to:  Develop and test an improved reservoir operation model that can be integrated into any CM and LSM at any scale, but in particular at large scales.Of interest is a simple but effective parametrization that can adjust to varying levels of data availability.
 Integrate the developed reservoir operation model into an LSM and evaluate its performance when working in combination of other processes in the model.Also of interest is to assess the potential conceptual and technical issues in this integration.This paper looks to improve the representation of dams and reservoirs within CMs and LSMs.Another approach would be to couple CMs and LSMs with WMMs, but that approach is not examined here due to the fact that WMMs generally require information on how water is managed within a basin, and we are particularly interested in the more generic case when this information is likely to be limited or unavailable.
The organization of the remainder of the paper is as follows.Section 2 reviews different existing approaches in the literature for the representation of reservoir operation in hydrologic models.Section 3 presents the proposed reservoir operation model and the metrics used to evaluate it, in comparison with other existing models.Section 4 provides a description of the reservoir dataset used for the developments and testing.Section 5 presents the assessment results and comparisons.Section 6 ends the paper with a summary of the main findings and conclusions.

Existing Reservoir Models in Catchment Models and Land Surface Models
An adequate representation of human interventions in Earth systems models is a major challenge.Systematic approaches towards full integration are needed as outlined in the recent studies of Nazemi andWheater (2015a and2015b), Wada et al. (2017), andPokhrel et al. (2016).In this work, our focus is on the representation of dam and reservoir operation in catchment models (CMs) and Land Surface Models (LSMs), particularly when used at large scales.While there has been tremendous progress in the last decades in modelling the operation and management of reservoir systems at local to regional scales (e.g., Razavi et al. 2012;Asadzadeh et al. 2013;Guo et al., 2013;Castelletti et al., 2010;Chang et al., 2010;Fraternali et al., 2012;Kasprzyk et al., 2013), a gap still exists between the methodologies applied for local/regional-scale reservoir operation and management and the representation of reservoir operation in Earth systems models, particularly in LSMs.This gap is due to a two-fold challenge.First, the upscaling of methodologies used at smaller scales to larger scales is non-trivial; and second, the availability of data on reservoir operation and water use is often limited in many parts of the world.For example, the reservoir purpose and operational details are not always known and large reservoirs typically serve several purposes (Wisser et al., 2010).As a result, most current hydrological modeling activities with CMs and LSMs, if not all, offer only a limited capability in simulating reservoir operations, whereas reservoir operation in practice involves a complex set of human-driven processes and decisions.
The existing reservoir operation methods in hydrologic models can be categorized roughly into three groups, based on their level of complexity in representing flow regulation; (I) natural lake methods, (II) inflow-and-demand based methods, III) target storage-and-release based methods.
Methods in the first category use formulations developed for the simulation of natural lakes or uncontrolled reservoirs.In these methods, the downstream release is calculated as a function of reservoir storage characterized by some empirical parameters (Meigh et al., 1999;Döll et al., 2003;Pietroniro et al., 2007;Rost et al., 2008).For instance, Meigh et al. (1999) calculate the release by   =   1.5 where   and   are release and reservoir storage, respectively.Their method was later modified by Döll et al. (2003) such that where  1 and  2 are release coefficients, and   and   are minimum and maximum allowable reservoir storages.The advantage of this method, as shown in Döll et al. (2003), is its minimal data requirement, which supports its global applicability to model lakes, reservoirs and wetlands.
However, it has limited functionality to adequately represent managed reservoirs due to not accounting for reservoir operation policies to constrain or increase releases at different phases of reservoir storage dynamics.Such simplistic methods ignore the fact that the operation of a reservoir depends on the reservoir purpose and the seasonal pattern of the mismatch between the demands it supports and the inflow it receives.
The inflow-and-demand based methods determine reservoir release using inflow or a combination of inflow and demands.
The simplest method assumes the downstream release is equal to inflow minus the sum of seepage and evaporation loss, i.e., there is no change in reservoir storage over time.This method is, for example, one of the reservoir operation options of SWAT (Arnold et al., 1998).Obviously, this method is handicapped as it neglects the reservoir storage dynamics.The method used in Wisser et al. (2010) is more elaborate as it estimates the release as a function of mean annual inflow and a set of empirical parameters that can be calibrated in the absence of information on the actual operation of a reservoir.Hanasaki et al. (2006) devised a method that accounts for water withdrawals for both irrigation and non-irrigation demands as well as inflows and reservoir maximum storage capacity.Hanasaki et al. (2006) showed that this method improved monthly discharge simulations compared with the results of a natural lake method.Haddeland et al. (2006a) further enhanced this method by constructing an objective function to be optimized for optimal release to satisfy demands of different sectors including flood control.
While the above methods have limitations in accounting for details of reservoir operation, they have found wide applicability in several global hydrological models (Hanasaki et al., 2008a, b;Pokhrel et al., 2012a;Döll et al., 2009;Biemans et al., 2010;Voisin et al., 2013).Nonetheless, for an adequate representation of reservoirs, particularly multi-purpose reservoirs and/or those with multi-year carry-over capacity, it is important to consider reservoir zoning and adjust reservoir release formulations for different storage levels; the absence of this consideration may limit the capability of this group of methods in representing complex reservoir operations (Wu and Chen 2012).
The target storage-and-release based methods (the third group mentioned above) aim to emulate actual rule curves (i.e., reservoir target storage and release for different times of the year) that guide reservoir operators to decide on downstream releases (Burek et al., 2013;Yates et al., 2005;Neitsch et al., 2011).The target levels of storage divide the total reservoir storage capacity into multiple zones.For example, in the SWAT model (Arnold et al., 1998), a reservoir model is available in which the total storage of a reservoir is divided into sediment, principal, flood control, and emergency flood control zones.
Wu and Chen (2012) modified this approach by changing its reservoir zoning model and developed a reservoir release simulation strategy that uses a decision-based parameterization that may better fit both storage and release of multi-purpose reservoirs.However, they reported only one application of this strategy to a local-scale reservoir, and its comprehensive evaluation needs to be performed on other reservoirs in other regions with different climates, levels of regulation, and allocation objectives.
For the LISFLOOD model, Burek et al. (2013) divided the total reservoir storage into conservative, normal and flood zones and defined release in accordance with these storage zones using multiple-linear regression.Zajac et al. (2017) showed the applicability of this method to capture the effects of lakes and reservoirs globally using a parameterization that depends on naturalized inflow and maximum storage.Their results showed that the inclusion of reservoirs and lakes in a hydrologic model through this method helped improve streamflow simulation for many stations, but the performance in replicating observed storage dynamics was not reported.Overall, the primary advantage of methods in the third group is that they allow approximation of reservoir-release policy and have the potential of making use of detailed data on a reservoir.Their main limitation, however, is their relatively high data demands.When data are available, methods under this category have the potential to enhance the representation of dams and reservoirs in terms of both reservoir storage and release, while adapting to the seasonality and change in operations on different time scales from daily to seasonal.Given the advances in the field and the growing availability of data sources, the target storage-and-release methods seem to be the most promising, as they can better simulate the reservoir operation dynamics (the dynamics of both storage and release).The data requirement includes data on observed inflow, observed release, observed storage (level) and reservoir physical characteristics.Reservoir level data are available for most lakes and reservoirs in the public domain, particularly in North America.These data can be converted to reservoir storage using reservoir elevation-area-volume relationships or by using area-volume relationships approximated by regular geometric shapes (Yigzaw et al., 2018;Liebe et al., 2005;Lehner et al., 2011).Inflows to and releases from a reservoir can be approximated by streamflow stations located upstream and downstream of the reservoir, respectively.Further, satellite missions such as MODIS (Savtchenko et al., 2004) and satellite radar altimetry are providing information on lake and reservoir surface area dynamics and reservoir water elevation for some large reservoirs.The combination of MODIS and satellite radar altimetry allows to derive storage-area-depth relationships (Gao et al., 2012;Andreadis et al., 2007;Zhang et al., 2014;Yoon and Beighley, 2015).The planned SWOT ( 2021) mission (Garambois and Monnier, 2015) will increase the availability of water level data for smaller rivers (with widths going down to 100m) that can be potentially converted to discharge to estimate reservoir inflows and downstream reservoir releases.
This study is an attempt to develop an improved method that better emulates reservoir operations for large-scale hydrologic modelling applications in terms of both reservoir storage and release, as outlined in section 2.

Material & Methods
In this section, we present the characteristics and formulation of our reservoir model.The reservoir water balance is maintained using the continuity equation, as shown in finite difference form in equation 1.The aim is to estimate unknown storage   and release   at the current time step based on the storage at the previous time step  −1 and precipitation (P), evaporation (E) and inflows (I) during the current time step.This equation is solved in conjunction with the parametrization equations presented in the next section for reservoir releases to compute   and   .

Proposed Reservoir Operation Model
A detailed description of our proposed target storage-and-release model (or target-release model for brevity) is provided here.This model is formulated in the form of parametric piecewise linear functions that approximate the reservoir release rules that may be used by reservoir operators.This model can be set up on any time scale; in the case studies reported here, we define the target levels to dynamically change on a monthly basis.We call the model the "Dynamically Zoned Target Release (DZTR)" Model.Piecewise linear function-based reservoir operation models have already been used to solve complex reservoir operation and water resources management problems (e.g., Razavi et al., 2013;Asadzadeh et al., 2014).A  2013), which may also be viewed as a modification to the model proposed by Burek et al. (2013) in terms of parametrization and reservoir zoning.Figure 1 shows the schematic representation of DZTR; Figure 1a shows the reservoir zoning and Figure 1b shows the piecewise-linear functions to estimate the release for each zone based on DZTR.
The DZTR model divides reservoir storage into five zones in a similar fashion to Wu and Chen ( 2012) and Burek et al. (2013), namely dead storage, critical storage, normal storage, flood storage and emergency storage.Whenever storage is below the emergency storage zone, release only occurs through the bottom outlet, but when the storage is within that zone, release happens through both of bottom outlet and the spillway.In the absence of data, the dead storage (Zone 0) is assumed to be 10% of the maximum storage after Döll et al. (2009).To estimate the remaining storage zones in cases where no operational information on a reservoir is available, we propose two strategies: (1) setting the zones based on some suggested exceedance probabilities based on historical reservoir storage time series, (2) optimizing these zones to reproduce the observed storage and release time series.Target releases for each zone can be obtained in a similar fashion.These target storages and releases are allowed to vary each month (or on any other arbitrarily selected time resolution) to allow a better representation of the seasonality of reservoir operation.
When reservoir storage is within the dead storage zone (Zone 0), the reservoir release is zero (equation 3).In Zone 1 (critical storage zone), the reservoir release is a function of storage at a given time step and the critical release target value (equation 4).In this zone, the reservoir operates to avoid storage depletion while trying to support environmental flow requirements defined as a critical (or minimum) release.In Zone 2 (normal storage zone), the reservoir release is purely governed by reservoir storage and varies between critical and normal release targets (equation 5).In this zone, the downstream release is greater for higher levels of storage.In Zone 3, the release decision considers both reservoir storage and inflow in that time step as well as the normal and maximum release targets (Equation 6).When in this zone, two scenarios may occur: (A) the amount of inflow in a time step is equal to or less than the normal release rate; (B) the amount of inflow in this time step is greater than the normal release rate.As formulated in Equation 6, in the case of scenario B, the inflow rate comes to play to augment the release in an attempt to keep the reservoir level within the normal storage zone.Scenario B is expected to occur more frequently in smaller reservoirs that only have "within-year" storage capacity, while scenario A should be more commonly seen with larger reservoirs that have "multi-year" carry over capacity.Hanasaki et al. (2006) suggested that reservoirs that have a ratio of storage capacity to mean annual inflow (referred to as c) of less than 0.5 be assumed as withinyear reservoirs and the ones with a ratio of 0.5 and above be considered as multi-year reservoirs.Other values for this threshold were also suggested in the literature; e.g.Wu and Chen (2012)  Zone 0 Zone 1 Zone 2

Evaluation Criteria
We evaluated the performance of the proposed reservoir operation model in emulating the outflow and storage data collected for many reservoirs around the world.As this model was intended to be integrated into large-scale H-LSMs, we further evaluated it when embedded in the MESH model (Modélisation Environmentale-Surface et Hydrologie) (Pietroniro et al., 2007).For all of these evaluations, we used Nash-Sutcliffe Efficiency (NSE) (Nash and Sutcliffe, 1970) and Kling-Gupta Efficiency (KGE) (Gupta et al., 2009), metrics to assess the goodness of fit of the model to observed reservoir outflow and storage data.

Identification of Reservoir Operation Model Parameters
As demonstrated in Section 3. every reservoir, we use cumulative distribution functions (CDFs) of historical storage and release values; see Figure 2 for example CDFs of the Lake Diefenbaker reservoir (Gardiner dam) in the Saskatchewan River Basin, Canada.Our preliminary analysis indicated that target storage and release values corresponding to 10%, 45%, and 85% non-exceedance probabilities generally perform reasonably well.We call these our 'generalized parameterization'.
However, optimal values of parameters for a given reservoir can be identified, when data are available, through parameter identification techniques.For this purpose, we use a bi-objective optimization approach that begins with the generalized parameter values as the starting point and optimizes the model fit to both storage and release data simultaneously.The first objective function measures the performance in reproducing releases and the second objective function measures the storage dynamics performance.Section 3.2 mentioned the performance metrics used.
where x is a vector of decision variables (parameter values) and Ω is decision space.For parameter identification on a monthly basis, a total of 72 decision variables are used in the optimization.We choose rather arbitrarily the storage and release target intervals that correspond to [5-35%], [35-75%], [75-95%] non-exceedance probabilities as the ranges of variation for critical, normal, and maximum (flood) storage and release, respectively.
For  1 () and  2 () in Equation 9, we used NSE(flow) that measures the goodness-of-fit to observed downstream flows and NSE(Storage) that measures the goodness-of-fit to observed reservoir storage.The multi-objective optimization problem to calibrate 72 reservoir target release and storage parameters was conducted using AMALGAM evolutionary multiobjective optimization algorithm (Vrugt and Robinson, 2007).AMALGAM was selected because it provides effective and reliable solutions for multi-objective optimization using multiple search operators (genetic algorithm, particle swarm optimization, adaptive metropolis search, and differential evolution) and self-adaptive offspring creation.Vrugt et al. (2009), Wöhling and Vrugt (2011);Zhang. (2011); Raad et al. (2009); Dane et al. (2010) and others showed that the performance of AMALGAM model parameter calibration was better than or equivalent to some other calibration algorithms across different complex response surfaces.Further detail on AMALGAM is provided in Vrugt and Robinson (2007).AMALGAM was run using an initial population size of 100, resulting in a total of 15,000 model evaluations to estimate final Pareto solutions for every single reservoir.

Comparison of reservoir operation models
We compared the performance of our DZTR model against those of Hanasaki et al. (2006) and Wisser et al. (2010) using NSE and KGE performance metrics defined on both storage and release simulations.The comparisons were made only for selected non-irrigation reservoirs because their irrigation reservoir formulation requires additional data on water demands.
For the method of Wisser et al. (2010), reservoir release was estimated under two conditions as shown in Equation 10and storage dynamics are estimated from the reservoir water balance (Equation 11).where  and  are empirical constants set to 0.16 and 0.6 respectively and   is the mean annual inflow (m 3 /s) and   is inflow to the reservoir (m 3 /s) at time t.
In the method of Hanasaki et al. (2006), the release of non-irrigation reservoirs was estimated by multiplying the mean annual inflow with release constraining coefficients (Equation 12).The release coefficient is estimated by dividing the initial storage of a given operational year by the maximum storage and the value is calculated for each simulation year (equation 13).The start of the operational year is the month when the mean monthly inflow shifts from being greater than the mean annual inflow to being lower than the mean annual inflow.
where  is the ratio of maximum reservoir storage to the mean total annual inflow; and  , is the release coefficient;  , ′ is the provisional monthly release (m 3 /s) which is equal to mean annual inflow (m 3 /s);  is a non-dimensional constant set to 0.85.Equation 12 differentiates between multi-year and single year storage reservoirs based on a threshold value of 0.5 for c.
Currently, the reservoir representation in MESH model is rudimentary.MESH offers two approaches to account for reservoir operation.In the first approach, the observed reservoir release rate at the reservoir location is provided as input to the model.
In this approach, the flow from the catchment upstream of the reservoir is discarded as the release is replaced by observations, a process referred to as "streamflow insertion", which limits the utility of the model to simulate future scenarios for which releases are not yet known.This approach violates the water conservation law in the model and also creates discontinuities with the model setup, especially if there are reservoir cascades.Nevertheless, streamflow insertion could be used when coupling water management models with MESH, and these coupled models could be used to formulate scenarios for reservoir operations.As mentioned in the objectives, however, model coupling is not the focus of this study as we are looking to examine the internal representation of reservoir operations within CMs and LSMs.The second approach is a natural lake or uncontrolled reservoir representation model similar to that of Döll et al., (2003), which was shown to be unsuitable for highly managed reservoirs.To improve the reservoir representation in MESH, this study aims to incorporate the DZTR model for controlled reservoirs into the MESH framework and evaluate its performance.

Case Studies and Data
The data set required to build and evaluate a reservoir operation model includes ( 1) reservoir physical characteristics such as the volume-level-area relationship and maximum capacity, which are static (in the absence of sedimentation or dam heightening), (2) time series of hydrologic variables such as inflow, release, and water level (or storage), and environmental flows.In this study, we assembled such a dataset for 37 reservoirs located in several regions across the globe (Figure 3).Most of these are located in the western US and western Canada; while some are located in Vietnam, central Asian countries and Egypt.Table 1 provides a summary of reservoir locations, construction years, main purposes, data periods, and other dam characteristics.Measured inflow, release, and storage time series were collected from different sources.Most of them were provided by the authors of previous studies (Hanasaki et al., (2006) and Coerver et al., (2018)).The inflow, release, and storage data for reservoirs located in Canada were acquired from Water Survey Canada, Alberta Environment and Parks, and the Saskatchewan Water Security Agency.Data for the High Aswan dam were acquired from the Nile Basin Encyclopedia via the Nile Basin Initiative.Additional information about the degree of regulation, dam height, and catchment area was obtained from GRanD database (Lehner et al., 2011).Reservoir operation simulations were performed on daily and monthly bases with simulation periods varying from 8 to 62 years.The choice of simulation period and time scale was based on data availability (Table 1).The first year of the reservoir simulations was used for spin-up, while the first half of the remaining data periods were used for calibration and the second half of the data periods were used for model validation.
We also evaluated the integration of our reservoir model into the MESH model on six reservoirs in two major bsins in Western Canada.Six of the test reservoirs (Gardiner, St Mary, Waterton, Oldman, Ghost and Dickson dams) are located in the heavily regulated Saskatchewan River basin (SaskRB) an done reservoir (Bennet dam) is located in the Mackenzie River Basin (MRB).For both of the basins, the MESH model was set up on a grid resolution of 0.125° and the data required to build the MESH model were obtained from different sources.The topographic data are based on the Canadian Digital Elevation (CDED) at a scale of 1:250,000 and were obtained from the GeoBase website (http://www.geobase.ca/).The data on seven climate forcing variables at a 30-min temporal resolution were obtained from Global Environmental Multi-scale (GEM) NWP model (Côté et al., 1998) and Canadian Precipitation Analysis (CaPA) (Mahfouf, et al., 2007).comparison against other reservoir operation approaches and a "base case" where the existence of a reservoir is ignored in a model (we refer to this base case as the "no-reservoir assumption").Under the no-reservoir assumption, the release was considered equal to inflow, and storage was considered constant, and as such, the performance metrics were estimated by directly comparing inflow with observed release.However, for reservoirs with c>0.5 such as Bhumibol, Flaming Gorge, Fort Peck, High Aswan, W.A.C. Bennett, the NSE (base-case) is negative, which indicates the significant influence of their regulations on the hydrograph shape.
Similarly, Figure 4b shows the evaluation of the different reservoir models based on the KGE metric (Gupta et al., (2009)).
The values of KGE (Flow) and KGE (Storage) are greater than 0.25 and 0.5 for 100% and 86% of the reservoirs, respectively.The KGE (base-case) values of 21% of reservoirs are less than 0, while those of 57% and 49% of the reservoirs are greater than 0.25 and 0.5, respectively.The NSE and KGE results show that the DZTR is capable of simulating flow and storage simulation well.Figure 6 compares the reservoir simulation and observation time series for the whole simulation period, while Figure 7 shows the long-term average of these simulations.Inflows are also included in Figure 6 and  (1984-1986, 1988-1989, 1999-2004) and wet years (1993,2005,(2010)(2011).Furthermore, as expected, Figure 7 shows that lowly-regulated reservoirs (c <0.5) have less impact on the flow regime, but with fairly significant storage seasonality (Oldman, E.B. Campbell, Palisades, Andijan).In general, the DZTR model with the generalized parameterization of reservoir zones and releases showed an improved performance and can be applied to any hydrological model (CM or H-LSM) that involves reservoir simulation.

A Comparison with previously developed reservoir operation models
To further illustrate the reliability of DZTR model in representing reservoir simulation, a comparison with the methods of Hanasaki et al, (2006) and of Wisser et al., (2010)  In addition, we compared the DZTR result shown in Figure 7 and 8 with the results reported in a study by Coerver et al. (2018) who applied a fuzzy logic concept to extract 11 reservoir operating rules.This comparison showed that the performance of our generalized parameterization is comparable with that Coerver et al. (2018) in simulating reservoir release; note that performance on storage is not reported in Coerver et al. (2018).This indicates that the simple parameterization applied in the DTZR model can provide effective solutions comparable to those of Coerver et al., (2018).

Initial storage and inflow sensitivity test
To test the effect of initial storage used in the reservoir simulation performance, two experiments were conducted on three reservoirs with different scale of regulations 1) Charvak (c=0.28), 2) Gardiner (c=1.46), and 3) High Aswan (c=2.84).In the first experiment, the initial storage was allowed to vary between within ten percent of maximum storage (0.1 *   ) to maximum storage (   ).In the second experiment, the initial storage range was narrowed to starting simulation month minimum and maximum historical observations.In both tests, 150 simulations were conducted by sampling the initial storage using random sampling from the defined storage range.Overall, as expected, the experiments suggest that the effect of initial storage on reservoir simulation performance depends on the regulation scale.Starting from observed storage values and using a oneyear warm-up period allows stabiliztion of the initial storage effect for low and medium regulated reservoirs.However, for highly regulated reservoirs, as in the case of High Aswan, longer spin-up periods are needed to stabilize the simulations.For example, a five-year spin-up period was required to fully stabilize the performance for the High Aswan dam simulations.
The existence of inflow bias is inevitable in any hydrological modeling practice.To understand the behavior of the DZTR model under biased inflow conditions, we conducted a sensitivity experiment on the Charvak, Gardiner and High Aswan reservoirs.To do so, the DZTR model performance was tested using five simulations in which the entire inflow time series was changed by -50%, -25%, 0%, +25%, and +50%.The sensitivity of simulations to bias in inflow was evaluated using the NSE (Flow) and NSE (Storage) performance metrics.
Figure 10 shows the results of the inflow bias test and that the reservoir simulation performance significantly changes as a result of this bias.Reducing the inflow by 50% considerably reduced the reservoir storage and release and led to negative values of NSE (Flow) and NSE (Storage) for all reservoirs.For such a large negative inflow bias, the reservoir operation tries to recover the storage to the target (observed) level by releasing a low as possible.Conversely, the positive inflow bias increased simulated storage and releases for all reservoirs, which led to negative performance metrics for all reservoirs 25% and +25% showed similar behavior as -50% and +50% bias for all reservoirs, but the simulation performance metrics during -25% and +25% provide significant positive NSE values for the Charvak and Gardener dams except for the Gardiner NSE(Flow) for +25% which resulted a negative NSE value.However, on the highly regulated High Aswan dam, the ±25% inflow bias significantly reduced the performance to negative values.

Parameter calibration and validation of the DTZR model
We tried to improve upon the generalized parameterization by calibrating the DTZR parameters via bi-objective optimization for two objective functions, Nash Sutcliffe on reservoir storage (NSE (Storage)) and Nash Sutcliffe on reservoir release (NSE (Flow)).This is an important step when the data and computational resources for optimization are available, to enhance reservoir simulation and consequently hydrological modeling of the region of interest.This general improvement of flow simulation when comparing a reservoir model to the no-reservoir assumption is, of course, not surprising.What is important to note, however, is that the improvement in NSE can be dramatic without calibration of the DZTR parameters.This is important for many LSM applications where calibration is generally not performed.Hanesaki et al. (2006) illustrate that their method is superior to the natural lake (or unregulated reservoir) method applied in many CMs and H-LSMs, and this paper shows that the DZTR model improves upon the results of Hanesaki et al. (2006).Therefore, it is natural to assume that the DZTR model would also be an improvement in uncalibrated H-LSM applications.
However, calibration is very common in CM or H-LSM applications in which the DZTR model would likely be employed.
A full comparison of calibrated results between a no-reservoir case, natural lake (or unregulated reservoir), and the DZTR model (and the other reservoir models) is beyond the scope of this paper.Again, given the improvements shown with the uncalibrated DZTR model when compared with other uncalibrated models, and the general improvements shown here when calibrating the DZTR model, it is assumed that calibrating the DZTR model within a CM or H-LSM would improve upon calibrating an unregulated reservoir model, or the other reservoir models compared in this paper.
The storage simulation showed low NSE (Storage) value for St. Mary and Waterton dams and negative NSE (Storage) for Oldman and Ghost dams.However, the simulation showed a reasonable representation of storage variability, but with considerable underestimation.This underestimation in storage in Figure 12 is attributable to the fact that the modelled inflow is underestimated.It is expected that calibration of the land-surface parameters in conjunction with the DZTR parameters in MESH would improve the modelled inflows and resulting modelled reservoir storage.
It is worth mentioning again that H-LSMs, such as MESH, can also be used for the original purpose of LSMs, which is to represent fluxes from the land-surface to the atmosphere.If the approach improves modelled flows where reservoirs operate, it could result in a better parameterization of the LSM, which should in-turn improve land-surface ET.Improving the land-  As expected, calibration significantly improved the performance of the DZTR model compared with the performance of the "generalized parameters".There often exists, however, a significant tradeoff between fitting reservoir storage versus release, signifying the importance of accounting for both storage and release in a multiobjective fashion.
 The integration of the DTZR reservoir model into the MESH hydrology-land surface modelling system was straightforward and improved the overall model performance compared with the traditional methods of accounting for reservoirs in H-LSMs.This integration can be viewed as a successful example for improving the representation of reservoir operation in CMs, LSMs and GWSMs.
Future research work may include (1) examining the applicability of the DTZR model for regions with severely-limited data by examining the utility of other data sources such those derived from satellite-based observations (Savtchenko et al., 2004;Garambois and Monnier, 2015;Gao et al., 2012) and using area-volume relationship approximated by regular geometric shapes (Yigzaw et al., 2018); and (2) an examination of a direct one-and/or two-way coupling of WMMs with CMs and LSMs towards developing a seamless coupled framework for the simulation of natural-engineered watershed systems.
Hydrol.Earth Syst.Sci.Discuss., https://doi.org/10.5194/hess-2019-7Manuscript under review for journal Hydrol.Earth Syst.Sci. Discussion started: 16 January 2019 c Author(s) 2019.CC BY 4.0 License.systematic integration of such models into large-scale hydrological modeling has been reported inBurek et al. (2013) as implemented in the LISFLOOD hydrological model, and inNeitsch et al. (2012)  as implemented in the SWAT model.Our DZTR model is a generalization of the method developed byRazavi et al. ( used a c value of 0.3.In this study, scenario A is used for reservoirs that have multi-year capacity and scenario B for reservoirs that have within-a-year capacity.Lastly, in Zone 4 (emergency storage zone) the reservoir algorithm operates to avoid reservoir overtopping by releasing the larger of the maximum release target or all excess storage above the maximum storage value (the flood storage) constrained to the Hydrol.Earth Syst.Sci.Discuss., https://doi.org/10.5194/hess-2019-7Manuscript under review for journal Hydrol.Earth Syst.Sci. Discussion started: 16 January 2019 c Author(s) 2019.CC BY 4.0 License.downstream channel capacity   .If not specified, a rough estimate of the downstream channel capacity value could be the 99 percentile of non-exceedance probabilities of discharges from historical data.

1
The land cover data are based on 2005 land-cover map from the Canada Centre for Remote Sensing (CCRS).Soil texture data were obtained from Soil Landscapes of Canada (SLC) data of Agriculture and Agri-Food Canada.The MESH parameter values were estimated using MESH model calibration on streamflow at major subbasins of SaskRB and MRB.Hydrol.Earth Syst.Sci.Discuss., https://doi.org/10.5194/hess-2019-7Manuscript under review for journal Hydrol.Earth Syst.Sci. Discussion started: 16 January 2019 c Author(s) 2019.CC BY 4.0 License.Evaluation of the dynamically zoned target release (DZTR) model with generalized parameters Individual reservoir simulations were conducted by the DZTR model with generalized monthly storage and release parameter values set at non-exceedance probabilities of 10%, 45%, and 85%, representing reservoir storage zones and their respective target releases.The evaluation of the DZTR model is based on the performance evaluation metrics and a

Figure 4
Figure4shows performance metrics results of the DZTR model in terms of NSE and KGE for storage and release simulations compared to those of the base case.As shown in Figure4a, both NSE (Flow) and NSE (Storage) results are greater than 0.25 and 0.5 for 90% and 50% of reservoirs, respectively.Although NSE (Flow) results are greater than zero for all reservoirs, 1% of reservoirs resulted in a negative NSE (Storage) values.The no-reservoir assumption resulted in NSE (base-case) values of greater than 0.25 and 0.5 for 45% and 30% of reservoirs respectively, but when compared to NSE (Flow) of the DZTR model, a few reservoirs showed comparable performance while the majority are much better represented by the DZTR model.Under the no-reservoir assumption, 48% of the reservoirs resulted in a negative NSE (basecase).Almost all positive NSE (base-case) results were observed on reservoirs with c<0.5 such as Dickson, E.B. Campbell, Kayrakkum, Oldman and Tyuyamuyun (as explained in Section 3, c is the ratio of storage capacity to annual inflow volume).

Figure 5
Figure5shows the scatter plot between KGE and NSE along with the regulation scale represented by c.The orientation of the scatter plot between NSE and KGE on flow and storage shows a strong positive correlation between the evaluation metrics which indicates that both metrics provide somewhat similar evaluation information.The relation between regulation scale and evaluation metrics shows that the performance in simulating the release improves in most cases as the scale of reservoir regulation reduces, as shown in Figure5aand 5b.However, as depicted in Figure5c, the regulation scale and the storage simulation performance (in terms of both KGE (Storage) and NSE (Storage)) do not show a strong correlation, which indicates that the performance for storage is not strongly affected by regulation scale.
Figure 7 to show the regulation pattern and changes caused by reservoir operation.Both figures indicate that the DZTR model captures well both release and storage dynamics in terms of reproducing the daily and monthly seasonality and in terms of both the magnitude and timing of storage and releases for almost all reservoirs, especially for reservoirs with high regulation (multipurpose, multiyear reservoirs) such as American Falls, Bhumibol, High Aswan, Sirikit, Trinity, and Bennett dams.However, the simulations also show some systematic over-and under-estimations; for example, the simulations of Bhumibol, Fort Peck, High Aswan, Int.Falcon, Navajo, Bennet, and Int.Amistad reservoirs show continuous underestimation and overestimation of reservoir storage.Some reservoirs such as Trinity, Palisades, Kayrakkum, Flaming Gorge, and Garrison show underestimation and overestimation of reservoir storage only for some seasons.A closer look at American Falls, Flaming Gorge, Fort Peck, Glen Canyon, Navajo dams in Figure 6 indicates that the DZTR model reliably captured storage and release seasonality, interannual trends, and release pattern shifts during consecutive wet years 1982-1986 followed by consecutive dry years 1987-1993.Similar patterns can be observed for the Gardiner dam with good simulation during both dry years was conducted as shown in Figure 8.The comparison shows that the DZTR model provides a considerable improvement according to all of the performance criteria, notably NSE (Storage) and NSE (Flow), except in the case of the E.B. Campbell dam where Hanasaki et al.'s method showed similar performance to DZTR.Also, the method of Hanasaki et al. (2006) outperforms that of Wisser et al. (2010).Out of the thirteen reservoirs compared, the DZTR resulted in positive values for both NSE (Storage) and NSE (Flow) in all except for E.B. Campbell storage.However, the method of Hanasaki et al. (2006) resulted in only five positive NSE (Flow) and NSE (Storage) values, whereas that of Wisser et al, (2010) generated only three positive NSE (Flow) values while producing negative values for NSE (Storage) for all the reservoirs compared.A similar performance pattern was observed for KGE metrics for flow and storage.

Figure 9
Figure 9 shows the results of these initial storage perturbation experiments.For both experiments the simulations on the Charvak dam showed a similar range for NSE (Flow) [0.79, 0.83] and NSE (storage) [0.61, 0.74].Using one year as a spinup period on Charvak dam simulations stabilized the initial storage effects, resulting in NSE (Flow) of 0.82 and NSE (Storage) of 0.74.The simulations on Gardiner dam in the first experiment showed a range of [0.35, 0.51] for NSE (Flow) and [-0.43, 0.88] range for NSE (storage), while in the second experiment the ranges were narrowed to [0.44, 0.49] for NSE (Flow) and [0.87, 0.88] for NSE (storage).For a one year spin-up period on the Gardiner dam this simulation converged the NSE (Flow) range to [0.49, 0.51] and the NSE (Storage) range to [0.76, 0.87] in the first experiment and to 0.49 NSE (Flow) and 0.87 NSE (Storage) for the second experiments.On the other hand, the simulation on the High Aswan dam showed a range of [-0.28, 0.85] for NSE (Flow) and [0.38 0.91] for NSE (storage) for the first experiment and [0.52, 0.85] for NSE (Flow) and [0.42 0.91] for NSE (storage) for the second experiment.Excluding a one year spin-up period from the metric calculation on the High Aswan dam simulation narrowed the NSE (Flow) range to [0.62, 0.85] and the NSE (Storage) range to [0.58 0.91] for both experiments.Overall, as expected, the experiments suggest that the effect of initial storage on Hydrol.Earth Syst.Sci.Discuss., https://doi.org/10.5194/hess-2019-7Manuscript under review for journal Hydrol.Earth Syst.Sci. Discussion started: 16 January 2019 c Author(s) 2019.CC BY 4.0 License.except on Gardiner NSE (Storage).As shown in Figure 10, with large positive inflow bias, storage quickly moves towards flood and maximum storage targets resulting in insufficient storage left to attenuate flood peaks and the operation model starts discharging large releases through the spillway to maintain the storage at the maximum storage target.Inflow bias of - Figure11shows the multicriteria reservoir calibration (yellow circles) and validation (red circles) Pareto solutions for all reservoirs.The Pareto solutions show strong tradeoffs between fitting observed reservoir storage versus downstream release, which also reflects the fact that the problem is multi-objective by nature and it is required to consider both storage and release, instead of fitting one at the cost of degrading the other.The generalized parameterization solution for the calibration (Yellow Square with blue border) and validation periods (red square with blue border) is also added in Figure11for each reservoir to show the improvement gained through parameter calibration.Relative to the generalized solution for the calibration period, reservoir parameter calibration improved both NSE (Flow) and NSE (Storage) for all reservoirs with a median improvement of 0.11 and 0.21, respectively.The NSE (Flow) improvement ranged from 0.017 to 0.575, and NSE (Storage) improvement ranged from 0.02 to 0.66.The Parameter calibration has shown significant improvement on reservoirs that have lower performance with generalized parameterization.The best examples of this case are Fort Randall, Int.Amistad, Trinity, Int.Falcon, and E.B. Campbell, as shown in Figure11.Small improvements in performance have also been observed on reservoirs that have greater performance with generalized parameterization such as American Falls, Andijan, Nurek, High Aswan, Waterton, and Charvak.The validation of calibrated solutions improved the NSE (Flow) and NSE (Storage) for 56% of the reservoirs with a median improvement of 0.035 and 0.092, respectively.The NSE (Flow) improvement in the validation period ranged from 0.001 to 0.335, and NSE (Storage) improvement ranged from 0.004 to 1.02.During validation, the remaining reservoirs (44% of them) resulted in NSE (Flow) and NSE (Storage) reductions with a median reduction of 0.032 and 0.089, respectively.The reductions of NSE (Flow) ranged from 0.001 to 0.073, and those of NSE (Storage) ranged from 0.001 to 0.257.Overall, considerable improvement was observed both in calibration and validation periods for several reservoirs such as the Dickson, Gardiner, Ghost, Int.Amistad, Int.Falcon, Kayrakkum, Sirikit, Yellowtail, and Glenmore.However, the improvements of DZTR model performance during calibration do not guarantee performance improvement in validation.Hydrol.Earth Syst.Sci.Discuss., https://doi.org/10.5194/hess-2019-7Manuscript under review for journal Hydrol.Earth Syst.Sci. Discussion started: 16 January 2019 c Author(s) 2019.CC BY 4.0 License.4.5 DZTR model test in MESH model Finally, the generalized parametrization of the DZTR model was integrated into the MESH model and tested to simulate six reservoirs in the Saskatchewan River Basin (Gardiner, St Mary, Waterton, Oldman, Ghost and Dickson dams) and one reservoir (Bennet dam) in the Mackenzie River Basin, both in Western Canada.The reservoir simulation was run using MESH modeled inflows at a half-hourly time step, the same as the MESH modeling time step, and the performance metrics were calculated at a daily time step.The MESH modelled inflows are considered the natural flows that would have occured if there had been no reservoirs.

Figure 12
Figure 12 illustrates that the uncalibrated DZTR model generally improves upon having no representation of the reservoirs in the model.This improvement is apparent in the NSE values of the flow, which increase with the DZTR model.The only exception is Dickson dam with a small reduction in NSE.The importance of integration of the DZTR model was predominant for the Gardiner and Bennett dams, which are highly regulated when compared to the other reservoirs tested in MESH.
Hydrol.Earth Syst.Sci.Discuss., https://doi.org/10.5194/hess-2019-7Manuscript under review for journal Hydrol.Earth Syst.Sci. Discussion started: 16 January 2019 c Author(s) 2019.CC BY 4.0 License.surface parameterization requires that the reservoir flows and/or storages are used to inform the land-surface parameterization through parameter estimation techniques.This is as an area requiring future research.5 Summary and Conclusions Human interventions in hydrologic systems through dams and reservoirs significantly change the flow regime.In this paper, we presented an improved reservoir operation model, called Dynamically Zoned Target Release (DZTR) that can be integrated into any large-scale hydrological model; here we integrated DZTR into the MESH hydrology-land surface model.The DZTR model is based on parametric piecewise linear functions that approximate reservoir release rules that may be used by reservoir operators.We proposed two strategies to identify the parameters of this model: one based on the distributions of historical storage and release to generate the so-called "generalized parameters" and the other one based on direct calibration to observed storage and release time series via multi-objective optimization.We first tested the DZTR model individually across a number of reservoirs around the globe, and then tested its performance when plugged into the MESH model.Our conclusions include: The DZTR reservoir operation model performed well in reproducing observed storage and release time series in (almost) all reservoirs tested and outperformed the existing reservoir models proposed byHanasaki et al. (2006) andWisser et al. (2010).The model was capable of capturing inter-and intra-annual variability of both reservoir storage and release.

Figure 1 :
Figure 1: The schematic representation of reservoir zoning and storage-release function: a) Four (active) reservoir zones with inflow and outflows; b) piecewise linear reservoir release function, m1, and m2 control the slope of the release curve and they change monthly.The upward blue arrow is to indicate that inflow to the reservoir may also be considered in determining the release in zone 3. 5

Figure 2 :
Figure 2: Cumulative Distribution Function (CDF) a) Storage CDF of Gardiner dam b) Reservoir release CDF of Gardiner dam.Double arrows on y-axis shows parametrizations ranges for each generalized parameters.

Figure 3 :
Figure 3: Locations of dams used to evaluate reservoir routing model.

Figure 4 :
Figure 4: Performance evaluation result of the DZTR model reservoir operation algorithm.

Figure 5 :
Figure 5: Scatter plot between KGE and NSE with regulation scale represented in terms of c a) KGE and NSE on no reservoir condition, b) KGE and NSE on DZTR release, and c) KGE and NSE on DZTR storage.

Figure 6 :
Figure 6: Daily and monthly reservoir simulations using DZTR model with a generalized parametrization, x-axis shows month/year, the primary y-axis shows release (m 3 /s) and the secondary y-axis shows storage (m 3 ).

Figure 7 :
Figure 7: Long-term average daily or monthly reservoir simulations with generalized parametrization, the x-axis shows days (1-365) or months, (1-12) the primary y-axis shows release (m 3 /s) and the secondary y-axis shows storage (m 3 ).

Figure 8 :
Figure 8: A comparison of our proposed reservoir operation model with generalized parameters with the models of Hanasaki, et al. (2006) and Wisser et al. (2010).

Figure 9 :
Figure 9: Reservoir initial storage effect on storage and release simulation.5

Figure 10 :
Figure 10: Inflow bias sensitivity test on storage and release simulation.

Figure 12 :
Figure 12: Reservoir simulation results within MESH model run for selected reservoirs.X-axis shows time (days), the primary yaxis shows release (m 3 /s) and the secondary y-axis shows storage (m 3 ).