## Abstract

Double-reservoir-and-double-pumping-station systems are commonly used for irrigation water supply in hilly regions of southern China. An optimization model for this water supply system is proposed to minimize water shortage. The model features few coupling constraints, including available water in the system and pumping volume limited by regional water rights. Dynamic programming was adopted to solve the subsystem and aggregation models. The results with the model and that with the standard operation policy were compared; the total water shortage was reduced by 87.7%, total water replenishment from outside was reduced by 2.2%, and total water spill was reduced by 60.6% for a system in Nanjing, China. The method may provide a reference for optimal operation of water supply systems comprising reservoirs and pumping stations.

## INTRODUCTION

With the worldwide construction of reservoirs and inter-basin water transfer projects, an optimal joint operation of multi-reservoir systems has recently become a research hotspot. It has been recognized that benefits derived from the joint operation of a multi-reservoir system may exceed the sum of that from the independent operation of each reservoir in it (Wang *et al.* 2016).

Generally, in an optimization model of reservoir operation, the objective function is to minimize the sum of the squared deviation of actual water supply from the desired target, in an attempt to offer the same deviation at each period as far as possible (Celeste & Billib 2009). Such models are subject to several constraints, such as lower and upper bounds of release and storage, water balance equation, and non-negativity constraints. Furthermore, hydraulic connections between reservoirs should be added into constraints in a multi-reservoir system.

Optimization models include two types of decision variables (water supplies and spills) that increase with the number of reservoirs; increased variables lead to the ‘curse of dimensionality’ in the solving process (Chen *et al.* 2016). At present, different types of meta-heuristic algorithms featuring high applicability and high efficiency have been applied to solving such problems. The genetic algorithm was applied to optimizing the operation policy of a four-reservoir system (Jothiprakash *et al.* 2011). Particle swarm optimization was adopted for the optimal water resources management of a multi-reservoir system (Nabinejad *et al.* 2017). Gu *et al.* (2017) used a global search algorithm to optimize the joint operation of a system comprising multiple donor reservoirs and one recipient reservoir. Ehteram *et al.* (2018) suggested the Spider Monkey Algorithm for optimizing operation rules of multiple reservoirs. However, meta-heuristic algorithms may fail to solve models with equality constraints because of their fatal weakness in that they cannot guarantee globally optimal final results, which can be attributed to the random sampling (Birhanu *et al.* 2014; Rani & Srivastava 2015). In addition, the application of meta-heuristic algorithms to reservoir operation problems may lead to spills when the reservoir storage is below its capacity, which is against the normal operation rule.

Dynamic programming based on Bellman's principle (Kennedy 1986) can provide reliable and optimal results. In addition, it has good applicability in multi-stage decision-making processes, regardless of the matter whether the model is continuous or linear. Therefore, it has been widely used in solving optimization problems of reservoir operation. Although the ‘curse of dimensionality’ (Rani & Moreira 2009; Ahmad *et al.* 2014; Wang *et al.* 2016) is induced under the condition of a large number of reservoirs, appropriate dimension reduction methods, such as the decomposition-aggregation method (Gong & Cheng 2018) and decomposition-coordination method (Mahey *et al.* 2017), can be used to obtain global optimal results by transforming a high-dimensional problem into a series of low-dimensional problems.

The purpose of optimizing the joint operation of a reservoir and a pumping station is to minimize the operation cost, mainly the operation cost of the pumping station. In addition to water supplies and spills of the reservoir, water volume pumped by the pumping station during each period is also a type of decision variable subject to the maximum pumping capacity and pumping volume limited by water rights.

Accordingly, the difficulty of solving the model will increase with increasing numbers of reservoirs and pumping stations in the system. Yu *et al.* (1994) and Pulido-Calvo & Gutiérrez-Estrada (2011) developed a nonlinear optimization model for a reservoir and its replenishment pumping station to minimize the operation cost. Reca *et al.* (2014) constructed a nonlinear model for the joint operation of reservoirs and pumping stations to obtain the minimum operation cost and improved it considering the evaporation of the reservoir. Ðurin (2016) derived an optimal operating policy with a regression analysis of the running time mode of pumps, design discharge of pumping station, and storage capacity of the reservoir, according to the deterministic water supply scheme. However, all the above-mentioned works optimized the operation schedule of pumping stations according to seasonal electricity price or time-of-use price and contributed toward saving energy costs. In addition, they did not limit the total water volume pumped by the systems.

In different countries, concepts and connotations of water rights are not the same under different social backgrounds (Molle 2004; Heikkila 2015). However, it especially refers to the water resources usufruct here. Since 2016, the agricultural water price reform has been promoted stage by stage in China (Shen & Wu 2017). As an important premise, the water resources usufruct is defined in the form of the annual total volume of water which is distributed to users. In agriculture, the Chinese government has established an upper bound for the annual amount of irrigation water. Furthermore, the total amount of water that can be obtained from each irrigation gate or pumping station within a year is strictly limited (Sun *et al.* 2016).

This study aimed to develop an optimization method based on the decomposition-aggregation theory for the joint operation of a double-reservoir-and-double-pumping-station system that can provide optimal results following normal operation regulations. For this purpose, an optimization model of this system in Nanjing, China, was established as a case study, and the total available water in the system was coupled into it, considering regional water rights.

## DOUBLE-RESERVOIR-AND-DOUBLE-PUMPING-STATION SYSTEM

### Research area

Liuhe district of Nanjing, the Jiangsu province is located in eastern China. The hilly regions of this district lie in the subtropical monsoon climate zone. Therefore, local residents build reservoirs and water transfer projects to deal with the uneven temporal and spatial distribution of water resources. The typical irrigation water supply system in this region is the double-reservoir-and-double-pumping-station system shown in Figure 1, which consists of the SH reservoir, HWB reservoir, and their replenishment pumping stations.

Reservoirs and pumping stations

The SH reservoir is the second largest reservoir in Liuhe district and its main function is to supply water for irrigation. During the period of operation, the XZ pumping station replenishes it with water from the BL river through a 2 m × 2 m square concrete culvert before the water level reaches below the lower boundary limit. Because of the higher topography, the irrigation area of the HWB reservoir is larger although the catchment area is relatively small, and thus, the probability of water shortage is relatively high. During water shortage, the HZ pumping station replenishes the HWB reservoir with water from the SH reservoir through a 1.5 m × 1.5 m square concrete culvert. The characteristics of reservoirs and pumping stations are shown in Tables 1 and 2, respectively.

2. Inflows and water demand

Monthly inflows and water demand at the 75% probability of exceedance are shown in Tables 3 and 4, respectively.

3. Evaporation

*E*

_{601}evaporator (Table 5), and the data required correction by coefficient

*k*. The average free water surface of each reservoir in a specific period is determined by its surface–volume relationship, which was provided by the Liuhe Water Authority. Finally, the reservoir evaporation for a specific period can be computed using Equation (1): where

*EF*(10

_{i}^{4}m

^{3}) is the evaporation loss of a reservoir in period

*i*,

*E*(mm) is the evaporation depth of

_{i}*E*

_{601}evaporator in period

*i*;

*k*is the correction coefficient for period

_{i}*i*;

*V*(10

_{i}^{4}m

^{3}) is the average water storage in period

*i*; and

*α*and

*β*are the reservoir coefficients. For the SH reservoir,

*α*

*=*1.194 × 10

^{−3},

*β*

*=*2.575, and for the HWB reservoir,

*α*

*=*1.657 × 10

^{−3},

*β*

*=*1.862.

Reservoir | Dead storage capacity (10^{4} m^{3}) | Utilizable capacity (10^{4} m^{3}) | Total storage capacity (10^{4} m^{3}) | Limited storage capacity in flood season (10^{4} m^{3}) | Catchment area (km^{2}) | Irrigation area (hm^{2}) |
---|---|---|---|---|---|---|

SH | 600 | 1,157 | 2,473 | 1,000 | 31.9 | 1,667 |

HWB | 457 | 936 | 1,891 | 1,393 | 20.8 | 2,067 |

Reservoir | Dead storage capacity (10^{4} m^{3}) | Utilizable capacity (10^{4} m^{3}) | Total storage capacity (10^{4} m^{3}) | Limited storage capacity in flood season (10^{4} m^{3}) | Catchment area (km^{2}) | Irrigation area (hm^{2}) |
---|---|---|---|---|---|---|

SH | 600 | 1,157 | 2,473 | 1,000 | 31.9 | 1,667 |

HWB | 457 | 936 | 1,891 | 1,393 | 20.8 | 2,067 |

Pumping station | Design discharge (m^{3}/h) | Design pumping head (m) | Daily operation duration (h) | Water rights (75%) (10^{4} m^{3}) |
---|---|---|---|---|

XZ | 10,000 | 19.4 | 20 | 360 |

HZ | 7,500 | 15.0 | 20 | – |

Pumping station | Design discharge (m^{3}/h) | Design pumping head (m) | Daily operation duration (h) | Water rights (75%) (10^{4} m^{3}) |
---|---|---|---|---|

XZ | 10,000 | 19.4 | 20 | 360 |

HZ | 7,500 | 15.0 | 20 | – |

*Notes*: The annual water rights of the XZ pumping station at 75% probability of exceedance is allocated by the Liuhe Water Authority. The HZ pumping station is an internal pumping station in the system; there is no definite limit on its water rights.

Probability | Reservoir | Period | Total | |||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|

Oct. | Nov. | Dec. | Jan. | Feb. | Mar. | Apr. | May | Jun. | Jul. | Aug. | Sep. | |||

75% | SH | 115 | 45 | 11 | 32 | 2 | 28 | 32 | 50 | 445 | 781 | 618 | 320 | 2,479 |

HWB | 33 | 55 | 1 | 36 | 10 | 9 | 37 | 56 | 218 | 133 | 157 | 20 | 766 |

Probability | Reservoir | Period | Total | |||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|

Oct. | Nov. | Dec. | Jan. | Feb. | Mar. | Apr. | May | Jun. | Jul. | Aug. | Sep. | |||

75% | SH | 115 | 45 | 11 | 32 | 2 | 28 | 32 | 50 | 445 | 781 | 618 | 320 | 2,479 |

HWB | 33 | 55 | 1 | 36 | 10 | 9 | 37 | 56 | 218 | 133 | 157 | 20 | 766 |

Probability | Reservoir | Period | Total | |||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|

Oct. | Nov. | Dec. | Jan. | Feb. | Mar. | Apr. | May | Jun. | Jul. | Aug. | Sep. | |||

75% | SH | 23 | 15 | 18 | 17 | 23 | 26 | 30 | 31 | 523 | 323 | 159 | 33 | 1,221 |

HWB | 59 | 79 | 17 | 16 | 23 | 27 | 31 | 35 | 550 | 559 | 345 | 77 | 1,818 |

Probability | Reservoir | Period | Total | |||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|

Oct. | Nov. | Dec. | Jan. | Feb. | Mar. | Apr. | May | Jun. | Jul. | Aug. | Sep. | |||

75% | SH | 23 | 15 | 18 | 17 | 23 | 26 | 30 | 31 | 523 | 323 | 159 | 33 | 1,221 |

HWB | 59 | 79 | 17 | 16 | 23 | 27 | 31 | 35 | 550 | 559 | 345 | 77 | 1,818 |

Period | Oct. | Nov. | Dec. | Jan. | Feb. | Mar. | Apr. | May | Jun. | Jul. | Aug. | Sep. |
---|---|---|---|---|---|---|---|---|---|---|---|---|

E (mm) _{i} | 56 | 42 | 38 | 23 | 28 | 41 | 58 | 86 | 98 | 115 | 110 | 71 |

k _{i} | 1.04 | 1.12 | 1.12 | 1.05 | 0.92 | 0.9 | 0.88 | 0.92 | 0.94 | 0.94 | 0.98 | 1.06 |

Period | Oct. | Nov. | Dec. | Jan. | Feb. | Mar. | Apr. | May | Jun. | Jul. | Aug. | Sep. |
---|---|---|---|---|---|---|---|---|---|---|---|---|

E (mm) _{i} | 56 | 42 | 38 | 23 | 28 | 41 | 58 | 86 | 98 | 115 | 110 | 71 |

k _{i} | 1.04 | 1.12 | 1.12 | 1.05 | 0.92 | 0.9 | 0.88 | 0.92 | 0.94 | 0.94 | 0.98 | 1.06 |

*Notes*: The above data are provided by the Liuhe Water Authority.

### Generalization of the system

Figure 2 shows a generalized schematic of the double-reservoir-and-double-pumping-station system comprising the SH reservoir, HWB reservoir, and their replenishment pumping stations. The SH reservoir is a donor reservoir, numbered as reservoir 1; the HWB reservoir is a recipient reservoir, numbered as reservoir 2. The XZ pumping station, numbered as pumping station 1, diverts water from the BL river to the SH reservoir; the HZ pumping station, numbered as pumping station 2, diverts water from the SH reservoir to the HWB reservoir.

It has been proved that if there is no appropriate joint operation scheme for such a system, water replenishment, water spills, and water shortages may occur simultaneously. Therefore, an optimization method is necessary for such a system to improve the utilization efficiency of inflows as well as to reduce water spills and shortages, considering limited water rights on the river.

In Figure 2, *YS*_{1,i}, *LS*_{1,i}, *X*_{1,i}, and *PS*_{1,i} are water demand, inflow, water supply, and water spill in the SH reservoir in period *i*, respectively; *YS*_{2,i}, *LS*_{2,i}, *X*_{2,i}, and *PS*_{2,i} are water demand, inflow, water supply, and water spill in the HWB reservoir in period *i*, respectively; *Y*_{1,i}, and *Y*_{2,i} are water replenishment by the XZ and HZ pumping stations in period *i*, respectively.

## METHODOLOGY

### Optimization model

#### Objective function

*F*is the annual sum of squared water shortage of each reservoir in each period and

*i*is the period number (

*i*

*=*1, 2, …, 20).

#### Constraints

Annual available water in the system

*j*is the reservoir number (

*j*

*=*1,2);

*SK*

_{1}(10

^{4}m

^{3}) and

*SK*

_{2}(10

^{4}m

^{3}) are annual available water of the SH and HWB reservoirs, respectively; and

*BZ*

_{1}(10

^{4}m

^{3}) is the maximum annual pumping volume of the XZ pumping station.

- 2. Maximum annual pumping volume expressed as Equation (4): where
*BZ*_{2}(10^{4}m^{3}) is the maximum annual pumping volume of the HZ pumping station.

*V*(10

_{j,i}^{4}m

^{3}) is the water storage of reservoir

*j*in period

*i*;

*V*

_{j,i}_{(min)}and

*V*

_{j,i}_{,(max)}are the lower and upper bounds of water storage

5. Maximum pumping capacity

*Q*

_{1}(m

^{3}/h) and

*Q*

_{2}(m

^{3}/h) are the design discharges of the XZ and HZ pumping stations and

*N*(h) is the maximum operation duration in period

_{i}*i*.

6. Initial and boundary conditions

A restriction was imposed to make the final storage *V _{j}*

_{,20}equal to the initial storage

*V*

_{j}_{,0}. If this limitation is not imposed, the solution would tend to empty the reservoir in the final period.

7. Non-negative constraints

It is also necessary to introduce non-negativity constraints to avoid negative values to the variables, which would be physically impossible.

### Solving method

#### Solving schematic

The above model is a multi-dimensional nonlinear model, which can be divided into multiple stages. There are six types of decision variables, including *X*_{1,i}, *X*_{2,i}, *PS*_{1,i}, *PS*_{2,i}, *Y*_{1,i}, and *Y*_{2,i}, each of which has 20 dimensions. Therefore, it is difficult to solve it directly. As shown in Figure 3, the double-reservoir-and-double-pumping-station system can be decomposed into two subsystems that consist of one reservoir and one pumping station, and hydraulic connections between them are formed through the HZ pumping station. In this manner, the dimensionality of the problem can be reduced.

#### Solving steps

Decomposition of the system

The maximum annual pumping volume of the HZ pumping station, *BZ*_{2}, is taken as a coordinating variable, and the model of the double-reservoir-and-double-pumping-station system can be decomposed into two subsystem models that consist of one reservoir and one pumping station as Equations (10)–(13).

*f*

_{1}is the annual sum of squared water shortage in each period of subsystem 1;

*f*

_{2}is the annual sum of squared water shortage in each period of subsystem 2;

*W*

_{1}(10

^{4}m

^{3}) is the annual available water of subsystem 1; and

*W*

_{2}(10

^{4}m

^{3}) is the annual available water of subsystem 2.

The lower and upper bounds of water storage, water balance equation, maximum pumping capacity, maximum annual pumping volume, initial and boundary conditions, and non-negative constraints should also be imposed on the subsystems.

2. Optimization of subsystems

The subsystem models are two nonlinear models with only one coupling constraint. Therefore, one-dimensional dynamic programming (DP) can be used to solve them, taking water supply *X _{i}* as a decision variable. The normal operation rule is integrated into recursive procedures of DP (Shi

*et al.*2015; Gong

*et al.*2019) to correct water storage and obtain water spill

*PS*and water replenishment

_{i}*Y*of each period simultaneously; in this manner, the defect of meta-heuristic algorithms, which leads to random occurrences of spills or replenishments against the normal operation rule, can be avoided. The specific procedure is as follows.

_{i}If

*V*<_{i}*V*_{i}_{, min}, then the pumping station should replenish the reservoir and the final water storage should be*V*_{i}_{, min}. The water replenishment and spill in period*i*are expressed as Equations (14) and (15):If

*V*>_{i}*V*_{i}_{, max}, then excess water should be drained through the spillway. The water replenishment and spill in period*i*are expressed as Equations (16) and (17):- If
*V*_{i}_{, min}<*V*<_{i}*V*_{i}_{, max}, then the water replenishment and spill in period*i*should both be zero as Equation (18):

Given a list of discrete values of *BZ*_{2} with a certain increment in the feasible region, discrete values of *W*_{1} and *W*_{2} can be determined and substituted into the two subsystem models. DP is used to derive a list of *f*_{2} ∼ *W*_{2} and *f*_{1} ∼ *W*_{1} relationships and a list of corresponding solutions, [*X*_{1,i}, *PS*_{1,i}, *Y*_{1,i}] and [*X*_{2,i}, *PS*_{2, i}, *Y*_{2,i}].

3. Aggregation of the system

According to the relationships between the objective function and annual available water of subsystem, *f*_{2} ∼ *W*_{2} and *f*_{1} ∼ *W*_{1}, an aggregation model can be derived as Equations (20) and (21):

The aggregation model can also be solved by DP, with the annual available water of subsystems *W _{j}* as decision variables. After optimal results of the aggregation model

*F**,

*W*

_{1}

^{*}, and

*W*

_{2}

^{*}are derived, the final optimal operation scheme [

*X*

_{1,i},

*PS*

_{1,i},

*Y*

_{1,i},

*X*

_{2,i},

*PS*

_{2,i},

*Y*

_{2,i}]* of the double-reservoir-and-double-pumping-station system can be obtained by searching optimal results of subsystems corresponding to

*W*

_{1}

^{*}and

*W*

_{2}

^{*}.

## RESULTS ANALYSIS

As shown in Figure 5(a), the objective function value *F* decreases with increasing *BZ*_{2} in the range of [1,100, 1,260] and begins to converge when *BZ*_{2} exceeds 1,220 (10^{4} m^{3}). This means that the benefits of increasing the water transfer volume of the HZ pumping station are no longer remarkable from this point.

Finally, the minimum value of *F* is 241, and the optimal annual pumping volume of the HZ pumping station is 1,245 (10^{4} m^{3}); the corresponding optimal operation results are shown in Table 6. Compared with the results using the standard operation policy (SOP), the total water shortage is reduced by 87.7%, total water replenishment from the river is reduced by 2.2%, and total water spill is reduced by 60.6% at the 75% probability of exceedance.

Method | Reservoir | Water supply | Water spill | Water replenishment | Water shortage | Evaporation |
---|---|---|---|---|---|---|

SOP | SH | 1,220 | 241 | 360 | 0 | 278 |

HWB | 1,664 | 0 | 1,100 | 154 | 202 | |

Optimization method | SH | 1,220 | 95 | 352 | 0 | 272 |

HWB | 1,799 | 0 | 1,245 | 19 | 211 |

Method | Reservoir | Water supply | Water spill | Water replenishment | Water shortage | Evaporation |
---|---|---|---|---|---|---|

SOP | SH | 1,220 | 241 | 360 | 0 | 278 |

HWB | 1,664 | 0 | 1,100 | 154 | 202 | |

Optimization method | SH | 1,220 | 95 | 352 | 0 | 272 |

HWB | 1,799 | 0 | 1,245 | 19 | 211 |

The optimal operation redistributes water resources between two reservoirs (Table 6). Furthermore, the optimal operation not only reduces water shortages of the HWB reservoir but also adjusts its water shortage distribution. As shown in Figure 5(b), the water shortage of the HWB reservoir using SOP reaches up to 1.34 MCM in June, accounting for 87% of the annual water shortage. This shortage can be attributed to the mass water consumption during the steeping period of paddy fields. After optimization, although the frequency of water shortage increases, any severe drought in a period could be avoided. In addition, the optimal operation also changes the evaporation of the system, corresponding to changes in the water storage during each period (Table 6). Nevertheless, the impacts of the changes on the local ecological environment can be ignored because the variation range is within 5%.

Operation curves of reservoirs derived using the optimization method and SOP are shown in Figure 6(a) and 6(b), respectively. As a donor reservoir, water storage of the SH reservoir before the flood season was lower with the optimization than with the SOP, which could help retain more flood water and reduce water spill from 2.41 MCM to 0.95 MCM. In contrast, as a recipient reservoir, the optimal water storage of the SH reservoir before the flood season was higher with the optimization than with the SOP, which could increase available water in the reservoir and avoid water shortages caused by the limited pumping capacity of the HZ pumping station during the peak period of irrigation.

According to the differences in the operation curves of the two reservoirs, operation modes of the two pumping stations in the system (XZ and HZ) were also adjusted. As a replenishment pumping station for the SH reservoir, water pumping patterns using different methods are shown in Figure 6(c). The frequency of water pumping at the XZ station did not change after optimization, but the optimal operation reduced the water pumping volume at the first stage, which was beneficial for relieving the pressure of flood control in the SH reservoir. As a water transfer pumping station for the SH reservoir as well as a replenishment pumping station for the HWB reservoir, different water pumping patterns of HZ station are shown in Figure 6(d). Under the SOP, water was transferred only from June to September, the peak period of irrigation. As a result, the pumping station was fully loaded, but it was still difficult to meet the water demand because of the limited maximum pumping capacity. After optimization, the HZ pumping station transferred 15.5% of the annual water pumping volume before the peak period; this could help the SH reservoir retain more flood water as well as increase available water in the HWB reservoir for irrigation.

In conclusion, the most appropriate approach for the optimal joint operation of the double-reservoir-and-double-pumping-station system is to form a joint operation policy. This would enable full play of the functions of all components of the system, making use of surplus water of the SH reservoir to appropriately compensate water shortages of the HWB reservoir. In this manner, the goal of simultaneously reducing water spills and water shortages can be achieved.

## DISCUSSION

At present, meta-heuristic algorithms have become mainstream methods to solve optimal operation models of water resources systems (Bozorgi *et al.* 2017; Mansouri *et al.* 2017). Meta-heuristic algorithms, for instance, genetic algorithm, particle swarm optimization, and ant colony algorithm, firstly generate initial solutions and modify them in each iteration according to certain rules until meeting the convergence (Bayat *et al.* 2011; Chang *et al.* 2013; Ehteram *et al.* 2017). However, they cannot guarantee globally optimal results due to the ineluctable random sampling during the iteration process. Comparatively, dynamic programing based on Bellman's principle is adopted to solve the sub models and the aggregation model after the decomposition and aggregation of the system which can ensure the global optimal solution.

In addition, the solving results of meta-heuristic algorithms may depend on the value of some algorithm parameters (Allawi *et al.* 2018; Ehteram *et al.* 2018), such as crossover rate and mutation rate in GA and inertia weight and acceleration factor in PSO. Therefore, it is necessary to carry out sensitivity analysis on these parameters for the final results. However, there is no need to input any algorithm parameter when using DP which can make the solving procedure more simple and reliable (Kennedy 1986; Peng *et al.* 2018).

These findings provide a reference for the optimal operation of similar systems that consist of reservoirs and pumping stations in hilly regions of southern China, northern Vietnam, Myanmar, and Laos. The mean annual rainfall is over 1,000 mm in these regions, but 70%–90% of that occurs in the flood season due to the monsoon climate, resulting in seasonal water shortages (Wang *et al.* 2017; Quinn *et al.* 2018). What is worse, in recent years, agricultural irrigation, as a traditional water user, has been greatly squeezed because of the rapid development of industry in south China and southeast Asia (Zhang *et al.* 2018). However, it is difficult for the original empirical operation method to handle these challenges. Under the limited agricultural water rights, the optimization model proposed in this paper can effectively adapt to the changes of water resources and alleviate seasonal water shortages by inflow regulation. In future research, the model should consider the stochastic nature of both inflows and water demands for real-time system optimization.

## CONCLUSION

An optimization model of a double-reservoir-and-double-pumping-station system was developed in this study to minimize water shortages. One such system in Nanjing, China, was selected for a case study, and the constraint of total available water in the system was coupled into it, considering regional water rights. The system was decomposed into two subsystems based on the decomposition-aggregation theory, and DP was adopted to solve the sub models and the aggregation model. A joint operation rule of the system is integrated into the optimization model to obtain an appropriate operation scheme of both reservoirs and pumping stations.

## ACKNOWLEDGEMENTS

This work was supported by the National Key Research and Development Program of China [grant number 2017YFC0403205] and the Postgraduate Research & Practice Innovation Program of Jiangsu Province [grant number KYCX19_2105]. The authors wish to acknowledge the managers of Liuhe Water Authority, Nanjing, China, for their support of providing reservoir and irrigation data.