Irrigation, damming, and streamflow fluctuations of the Yellow River

The streamflow of the Yellow River (YR) is strongly affected by human activities like irrigation and dam operation. Many attribution studies have focused on the long-term trends of streamflows, yet the contributions of these anthropogenic factors to streamflow fluctuations have not been well quantified with fully mechanistic models. This study aims to (1) demonstrate whether the mechanistic global land surface model ORCHIDEE (ORganizing Carbon and Hydrology in Dynamic EcosystEms) is able to simulate the streamflows of this complex rivers with human activities using a generic parameterization for human activities and (2) preliminarily quantify the roles of irrigation and dam operation in monthly streamflow fluctuations of the YR from 1982 to 2014 with a newly developed irrigation module and an offline dam operation model. Validations with observed streamflows near the outlet of the YR demonstrated that model performances improved notably with incrementally considering irrigation (mean square error (MSE) decreased by 56.9 %) and dam operation (MSE decreased by another 30.5 %). Irrigation withdrawals were found to substantially reduce the river streamflows by approximately $\mathrm{242.8}\pm \mathrm{27.8}\times {\mathrm{10}}^{\mathrm{8}}$ m³ yr⁻¹ in line with independent census data ($\mathrm{231.4}\pm \mathrm{31.6}\times {\mathrm{10}}^{\mathrm{8}}$ m³ yr⁻¹). Dam operation does not change the mean streamflows in our model, but it impacts streamflow seasonality, more than the seasonal change of precipitation. By only considering generic operation schemes, our dam model is able to reproduce the water storage changes of the two large reservoirs, LongYangXia and LiuJiaXia (correlation coefficient of ∼ 0.9). Moreover, other commonly neglected factors, such as the large operation contribution from multiple medium/small reservoirs, the dominance of large irrigation districts for streamflows (e.g., the Hetao Plateau), and special management policies during extreme years, are highlighted in this study. Related processes should be integrated into models to better project future YR water resources under climate change and optimize adaption strategies.

1 Introduction

More than 60 % of all rivers in the world are disturbed by human activities (Grill et al., 2019), contributing altogether to approximately 63 % of surface water withdrawal (Hanasaki et al., 2018). River water is used for agriculture, industry, drinking water supply, and electricity generation (Hanasaki et al., 2018; Wada et al., 2014), these usages being influenced by direct anthropogenic drivers and by climate change (Haddeland et al., 2014; Piao et al., 2007, 2010; Yin et al., 2020; Zhou et al., 2020). In order to meet the fast-growing water demand in populated areas and to control floods (Wada et al., 2014), reservoirs have been built for regulating the temporal distribution of river water (Biemans et al., 2011; Hanasaki et al., 2006), leading to a massive perturbation of the seasonality and year-to-year variations of streamflows. In the midlatitude–northern latitude regions, where a decrease of rainfall is observed historically and projected by climate models (IPCC, 2014), water scarcity will be further exacerbated by the growth of water demand (Hanasaki et al., 2013) and by the occurrence of more frequent extreme droughts (Seneviratne et al., 2014; Sherwood and Fu, 2014; Zscheischler et al., 2018). Thus, adapting river management is a crucial question for sustainable development, which requires comprehensive understanding of the impacts of human activities on river flow dynamics,, particularly in regions under high water stress (Liu et al., 2017; Wada et al., 2016).

The Yellow River (YR) is the second longest river in China. It flows across arid, semi-arid, and semi-humid regions, and the catchment contains intensive agricultural zones and has a population of 107 million inhabitants (Piao et al., 2010). With 2.6 % of total water resources in China, the Yellow River basin (YRB) irrigates 9.7 % of the croplands (http://www.yrcc.gov.cn, last access: 28 February 2021). Underground water resources are used in the YRB, but they only account for 10.3 % of total water resources, outlining the importance of streamflow water for regional water use. A special feature of the YRB is the huge spatiotemporal variation of its water balance. Precipitation is concentrated in the flooding season (from July to October) which constitutes ∼ 60 % of the annual discharge, whereas the dry season (from March to June) represents only ∼ 10 %–20 %. Numerous dams have been built to regulate the streamflows intra- and inter-annually in order to control floods and alleviate water scarcity (Liu et al., 2015; Zhuo et al., 2019). The YRB streamflows are thus highly controlled by human water withdrawals and dam operations, making it difficult to separate the impacts of human and natural factors on the variability and trends.

Numerous studies documented the effects of anthropogenic factors on streamflows and water resources in the YRB by statistical approaches (e.g., Liu and Zhang, 2002; Jin et al., 2017; Miao et al., 2011; Wang et al., 2006, 2018; Zhuo et al., 2019). To further elucidate the mechanisms, physical-based land surface hydrology models including natural and anthropogenic factors are required. Many previous model studies only considered natural processes, and YRB simulations were evaluated against naturalized streamflows (Liu et al., 2020; Xi et al., 2018; Yuan et al., 2018; Zhang and Yuan, 2020). YRB modeling studies simulating real streamflows and comparing their values to observed streamflows are scarce, the most important being from Jia et al. (2006) and Tang et al. (2008). Yet, Jia et al. (2006) prescribed census irrigation and dam operation data as input of their model. Tang et al. (2008) included irrigation as a mechanism in their DBH (distributed biosphere hydrological) model and investigated the long-term trends of streamflows, but they described the irrigation demand simply from satellite leaf area data, so that crop plant water requirements and phenology were not represented by physical laws. Several global hydrological models (GHMs) simulated both irrigation and dam operation processes and were applied for future projection of water resources regionally (Liu et al., 2019) or globally (Hanasaki et al., 2018; Wada et al., 2014, 2016). Those global GHM studies acknowledged the complex situation of the YRB where models' performances are limited, but none has focused on the sources of error or potential overlooked mechanisms in this catchment.

To model present water resources in the YRB and make future projections, not only natural mechanisms, but also anthropogenic ones must be represented in a model. If a key mechanism is missing in a model, a calibration of its parameters to match observations can compensate for structural biases, and projections may be erroneous. For example, the HBV model (Hydrologiska Byråns Vattenbalansavdelning) was well calibrated with different approaches in 156 catchments in Austria but failed in predicting streamflow changes due to climate warming (Duethmann et al., 2020), one of the key reasons being that the response of vegetation to climate change was missing in the model. In this study, we integrate two key anthropogenic processes (irrigation and dam operation) in the land surface model ORCHIDEE (ORganizing Carbon and Hydrology in Dynamic EcosystEms), which has a mechanistic description of plant–climate and soil water availability interactions and of river streamflows. Through a set of simulations with generic parameter values, we aim to preliminarily diagnose how irrigation and dam operation improve the simulations of observed YRB streamflows. After making sure we understand the impact of adding these two new and crucial processes, the model will be calibrated against a suite of observations so that it can be applied for future projections.

Using a standard version of ORCHIDEE without irrigation or dams, Xi et al. (2018) performed simulations with a 0.1^∘ hypo-resolution atmospheric forcing over China (Chen et al., 2011). They attributed the trends of several river streamflows to natural drivers from increased CO₂ and climate change and to land use change. Lacking irrigation and other human removals, their simulated results were higher than the observed streamflows for the YRB. By developing a crop module in ORCHIDEE (Wang et al., 2016; Wang, 2016; Wu et al., 2016), ORCHIDEE was able to provide precise estimation of crop physiology, phenology, and yield at both local and national scales, as well as other site-based crop models, e.g., EPICs (Folberth et al., 2012; Izaurralde et al., 2006; Liu et al., 2007, 2016; Williams, 1995), CGMS-WOFOST (de Wit and van Diepen, 2008), APSIM (Elliott et al., 2014; Keating et al., 2003), and DSSAT (Jones et al., 2003), and land surface models, e.g., CLM-CROP (Drewniak et al., 2013), LPJ-GUESS (Smith et al., 2001; Lindeskog et al., 2013), LPJmL (Waha et al., 2012; Bondeau et al., 2007), and PEGASUS (Deryng et al., 2011, 2014; Wang et al., 2017; Müller et al., 2017). ORCHIDEE-estimated irrigation accounts for potential ecological and hydrological impacts (e.g., physiological response of plants to climate change and short-term drought episodes on soil hydrology) with respect to other land surface models (LSMs) and GHMs (Hanasaki et al., 2008; Leng et al., 2015; Thiery et al., 2017; Nazemi and Wheater, 2015; Voisin et al., 2013). In a study focusing on China (Yin et al., 2020), ORCHIDEE was able to simulate irrigation withdrawals across China and to evaluate them against census data with a provincial-based spatial correlations of ∼ 0.68. It successfully explained the decline of total water storage in the YRB. In this study, we add a simple module describing the dam operations to further improve the model over the YRB.

A simple dam operation model is developed and firstly coupled to ORCHIDEE to simulate the real streamflows in this study. Similar to other GHMs and LSMs, our dam operation model is based on generic operation principles due to a lack of related data. Recent dam models are developed from different perspectives, such as the agent-based model River Wave (Humphries et al., 2014), the basin-specific model Colorado River Simulation System (Bureau of Reclamation, 2012), and the original dam module in the Variable Infiltration Capacity (VIC) model (Lohmann et al., 1998). However, the representation of dam operations in many global hydrological studies (e.g., Droppers et al., 2020; Haddeland et al., 2006, 2014; Hanasaki et al., 2018; Zhao et al., 2016; Yassin et al., 2019; Wada et al., 2014, 2016) is based on the ideas of Hanasaki et al. (2006). They categorized dams based on their regulation purposes (irrigation and non-irrigation). Irrigation-oriented rules adjust the dam retention to meet the irrigation demand downstream, while non-irrigation-oriented rules buffer floods and thus dampen the variability (Hanasaki et al., 2006). However, the water release target of a dam in the model of Hanasaki et al. (2006) is fixed at the beginning of the year and cannot adjust interactively to large intra- and inter-annual climate variations, which is a key feature of the YRB. To overcome this limitation, we propose a new dam operation model based on a target operation plan, constrained by the regulation capacity of a dam and historical simulated streamflows, with flexibility to adjust to climate variation. The effects of dams on streamflows could then be studied with ORCHIDEE and isolated from the effect of climate factors and irrigation demands. Different from classical approaches separating the YRB into an upper, middle, and lower stream (Tang et al., 2008; Zhuo et al., 2019), here we further divide both the upper and middle streams into sub-catchments based on the locations of five key gauging stations (Fig. 1). This approach splits regions with and without big dams (or large irrigation areas) in the upper and middle streams, which simplifies the assessment of the effects of irrigation and damming on streamflows.

Figure 1 (a) Map of YRB. Gray and blue lines indicate the catchment and network of YR based on GIS data, respectively. Dark circles are main artificial reservoirs on the YR. Triangles are gauging stations. Red triangles are main stations used for classifying sub-catchment and simulation comparison, and teal triangles are stations used to assess the impacts of XiaoLangDi Reservoir on river streamflows. Colored patterns are sub-catchments between two neighboring gauging stations based on the ORCHIDEE routing map. The water balances of specific sub-catchments are shown in the top left. (b) Conceptual figure of YR main stream, gauging stations, and artificial reservoirs. The sizes of circles indicate the regulation capacities of these reservoirs (Table 2).

In this study, ORCHIDEE, with the novel crop–irrigation module (Wang, 2016; Yin et al., 2020) and the new dam operation model, was applied in the YRB from 1982 to 2014 in order to (1) demonstrate whether ORCHIDEE and the dam model, with generic parameterizations, are able to improve the simulation of streamflow fluctuations and (2) attempt to separate the effect of irrigation and dams on the fluctuations of monthly streamflows. We first describe the ORCHIDEE model and our new dam model in Sect 2.1. Then we present the algorithm used for estimating sub-catchment water balances in Sect. 2.2, followed by the input and evaluation datasets, the simulation protocol, and metrics for evaluation in Sect. 2.3 to 2.5. Results are presented in Section 3, and limitations are discussed in Sect. 4.

2 Methodology

2.1 ORCHIDEE land surface model used in this study

2.1.1 Irrigation and crop modules

ORCHIDEE is a physical process-based land surface model that integrates the hydrological cycle, surface energy balances, the carbon cycle, and vegetation dynamics using two main modules. The SECHIBA (surface–vegetation–atmosphere transfer scheme) module simulates the dynamics of the water cycle, energy fluxes, and photosynthesis at a 0.5 h time interval, which are used by the STOMATE (Saclay Toulouse Orsay Model for the Analysis of Terrestrial Ecosystems) to estimate vegetation and soil carbon cycle at daily time step. The ORCHIDEE model used in this study is a special version with a newly developed crop and irrigation module (Wang et al., 2017; Wu et al., 2016; Yin et al., 2020). The crop module includes specific parameterizations for wheat, maize, and rice, calibrated over China using observations (Wang, 2016; Wang et al., 2017). It is able to simulate crop carbon allocation and different phenological stages, as well as related management (e.g., planting date, rotation, multi-cropping, and irrigation).

Irrigation amount is simulated in the land surface model ORCHIDEE (Wang, 2016; Wang et al., 2017) as the minimum between crop water requirements and the water supply. The crop water requirements are defined according to the choice of irrigation technique, namely minimizing soil moisture stress for the flooding technique, sustaining plant potential evapotranspiration for the dripping technique, and maintaining the water level above the soil surface during specific months for the paddy irrigation technique. Each crop is grown on a specific soil column (in each model grid cell), where the water and energy budgets are independently resolved. The water resources in the hydrological routing scheme are from three water reservoirs: (1) a streamflow component, (2) a fast reservoir with surface runoff, and (3) a slow reservoir with deep drainage, used in this order for defining the priorities of water use for irrigation. As long-distance water transfer is not modeled, streams only supply water to the crops growing in the grid cell they cross, according to the river routing scheme of the ORCHIDEE model (Ngo-Duc et al., 2007). Without dams, irrigation can be underestimated where dams store water to supply the crop demand. Transfer from reservoirs, lakes, or local ponds to adjacent cells is not considered, which should further lead to an underestimation of the irrigation supply, dependent on the cell size. Details of the coupled crop–irrigation module of ORCHIDEE are described in Yin et al. (2020).

2.1.2 New dam operation model

To account for the impacts of dam regulation on streamflow (Q) seasonality, we developed a dynamic dam water storage module based on only two generic rules: reducing flood peaks and guaranteeing base flow. This model depends on simulated inflows and is thus independent from irrigation demands. It has been developed for the main reservoirs of the YRB (e.g., the LongYangXia, LiuJiaXia, and XiaoLangDi in Fig. 1). Different from Biemans et al. (2011) and Hanasaki et al. (2006), we primarily consider the ability of reservoirs to regulate river flow seasonality. This means that the target base flow and flood control of our dam model are not fixed proportions of the mean annual streamflow but depend on the regulation capacity of the reservoir (C_max). Firstly, similar to Voisin et al. (2013), a multi-year averaged monthly streamflow (Q_s) is calculated based on ORCHIDEE simulations. To include the potential impacts of recent climate change on dam operation, here we only consider the latest past 10-year simulations, as

$$\begin{array}{}\text{(1)}& Q_{\text{s},i}={\displaystyle \frac{\mathrm{1}}{N}}{\displaystyle }\sum _j^{j\in N}{Q}_i^j.\end{array}$$

Here Q_s,i (m³ s⁻¹) is the multi-year averaged monthly streamflow of month i; j is the year index; and N is the number of years accounted for. For an upcoming year j, we only use the historical simulations (maximum latest 10 years) to calculate Q_s.

Secondly, we evaluate the target water storage change ΔW_t and monthly streamflow Q_t considering the regulation capacity of each reservoir. As shown in Fig. S1 in the Supplement, 1 year can be divided into two periods by comparing Q_s with ${\stackrel{\mathrm{‾}}{Q}}_{\text{s}}$. The longest continuous period of months with $Q_{\text{s}}>{\stackrel{\mathrm{‾}}{Q}}_{\text{s}}$ is the recharging season for reservoirs, and the rest is the releasing season. The amount of water stored during the recharging season (blue region in Fig. S1) is determined by C_max and is used during the releasing season (red regions in Fig. S1). The values of ΔW_t and Q_t can be estimated by

$$\begin{array}{}\text{(2)}& {\displaystyle }& {\displaystyle }k=min({\displaystyle \frac{C_{\text{max}}}{\mathit{\alpha }{\displaystyle }\sum _i^{i\in \text{Recharge}}{Q}_{\text{s},i}}},k_{\text{max}}),\text{(3)}& {\displaystyle }& {\displaystyle }\mathrm{\Delta }{W}_{\text{t},i}=\mathit{\alpha }\left[k\left(Q_{\text{s},i}-{\stackrel{\mathrm{‾}}{Q}}_{\text{s}}\right)+{\stackrel{\mathrm{‾}}{Q}}_{\text{s}}\right],\text{(4)}& {\displaystyle }& {\displaystyle }{Q}_{\text{t},i}=Q_{\text{s},i}-\mathrm{\Delta }{W}_{\text{t},i}/\mathit{\alpha }.\end{array}$$

Here k (–), varying between 0 and k_max (=0.7), indicates the ability of the reservoir to disturb streamflow seasonality. It is a ratio of the maximum regulation capacity of the reservoir C_max (10⁸ m³) over the streamflow amount throughout the recharging season. α (0.0263) converts monthly streamflow to water volume. Assuming that the water storage of the reservoir reaches C_max at the end of the recharging season, we can calculate target water storage W_t by using ΔW_t.

Finally, the variation of the actual water storage of the reservoir ΔW is a decision regarding actual monthly streamflow, current water storage, Q_t, ΔW_t, and W_t. During the releasing season, ΔW is calculated as

$$\begin{array}{}\text{(5)}& \begin{array}{ll}-W_i\frac{\left(-\mathrm{\Delta }{W}_{\text{t},i}\right)}{W_{\text{t},i}}& \text{if }{W}_i\le W_{\text{t},i};\\ \mathrm{\Delta }{\stackrel{\mathrm{̃}}{W}}_i-\left[\left(W_i+\mathrm{\Delta }{\stackrel{\mathrm{̃}}{W}}_i\right)-\left(W_{\text{t},i}+\mathrm{\Delta }{W}_{\text{t},i}\right)\right]& \text{if }{W}_i>W_{\text{t},i}\text{ and }\mathrm{\Delta }{\stackrel{\mathrm{̃}}{W}}_i>\mathrm{\Delta }{W}_{\text{t},i};\\ \mathrm{\Delta }{W}_{\text{t},i}-\left(W_i-W_{\text{t},i}\right)& \text{if }{W}_i>W_{\text{t},i}\text{ and }\mathrm{\Delta }{\stackrel{\mathrm{̃}}{W}}_i\le \mathrm{\Delta }{W}_{\text{t},i}.\end{array}\end{array}$$

$\mathrm{\Delta }{\stackrel{\mathrm{̃}}{W}}_i=\mathit{\alpha }{Q}_i-(\mathit{\alpha }{Q}_{\text{t},i}-\mathrm{\Delta }{W}_{\text{t},i})$. It is the expected release amount to make river streamflows equal to the target streamflows after reservoir regulation. If current water storage is less than the target value (the case of Eq. 5), the ΔW_i is calculated by the W_i with a proportion of ΔW_t,i over W_t,i. If the current storage is more than the target value (the cases of Eq. 5), the reservoir can release more water based on a balance between the target water storage change ΔW_t,i and the target water storage at the next time step W_t,i (represented by $\mathrm{\Delta }{\stackrel{\mathrm{̃}}{W}}_i$). Note that all water storage change variables are negative throughout the releasing season.

During the recharging season, we can calculate the ΔW_i as

$$\begin{array}{}\text{(6)}& \begin{array}{ll}max\left(min\left(W_{\text{t},i}+\mathrm{\Delta }{W}_{\text{t},i}-W_i,\mathit{\alpha }{Q}_i\right),\mathrm{0}\right)& \text{if }{W}_i>W_{\text{t},i};\\ min\left(\mathrm{\Delta }{W}_{\text{t},i}+\left(W_{\text{t},i}-W_i\right),\mathit{\alpha }{Q}_i\right)& \text{if }{W}_i\le W_{\text{t},i}.\end{array}\end{array}$$

If current water storage is larger than the target value (Eq. 6), we will try to recharge a volume of water to make $W_{i+\mathrm{1}}=W_{\text{t},i+\mathrm{1}}$. If current water storage is less than the target value (Eq. 6), we decide to recharge additional water volume besides the ΔW_t,i.

ΔW is then applied as a correction of simulated streamflows to generate actual monthly streamflows using the following equation:

$$\begin{array}{}\text{(7)}& {\widehat{Q}}_{\text{sim},i}=Q_{\text{sim},i}-{\displaystyle \frac{\mathrm{1}}{\mathit{\alpha }}}\mathrm{\Delta }{W}_i.\end{array}$$

Here ${\widehat{Q}}_{\text{sim}}$ (m³ s⁻¹) is the simulated regulated streamflows, and Q_sim (m³ s⁻¹) is the simulated monthly streamflows. Note that this model is a simplified representation of dam management because it ignores the direct coupling between water demand and irrigation water supply from the cascade of upstream reservoirs. This approach implies that, with a regulated flow, demands will be able to be satisfied and floods avoided without being more explicit. A complete coupling of demand, flood, and reservoir management is difficult to implement in the land surface model in the absence of data about the purpose and management strategy of each dam, given different possibly conflicting demands of water for industry and drinking versus cropland irrigation.

Before performing the simulation, we estimate the maximum regulation capacity of each studied reservoir in each river sub-catchment shown in Fig. 1. Table 1 lists collected information of the main reservoirs on the YR. Only large reservoirs like LongYangXia (LYX), LiuJiaXia (LJX), and XiaoLangDi (XLD) are considered in our model because of their huge C_max.

Table 1 Information of artificial reservoirs on the YR with considerable total capacity. Data are mainly from the YR Conservancy Commission of the Ministry of Water Resources (http://www.yrcc.gov.cn, last access: 28 February 2021). The “regulation purposes” follow the style of Hanasaki et al. (2006). “H” indicates hydropower; “C” indicates flood control; “I” indicates irrigation; “W” indicates water supply; and “S” indicates scouring sediment.

^* The total capacity shrink is due to sedimentation.

2.2 Sub-catchment diagnosis

Figure 1 shows the YRB and main gauging stations used in this study. To effectively use Q_obs for investigating impacts of irrigation and dam regulations on the streamflows of different river sub-catchments, we divided the YRB into five sub-catchments (R_i, $i\in [\mathrm{1},\mathrm{5}]$; Fig. 1) with an outlet at each gauging station. Thus we can evaluate the water balance in R_i by

$$\begin{array}{}\text{(8)}& {\displaystyle \frac{\mathrm{\Delta }{\text{TWS}}_i}{\mathrm{\Delta }t}}=P_i-{\text{ET}}_i+{\displaystyle \frac{Q_{\text{in},i}-Q_{\text{out},i}}{A_i}},\end{array}$$

where Δt is the time interval, ΔTWS_i (mm) is the change of total water storage in specific R_i, P_i (mm Δt⁻¹) is precipitation in R_i, ET_i (mm Δt⁻¹) is evapotranspiration in R_i, and A_i (m²) is the area of R_i. Q_in,i and Q_out,i (m³ Δt⁻¹) are inflow and outflow respectively. In addition, $q_i=Q_{\text{out},i}-Q_{\text{in},i}$ indicates the contribution of R_i to the river streamflows, that is the sub-catchment streamflows. This term can be negative if local water supply (e.g., precipitation and groundwater) cannot meet water demand. A conceptual figure of the water balance of a sub-catchment is shown in the top left of Fig. 1.

2.3 Evaluation datasets

Observed monthly streamflows (Q_obs) from the gauging stations shown in Fig. 1 are used to evaluate the simulations. Several precipitation (P) and evapotranspiration (ET) datasets were selected to evaluate the simulated water budgets in each sub-catchment of the YRB. The 0.5^∘ 3-hourly precipitation data from GSWP3 (Global Soil Wetness Project Phase 3) used as model input are based on GPCC v6 (Global Precipitation Climatology Centre; Becker et al., 2013) after bias correction with observations. The MSWEP (Multi-Source Weighted-Ensemble Precipitation) is a 0.25^∘ 3-hourly P product integrating numerous in situ measurements, satellite observations, and meteorological reanalysis (Beck et al., 2017). Three ET datasets are chosen for their potential ability to capture the effect of irrigation disturbance on ET (Yin et al., 2020) (noted as ET_obs). GLEAM v3.2a (Global Land Evaporation Amsterdam Model; Martens et al., 2017) provides 0.25^∘ daily ET estimations based on a two-soil layer model in which the top soil moisture is constrained by the ESA CCI (European Space Agency Climate Change Initiative) Soil Moisture observations. The FLUXCOM model (Jung et al., 2009) upscales ET data from a global network of eddy covariance tower measurements into a global 0.5^∘ monthly ET product. Since these towers do not cover irrigated systems, ET from irrigation simulated by the LPJmL (Lund–Postam–Jena managed Land) is added to ET from non-irrigated systems. The PKU ET product estimates 0.5^∘ monthly ET using water balances at basin scale, integrating FLUXNET observations to diagnose sub-basin patterns using a model tree ensemble approach (Zeng et al., 2014).

2.4 Simulation protocol

The 0.5^∘ half-hourly GSWP3 atmospheric forcing (Kim, 2017) was used to drive ORCHIDEE simulations. Yin et al. (2018) used four atmospheric forcing datasets to drive ORCHIDEE to simulate soil moisture dynamics over China, and they found that the GSWP3 provided the best performances; hence we chose this forcing for this study. A 0.5^∘ map with 15 different plant functional types (PFTs) containing crop sowing area information for the three PFTs corresponding to the modeled crop (wheat, maize, and rice) is used, based on a 1:1 million vegetation map and provincial-scale census data of China. Crop planting dates for wheat, maize, and rice are derived from the spatial interpolation of phenological observations from the Chinese Meteorological Administration (Wang et al., 2017). The soil texture map is from Zobler (1986). Two simulation experiments were performed to assess the impacts of irrigation on streamflows: (1) NI, no irrigation, and (2) IR, irrigated by available water resources. In IR, only surface irrigation is considered, that is, water applied on the cropland soil without interception by canopies. The soil water stress, a function of soil moisture and crop root density up to 2 m depth (Yin et al., 2020), is checked every half an hour. When it is less than a target threshold, irrigation is triggered with an amount equal to the difference between saturated and current soil moisture. To precisely estimate irrigation water consumption (direct water loss from the surface water pool excluding return flow), deep drainage of the three crop soil columns is turned off in the IR simulation. Simulations start from a 20-year spin-up in 1982 to initialize the thermal and hydrological variables then continued from 1982 to 2014. A validation against naturalized streamflows is shown in the Supplement Table S1.

Table 2 Definitions of sub-catchments and values of C_dam used in the dam regulation simulation.

The dam operation simulation starts from 1982 as an offline model applied to the simulated streamflows from the IR simulation (Q_IR) as input. The initial values of W were set to half of the C_max. Considering the potential joint regulation of reservoirs, we firstly estimate the total ΔW of all considered reservoirs by using Q_IR at HuaYuanKou (outlet of R₄; Fig. 1). Then we estimate the ΔW of LYX and LJX reservoir by using Q_IR at LanZhou. The difference between these two ΔW is assumed to be the ΔW of the XLD reservoir in between. Offline simulated ΔW values are used to estimate regulated monthly streamflows (${\widehat{Q}}_{\text{IR}}$) as in Eq. (7). As huge irrigation water withdrawals occur in R₃ and R₅ (YRCC, 2014), the water recharge of reservoirs may result in negative ${\widehat{Q}}_{\text{IR}}$ at TouDaoGuai and LiJin. To avoid this numerical artifact due to the offline nature of our dam model, we corrected all negative ${\widehat{Q}}_{\text{IR}}$ to zero by assuming that the streamflows cannot further drop when all stream water is consumed upstream. The impact of these corrections is accounted for at gauging stations downstream to ensure mass conversation.

2.5 Evaluation metrics

Three metrics are used to evaluate the performances of simulated monthly Q. The mean square error (MSE) evaluates the magnitude of errors between simulation and observations. It can be decomposed into three components (Kobayashi and Salam, 2000):

$$\begin{array}{}\text{(9)}& \text{MSE}={\displaystyle \frac{\mathrm{1}}{n}}\sum _{i=\mathrm{1}}^n{\left(S_i-O_i\right)}^{\mathrm{2}}=\text{SB}+\text{SDSD}+\text{LCS},\end{array}$$

where S_i and O_i are simulated and observed values, respectively, and n is the number of samples. SB (squared bias) is the bias between simulated and observed values. In this study, SB represents the difference between simulated and observed multi-year mean annual Q. SDSD (the squared difference between standard deviation) relates to the mismatch of variation amplitudes between simulated and measured values. It can reflect whether our simulation can capture the seasonality of Q_obs. LCS (the lack of correlation weighted by the standard deviation) indicates the mismatch of fluctuation patterns between simulated and observed values, which is equivalent to inter-annual variation of Q in this study. The formulas of these three components and a detailed explanation can be found in Kobayashi and Salam (2000).

The index of agreement ($d\in [\mathrm{0},\mathrm{1}]$) is defined as the ratio of MSE and potential error. It is calculated as

$$\begin{array}{}\text{(10)}& d=\mathrm{1}-{\displaystyle \frac{{\sum }_{i=\mathrm{1}}^n{\left(O_i-S_i\right)}^{\mathrm{2}}}{{\sum }_{i=\mathrm{1}}^n{\left(|S_i-\stackrel{\mathrm{‾}}{O}|+|O_i-\stackrel{\mathrm{‾}}{O}|\right)}^{\mathrm{2}}}},\end{array}$$

where d=1 indicates perfect fit, and d=0 denotes poor agreement.

The modified Kling–Gupta efficiency ($\text{mKGE}\in (-\mathrm{\infty },\mathrm{1}]$) is defined as the Euclidean distance of three independent criteria: correlation coefficient r, bias ratio β, and variability ratio γ (Gupta et al., 2009; Kling et al., 2012). It is an improved indicator from the Nash–Sutcliffe efficiency avoiding heterogeneous sensitivities to peak and low flows, which is crucial for this study that is not only interested in simulating peak flows but also concentrates on base flows regulated by dams for human usage. The mKGE is calculated as

$$\begin{array}{}\text{(11)}& {\displaystyle }\text{mKGE}=\mathrm{1}-\sqrt{(\mathrm{1}-r{)}^{\mathrm{2}}+(\mathrm{1}-\mathit{\beta }{)}^{\mathrm{2}}+(\mathrm{1}-\mathit{\gamma }{)}^{\mathrm{2}}},\text{(12)}& {\displaystyle }\mathit{\beta }={\displaystyle \frac{{\mathit{\mu }}_S}{{\mathit{\mu }}_O}};\mathit{\gamma }={\displaystyle \frac{{\text{CV}}_S}{{\text{CV}}_O}},\end{array}$$

where r is the correlation coefficient between observed and simulated streamflows, and μ (m³ s⁻¹) and CV (–) are the mean and the coefficient of variation of Q, respectively. These indicators are used for three comparisons: (1) Q_NI and Q_obs, (2) Q_IR and Q_obs, and (3) ${\widehat{Q}}_{\text{IR}}$ and Q_obs.

Figure 2 Top panels: annual mean of hydrological elements in each sub-catchment of the YR basin from both simulation (plain colors) and observation (hatched colors). Error bars represent standard deviation. Bottom panels: trends of these elements in each sub-catchment. Dark borders indicate the trend is statistically significant (p-value <0.05) according to the Mann–Kendall test.

3 Results

3.1 Water budgets at sub-catchment scale

Figure 2 displays water budgets and trends in R_i based on simulation and observations. Going from upstream to downstream, precipitation in P_GSWP3, which is consistent with P_MSWEP, decreases from 543.6 mm yr⁻¹ (R₁) to 254.2 mm yr⁻¹ (R₃) and then rises again to 652.1 mm yr⁻¹ (R₅). The magnitudes of simulated ET (both ET_NI and ET_IR) have no significant differences with ET_obs aggregated over sub-catchments R₁ to R₅. Grid-cell-based validation shows high agreement between simulated and observed ET across all sub-catchments. The lowest mean of correlation coefficients is 0.79, and the highest mean of relative root mean square error (RMSE) is 4.9 % (Table S2). Except for R₁ where cropland is rare, ET_IR accounts for an amount representing more than 80 % of P_GSWP3 in the YRB, with a maximum value of 96.5 % in R₃. The difference between ET_IR and ET_NI is due to the irrigation process, which accounts for 9.1 % and 8.2 % of ET_NI in R₃ and R₅ respectively as caused by the irrigation demand. The impact of irrigation can be detected from sub-catchment streamflows ($q_i=(Q_{\text{out},i}-Q_{\text{in},i})/A_i$) as well. For instance, both q_obs and q_IR are negative in R₃ and R₅, suggesting that local surface water resources cannot meet water demand for irrigation. As irrigation water transfers between grid cells are not represented in our simulations, the non-availability of water locally results in an underestimation of the irrigation withdrawals, likely explaining why q_IR>q_obs in R₃ to R₅.

The trends of P and ET are positive but not significant in most R_i during the period 1982–2014 (bottom panels of Fig. 2). However, significant trends can be found in simulated and observed q in some R_i. The decrease of q_obs in R₁ is not captured by the model, neither in q_NI nor q_IR. This underestimated decrease of river streamflows might be linked to decreased glacier melt or increased non-irrigation-related human water withdrawals, which are ignored in our simulations. In R₂ and R₃, the q_obs trends are determined by the joint effects of climate change (e.g., the P increase) and human water withdrawals. The trends of q_IR show the same direction as those of q_obs. In R₅, however, q_obs increased by 1.67 mm yr⁻¹, which was not captured by our simulation of q_IR. Besides the increase of P, another possible driver of increasing q_obs in R₅ is a decrease of water withdrawal due to the improvement of irrigation efficiency (Yin et al., 2020), which is not accounted for in our simulations. Moreover, water use management may play an important role in the observed positive trends of q_obs as well, with the aim to increase the streamflows downstream of the YR to avoid streamflow cutoff (Q_obs<1 m³ s⁻¹) that occurred in the 1990s (Wang et al., 2006).

Irrigation not only influences annual streamflows in the YR but also affects its intra-annual variation. In general, the discharge yield Y_Q, defined by the sum of surface runoff and drainage, of all grid cells in NI should be higher than in IR because our irrigation model can remove water from the stream reservoirs which is a fraction of drainage and runoff. However, our simulations show that Y_Q,NI can be less than Y_Q,IR (Fig. S2) at the beginning of the monsoon season. This is because irrigation keeps soil moisture higher than without irrigation in July in R₄ and R₅ (Fig. S2d and e), which in turn promotes Y_Q because the soil water-holding capacity is lower, and a larger fraction of P can go to runoff. This mechanism highlighted that irrigation could enhance the heterogeneity of water temporal distribution and may reinforce floods after a dry season.

Figure 3 (a–e) Time series of annual streamflows from observations and simulations at each gauging station. (f–j) Seasonality of observed and simulated streamflows at each gauging station. Q_obs is the observed annual mean streamflows. Q_NI and Q_IR are the simulated annual mean streamflows based on the NI and IR simulations (Sect. 2.4), respectively. These simulations do not account for dams, and therefore the seasonality has a higher amplitude than observed in the plots on the right-hand side.

3.2 Comparison between observed and simulated Q

Figure 3 shows time series of annual streamflows and of the seasonality of monthly streamflows. Our simulations underestimate Q_obs at TangNaiHai in R₁, likely because we miss glacier melt. After LanZhou, the values of Q_IR coincide very well with those of Q_obs, indicating that irrigation strongly reduces the annual streamflows of the YR. However, the seasonality of monthly Q_IR is different from Q_obs (Fig. 3f–j). Despite the good match of annual values, the model without dams (shown in Fig. 3) produces an underestimation of Q in the dry season and an overestimation of Q in the flood season. Such a mismatch of Q seasonality is likely caused primarily by dam regulation ignored in the model. The locations of several big reservoirs are shown in Fig. 1b, and their characteristics are listed in Table 1.

Figure 4 (a–b) Monthly observed (Q_obs) and simulated (Q_NI) streamflows at TangNaiHai and LanZhou stations. Green bar in (b) indicates the start of the LongYangXia dam regulation. (c) Observed annual streamflows at TongGuan and XiaoLangDi gauging stations, which are located up- and downstream of the XiaoLangDi reservoir, respectively (see Fig. 1). Blue bar in (c) indicates the start of the XiaoLangDi dam regulation.

Regarding dams, before the operation of the LongYangXia dam in 1986, which brought a regulation capacity of 193.5×10⁸ m³ (green bar in Fig. 4b), the peaks of monthly Q_NI at LanZhou were slightly lower than the peaks of Q_obs in R₂ (Fig. 4b), as well as at TangNaiHai (Fig. 4a). But after the construction of the LongYangXia dam, modeled peak Q_NI became systematically higher than the peak of Q_obs each year, suggesting that the construction of this dam caused the observed peak reduction (Fig. 4b). Moreover, the seasonality of Q_obs changed dramatically in the period (1982–2014), but no similar trend was found in monthly P (Fig. S3), suggesting that dam operation was the primary driver of the observed shift in seasonal streamflow variations of the YRB from 1982 to 2014. Dams can affect inter-annual variations of Q as well, although less than the seasonal variation. For instance, TongGuan and XiaoLangDi are two consecutive gauging stations upstream and downstream of the reservoir of XiaoLangDi in R₄ (Fig. 1). The annual Q_obs at the two stations shows different features after the construction of the XiaoLangDi reservoir in 1999.

Figure 5 Comparison between observed and simulated actual monthly streamflows at gauging stations. Q_obs (dark lines) is observed monthly streamflows. Q_IR (blue lines) is simulated monthly streamflows from the IR experiment (Sect. 2.4). ${\widehat{Q}}_{\text{IR}}$ (orange lines) is simulated monthly streamflows including impacts of dam regulation (Sect. 2.4).

Figure 5 shows monthly time series of Q_obs, Q_IR, and ${\widehat{Q}}_{\text{IR}}$ (see Sect. 2.1.2) at each gauging station. Discharge fluctuations are successfully improved in ${\widehat{Q}}_{\text{IR}}$. Especially the base flow of ${\widehat{Q}}_{\text{IR}}$ coincides well with that of Q_obs during winter and spring. The only exception occurs at LiJin, where ${\widehat{Q}}_{\text{IR}}$ overestimates the streamflows from January to May. In fact, the water release from XLD during this period would be withdrawn for irrigation and industry in R₅. However, our offline dam model is not able to simulate the interactions, leading to the overestimation.

The dam model is successful in reproducing flood control as well. At LanZhou, although ${\widehat{Q}}_{\text{IR}}$ underestimates the peak flow due to the bias of the simulated mean annual streamflows (Fig. 3b), its seasonality is much smoother than that of Q_IR. The underestimation of ${\widehat{Q}}_{\text{IR}}$ can reflect special water management during extreme years. From 2000 to 2002, the YRB experienced severe droughts, with 10 %–15 % precipitation less than usual, leading to a decrease of surface water resources by as much as 45 % (Water Resources Bulletin of China; http://www.mwr.gov.cn/sj/tjgb/szygb/, last access: 28 February 2021). To guarantee base flow, a set of policies were applied (e.g., reducing water withdrawn, increasing water price, and releasing more water from reservoirs). These policies are not accounted for in the model, which will produce a higher irrigation demand during dry years and promote the underestimation of the Q_obs. From TouDaoGuai to LiJin, the floods from August to October are dramatically reduced by our dam model. Nevertheless, the peaks are still overestimated in ${\widehat{Q}}_{\text{IR}}$, which might be due to numerous non-modeled medium/small reservoirs that were ignored by our model; no fewer than 203 medium reservoirs were documented at the end of 2014 (YRCC, 2014). In our simulation, a 326.5 ×10⁸ m³ regulation capacity is considered, which only accounts for 45 % of the total storage capacity of 720 ×10⁸ m³ (Ran and Lu, 2012). Moreover, in the five irrigation districts (http://www.yrcc.gov.cn/hhyl/yhgq/, last access: 28 February 2021; Tang et al., 2008), special irrigation systems in the YRB could contribute to flood reduction. For instance, the Hetao Plateau is the traditional irrigation district and is equipped with a hydraulic system that can divert river water into a complex irrigation network (bounded by 106.5–109^∘ E and 40.5–41.5^∘ N, shown in Fig. 1), by adjusting water level differences during the flood season. This no-dam diversion system of the Hetao Plateau can take 50 ×10⁸ m³ as an extra regulation capacity per year, equivalent to 14 % of the annual streamflow in R₃.

Simulated ΔW in R₂ is compared to observations (Jin et al., 2017) in the left panel of Fig. 6, and the good agreement suggests that our dam model is able to capture the seasonal variation of ΔW (r=0.9, p<0.001) and rectify the simulated streamflows. In the case of XiaoLangDi in the right panel of Fig. 6, where the correlation is smaller (r=0.34, p=0.28), the mismatch could be explained by sediment regulation procedures of that dam, given that it releases a huge amount of water in June for reservoir cleaning and sediment flushing downstream (Baoligao et al., 2016; Kong et al., 2017; Zhuo et al., 2019), a process not represented in our simple dam model. Moreover, because we ignored the buffering effect of numerous medium reservoirs, the simulated water recharge during the flood season could be overestimated.

Figure 6 The changes of water storage of dams (ΔW) in R₂ and R₄. The dark line represents the ΔW from literature. The multi-year mean of ΔW of LongYangXia and LiuJiaXia is from Jin et al. (2017). The ΔW of XiaoLangDi is from the 1-year record reported in Kong et al. (2017). Red lines represent corresponding simulated ΔW from our dam regulation model.

Figure 7 Indicators of Q comparisons in each sub-catchment of the YRB. Colors indicate different comparisons. The MSE is decomposed to SB, SDSD, and LCS, which are distinguished by different transparencies.

Figure 7 presents the model performances with different metrics in different R_i. The results show that MSE increases considerably from R₁ to R₅, implying accumulated impacts of error sources in increasing the error of modeled Q when going downstream in the entire catchment. Most likely, those error sources are omission errors of anthropogenic factors such as drinking and industrial water removals, but also of natural factors such as riparian wetlands and floodplains (e.g., the SanShengGong water conservancy hub) and non-represented small streams in the routing of ORCHIDEE (e.g., the irrigation system at the Hetao Plateau). From the decomposition of MSE, we found that adding irrigation to the model removes most of the bias in the average magnitude Q by reducing the SB bias error term of the MSE. The only exception occurs at LanZhou, where SB increases in IR, consequently leading to higher MSE. This misfit is due to the underestimation of Q upstream (Fig. 3a). Thus, modeled Q_IR is lower at LanZhou, which enlarges the SB. On the other hand, adding the dam operation contributes to improve the phase variations of Q which are dominated by the phase of the seasonal cycle, by reducing the SDSD error term. Nevertheless, the LCS error term, indicating the magnitude of the variability, mainly at inter-annual timescales, has no significant improvement with the representation of irrigation and dam regulations. It is because some of reservoirs are able to regulate Q inter-annually (Table 1), which can be observed from Fig. 4c. However, related operation rules are unclear and are not implemented in our dam model. Improvements were found in d as well, which demonstrates that the way human effects on Q of the YR were modeled brings more realistic results, despite ignoring the direct effect of irrigation demand on reservoir release and ignoring industrial and domestic water demands. The mKGE reveals a significant increase after considering dam operations (Fig. 7c). Particularly at LanZhou and HuaYuanKou, the mKGE of ${\widehat{Q}}_{\text{IR}}\sim Q_{\text{obs}}$ increases by 0.86 and 1.11 compared to Q_IR∼Q_obs, respectively. Note that the mKGEs of Q_IR∼Q_obs are smaller than those of Q_NI∼Q_obs from R₂ to R₄ because irrigation decreases the mean annual streamflow of Q_IR, which further increases the CV_S, leading to worse γ in mKGE (Eq. 11).

4 Discussion

This study shows that ORCHIDEE land surface model with crops, irrigation, and our simple dam operation model can reproduce streamflow mean levels, inter-annual variations, and the seasonal cycle in different sub-catchments of the YRB correctly. We preliminarily quantified the impacts of irrigation and dams on the fluctuations of streamflows. Simulated water balance components were compared to observations in different sub-catchments with a good agreement (e.g., 4.5±6.9% for ET). We found that irrigation mainly affects the magnitude of annual streamflows by consuming $\mathrm{242.8}\pm \mathrm{27.8}\times {\mathrm{10}}^{\mathrm{8}}$ m³ yr⁻¹ of water, consistent with census data giving a consumption of $\mathrm{231.4}\pm \mathrm{31.6}\times {\mathrm{10}}^{\mathrm{8}}$ m³ yr⁻¹ (YRCC, 2014). As the water of the YRB is reaching the limit of usage (Feng et al., 2016), we did not find any significant effect of irrigation on streamflow trends. Instead of increasing river water withdrawals, the growing water demand appeared to have been balanced by improving water use efficiency during the study period (Yin et al., 2020; Zhou et al., 2020). Our simulation reveals that the impact of irrigation on streamflows may even be positive under special situations, which was also shown in Kustu et al. (2011). However, our mechanisms are different from the irrigation–ET–precipitation atmospheric feedback mechanisms found by Kustu et al. (2011); we demonstrated that irrigation may significantly increase soil moisture and promote runoff yield during the following wet season. It implies that irrigation in such landscapes may reinforce the magnitude of floods during the rainy season by a higher legacy soil moisture.

We found that dams strongly regulate the temporal variation of streamflows (Chen et al., 2016; Li et al., 2016; Yaghmaei et al., 2018). By including simple regulation rules depending on inflows, our dam model reduced the simulation error by 48 %–77 % (MSE in Fig. 7), especially the variability component (SDSD) of the total error, which is dominated by seasonal misfit reduction from dams. Moreover, we confirmed that the change of Q_obs seasonality during the study period is not due to climate change (Fig. S3) but is determined by dam operations (Wang et al., 2006). Big dams, like the LongYangXia, LiuJiaXia, and XiaoLangDi, are able to regulate streamflows inter-annually (Wang et al., 2018) and smooth the inter-annual distribution of water resources in the YRB (Piao et al., 2010; Wang et al., 2006; YRCC, 2014). However, their detailed operation rules are unclear and were not implemented explicitly in our dam model. The error corresponding to inter-annual variation (LCS in MSE in Fig. 7) was thus not reduced by including dams. In the dam model, some functions of reservoirs, such as providing irrigation supply, industrial and domestic water, electricity generation, and flood control (Basheer and Elagib, 2018), are not explicitly represented. Particularly the XiaoLangDi dam carries a distinctive water-sediment mission, which scours sediments downstream by creating artificial floods in June (Kong et al., 2017; Zhuo et al., 2019). These functions are associated with many socioeconomic factors and drivers, leading to competing demands for water (e.g., policies, electricity price, water price, land use change, irrigation techniques, water management techniques, and dams inter-connection), which are difficult to model well due to a lack of data. However, with the upcoming Surface Water and Ocean Topography (SWOT) mission, it will be possible to monitor the water level and surface extent of more reservoirs, which will be helpful to improve and validate the dam operation simulations (Ottlé et al., 2020).

Our simulations ignored potential impacts of reservoirs on local climate, e.g., through their evaporation (Degu et al., 2011). The water area of several artificial reservoirs (LongYangXia, LiuJiaXia, BoHaiWan, SanShengGong, and XiaoLangDi) is approximately 1056 km², which is larger than the 10th largest natural lake in China. These water bodies can also significantly influence local energy budgets, and evaporative water loss may be considerable, especially in arid and semi-arid areas (Friedrich et al., 2018; Shiklomanov, 1999). In addition, the five large irrigation districts (http://www.yrcc.gov.cn/hhyl/yhgq/, last access: 27 February 2021) could dramatically alter the local climate through atmospheric feedbacks. For instance, the Hetao Plateau can take about 50 ×10⁸ m³ from streamflows every year during the flood season. Its irrigation area is 5740 km², with an evapotranspiration rate ranging between 1200–1600 mm yr⁻¹. However, as these irrigation districts divert river water without big dams, they are not taken into account in most YR studies. Another non-negligible factor in the case of YR is sedimentation, which reduces the regulation capacities of reservoirs and weakens streamflow regulation by humans. For instance, the total capacity of QingTongXia has declined from 6.06 to 0.4×10⁸ m³ since 1978 due to sedimentation. Therefore, how land use change and the evolution of natural ecosystems affect sediment load and deposition is another key factor to project dams' disturbances on streamflows in the YRB.

Simulating anthropogenic impact on river streamflows is challenging. In the case of the YR, well-calibrated models can provide accurate naturalized streamflow simulations, with Nash–Sutcliffe efficiency (NSE) as high as 0.9 (Yuan et al., 2016). However, when considering the impacts of irrigation and dams, the NSE values of simulations are generally worse. For instance, simulations with anthropogenic effects by Hanasaki et al. (2018) had lower NSE than the simulation with only natural processes. Similarly, Wada et al. (2014) showed NSE decrease after considering anthropogenic factors in the YRB, which was interpreted as complexity of the YRB under the impacts of human activities and climate variation. However, the NSE of naturalized streamflows cannot really be compared to the one of regulated streamflows. Even if a model can perfectly simulate the dam operations, the NSE of naturalized streamflows will be larger than that of regulated streamflows, given that dam operations automatically reduce the variation of river streamflows (a simple proof is available in Sect. A of the Supplement). In fact, our model performances are very similar to GHM simulations (Fig. S2 from Liu et al., 2019). By gradually considering anthropogenic factors (irrigation and dam operations), the performances of our simulations increase dramatically according to all three metrics.

Intensive calibrations using a suite of observations can allow catchment-scale studies to provide highly accurate simulated streamflows for short-term flood forecasts. However, for long-term projections, a model should include all key processes of the system studied. If key processes are missing in the model, a calibration will cover up the shortcomings, which lead to a lack of predictive capacities for long timescales, as shown by Duethmann et al. (2020). Therefore, by developing crop physiology and phenology, irrigation, and a (offline) dam operation model, we have tried to demonstrate that streamflow fluctuations of the YR can be reasonably reproduced by a generic land surface model. Although mismatches exist in the simulated streamflows, they are more likely caused by missing processes (joint impact of multiple medium reservoirs, special mission of dams, irrigation system characteristics) than by poor calibration of existing processes because other simulated hydrological variables coincide well with observations in the YRB, such as soil moisture dynamics (Yin et al., 2018), naturalized river streamflows (Table S1 in Xi et al., 2018), leaf area index (Sect. S2 in Xi et al., 2018), amount and trend of irrigation withdrawals (Yin et al., 2020), trends of total water storage (Sect. 3.4 in Yin et al., 2020), and ET (Table S2). On the contrary, these mismatches draw our attention to some key mechanisms overlooked in most models. For instance, our model underestimates the annual streamflow at LanZhou in the period 2000–2002 (Fig. 3b), during which ${\widehat{Q}}_{\text{IR}}$ was almost negatively correlated to the Q_obs (Fig. 5a). In summary, our results show that the errors of simulated streamflows decreased dramatically after considering crops, irrigation, and dam operations, suggesting that these are first-order mechanisms controlling streamflow fluctuations. Future work can be focused on completing the model by linking dam operation to the variable crop water demand.

5 Conclusions

A land surface model ORCHIDEE and a newly developed dam model are utilized to simulate the streamflow fluctuations and dam operations in the Yellow River basin. The impacts of irrigation and dam regulation on streamflow fluctuations of the Yellow River were preliminarily qualified and quantified in this study by using a process-based land surface model and a dam operation model. Irrigation mainly contributes to the reduction of annual streamflow by as much as $\mathrm{242.8}\pm \mathrm{27.8}\times {\mathrm{10}}^{\mathrm{8}}$ m³ yr⁻¹. The shifts of intra-annual variation of the Yellow River streamflows appear not to be caused by climate change, at least not by significant changes of precipitation patterns and land use during the study period, but by the construction of dams and their operation. After considering the impacts of dams, we found that dam regulation can explain about 48 %–77 % of the fluctuations of streamflows. The effect of dams may be still underestimated because we only considered simple regulation rules based on inflows but ignored its interactions with irrigation demand downstream. Moreover, our analysis showed that several reservoirs on the Yellow River are able to influence streamflows inter-annually. However, such effects are not quantified due to a lack of knowledge of the regulation rules across our study period.