Modeling and Solving of Joint Flood Control Operation of Large-Scale Reservoirs: A Case Study in the Middle and Upper Yangtze River in China

Flood disasters are the most frequent and most severe natural disasters in most countries around the world. Reservoir flood operation is an important method to reduce flood losses. When there are multiple reservoirs and flood control points in the basin, it is difficult to use reservoirs separately to fully realize their flood control potential. However, the multi-reservoir joint flood control operation is a multi-objective, multi-constrained, multi-dimensional, nonlinear, and strong-transition feature decision-making problem, and these characteristics make modeling and solving very difficult. Therefore, a large-scale reservoirs flood control operation modeling method is innovatively proposed, and Dynamic Programming (DP) combined with the Progressive Optimality Algorithm (POA) and Particle Swarm Optimization (PSO) methods, DP-POA-PSO, are designed to efficiently solve the optimal operation model. The middle and upper Yangtze River was chosen as a case study. Six key reservoirs in the basin were considered, including Xiluodu (XLD), Xiangjiaba (XJB), Pubugou (PBG), Tingzikou (TZK), Goupitan (GPT), and Three Gorges (TG). Studies have shown that DP-POA-PSO can effectively solve the optimal operation model. Compared with the current operation method, the joint flood control optimal operation makes the flood control point reach the flood control standard, moreover, in the event of the flood with a return period of 1000 years, Jingjiang, the most critical flood control point of the Yangtze River, does not require flood diversion, and the volume of flood diversion in Chenglingji is also greatly reduced.


Introduction
Since ancient times, flood disasters have been the most frequent natural disasters, with the most severe losses, the most affected people, and the most extensive impacts faced by humans. In order to avoid or reduce flood disasters, people have built a large number of flood control projects such as dikes, reservoirs, flood storage and detention areas, among which reservoirs are a very important part [1][2][3]. Generally, the purpose of reservoir operation includes flood control, water supply, hydropower generation, shipping, ecology, tourism, fishery, etc., among which flood control is the most important function of most reservoirs [4,5]. The purpose of flood control is to minimize the flood disaster loss of flood control points protected by the reservoir under the premise of ensuring the safety of the dam itself. The key variables to control the operation of the reservoir are the remaining flood control capacity of the reservoir and the forecasted flood inflow [6]. Some studies have established a multi-reservoir joint flood control operation method with the minimum flood disaster loss as the objective function [7]. At present, most of the research of multi-reservoir flood control operation focuses on minimizing the maximum releases If the reservoirs are operated alone, it is difficult to realize the flood control potential of the reservoirs due to the lack of cooperation between the reservoirs [16][17][18]. Therefore, it is necessary to establish a joint flood control operation method for reservoirs [19]. However, the multi-reservoir joint flood control operation is a multi-objective, multi-constrained, multi-dimensional, nonlinear [20], and strongly coupled decision-making problem. These characteristics cause modeling and solving are very difficult [7]. Zhou et al. [21,22] proposed a virtual reservoir method, which simplified multiple reservoirs into one virtual reservoir for decision-making. Although the decision-making dimension is reduced, the hydraulic connection between reservoirs cannot be considered, which cannot be ignored, especially in a large basin. He et al. [8] designed a flood control operation method which is based on equal water storage, and successfully applied to flood control in the middle and upper Yangtze River. But this method does not attain the optimal allocation of reservoir flood control storage capacity. Due to this defect, the flood control effect of the reservoir has not reached the ideal goal. Therefore, this paper proposes an optimal flood control operation model with multiple reservoirs and multiple flood control points considering hydraulic connections.
Over the years, many optimization algorithms have been put into application for the solution of the reservoir optimal operation model. These applied algorithms are divided into traditional optimization algorithms and artificial intelligence algorithms. Traditional optimization algorithms include Dynamic Programming (DP) [23,24], Progressive Optimality Algorithm (POA) [25][26][27], Linear Programming (LP) [28], Nonlinear Programming (NLP) [29], artificial intelligence algorithms include Particle Swarm Optimization (PSO) [30][31][32], Genetic Algorithm (GA) [33][34][35], and Differential Evolution Algorithm (DE) [12,36,37], Non-dominated Sorting Genetic Algorithm III (NSGA-III) [1,38], etc. When the number of reservoirs and flood control points is small, these algorithms can achieve an efficient solution to the model, but as the scale grows, the model goals, constraints, and dimensions increase exponentially, and it is difficult to solve the model directly using these algorithms. For example, when the number of reservoirs is large, DP suffers from the "dimensional disaster" [39], POA is greatly affected by the initial solution and easily falls into the local optimal solution [40], PSO locally converges quickly and leads to premature convergence [41], and LP is difficult to solve strong nonlinear problems through linear fitting [42]. For the DP algorithm, it can obtain the global optimal solution but has the "dimensional disaster", the POA algorithm has a small amount of calculation but is easy to fall into the local optimal solution, and the PSO search ability is strong but the convergence is difficult when the dimensionality is large. Based on the advantages and disadvantages of each algorithm, this paper proposes an algorithm that combines DP, POA, and PSO to solve a large-scale reservoirs joint flood control optimal operation model. This algorithm reduces the scale of DP calculation and avoids the "dimensional disaster". At the same time, it combines POA and PSO to improve the global search capability while avoiding falling into the local optimal solution.
The purpose of this research is to establish a flood control operation model for a flood control system composed of multiple reservoirs and flood control points, considering the flood routing, to minimize the system flood loss and balance the flood control requirements of each point and reservoir. A combination optimization algorithm is proposed to complete the efficient solution of the model, and the flood control system composed of six reservoirs and eight control points in the middle and upper reaches of the Yangtze River is taken as a case study.

Reservoirs Joint Flood Control Optimal Operation Model
The primary objective of reservoirs joint flood control operation is to minimize the total flood loss and risk of the system during the entire flood season while subjecting to several equality and inequality constraints. The objective function and constraints are described as follows.

Objective Function
Usually, reservoir joint flood control optimal operation has the following objectives: (1) Minimize the amount of flood diversion at the flood control point, this objective is to reduce the flood disaster loss of flood control points; (2) The minimum sum of squares of flow at the flood control point. This objective is to reduce the pressure of flood control points, which is used in the situation where the flood control point does not require flood diversion; (3) The highest operating water level of the reservoir is the lowest. This objective is to reduce the flood control risk of the reservoir itself; (4) The sum of squares of the inflow of the reservoir is the smallest. This objective is to avoid drastic changes in the inflow of the downstream reservoir caused by the operation of the upstream reservoir.
(5) The sum of the square of the outflow of the reservoir is the smallest. This objective is to avoid drastic changes in the outflow flow of the reservoir from affecting power generation, water supply, etc.
In this model, the above five objectives are combined to be one objective by using the weight method, the final objective function is shown as follows: where N is the total number of reservoirs; Z max n is the maximum allowable operating water level of the n-th reservoir; Z min n is the minimum allowable operating water level of the n-th reservoir; Z up n is the highest operating water level of the n-th reservoir; T is the total number of operation periods; Q n,t,in is the inflow of the n-th reservoir during the t-th period; Q n,t,out is the outflow of the n-th reservoir during the t-th period; Q n,max is the maximum inflow of the n-th reservoir when reservoirs are not in operation; M is the total number of flood control points; Q m,t is the inflow of the m-th point during the t-th period; Q m,max is the maximum inflow of the m-th point when reservoirs are not in operation; Q m,t,E is the flood diversion flow of the m-th point during the t-th period, which is equal to zero when the flow is less than the safe flow of the point; Q k m,t,E is the flood diversion flow of the m-th point during the t-th period when reservoirs are not in operation; ∆t is the length of operation period, which usually depends on the available data for operation; a n , b n , and c n is the weight of the highest operating water level, the square sum of inflow and square the sum of outflow of the n-th reservoir, respectively; d m and e m is the weight of the square sum of flow and the square sum of flood diversion flow of the m-th point, respectively. The weight value is determined according to the importance of the objective and the requirements of the project, it is given in advance based on the experience of the decision maker, and the weight is appropriately modified according to the operation result until a satisfactory result is achieved. Generally, engineers give multiple sets of weights and their corresponding operation results based on experience, and then the decision maker selects a set of better weight values from them. When the decision maker is dissatisfied with all the operation results, the weights are appropriately modified after consultation with the engineer. A satisfactory operation result needs to meet the requirements of all reservoirs and flood control points as much as possible. When the requirements of all reservoirs and flood control points cannot be guaranteed at the same time, the weight should be adjusted to meet the requirements of high-level reservoirs and flood control points first.

Constraints
In the process of flood control operation, various complex equality and inequality constraints, such as water balance, water level, output, and hydraulic connection, should be considered for restricting objective optimization. The constraints of the model are described as follows: (1) Water balance constraint: where V n,t is the reservoir storage of the n-th reservoir at the end of the t-th period.
(2) Water level constraint: where Z n,t is the water level of the n-th reservoir at the end of the t-th period; Z min i,t and Z max i,t are the maximum and minimum allowable operating water level; ∆Z n is the maximum amplitude of water level variation.
(3) Reservoir discharge constraint: |Q n,t,out − Q n,t−1,out | ≤ ∆Q n,out (6) where Q min n is the minimum discharge of the n-th reservoir; Q max n (Z n,t ) is the maximum discharge ability of the n-th reservoir at t-th period. The maximum discharge ability is related to the water level; ∆Q n,out is the maximum amplitude of discharge of the n-th reservoir.
(4) Constraint of discharge for downstream flood control points: where Q m,t is the flow of the m-th point at the t-th period; H is the number of reservoirs upstream of the m-th point; Q h,m,t is the discharge of the h-th reservoir considering the flood routing at the m-th point at the t-th period; ∆q m,t is the inflow in the interval zone between the reservoirs upstream and the m-th point; and Q max m is the allowed maximum discharge for the m-th point. This constraint is often used to turn the objective of non-critical points into constraints to reduce the difficulty of solving the model.
(5) Constraint of flood routing: The Muskingum method [43,44] is the most commonly used method of river channel flood routing and is applied to this model. Where O t and I t are the outflow and inflow of the river channel in t-th period, respectively; O t−1 and I t−1 are the outflow and inflow of the river channel in t −1 -th period, respectively; K is the storage time constant for the river channel, which has a value reasonably close to the flow travel time within the river channel; and x is a weighting factor varying between 0 and 0.5; C 0 , C 1 , and C 2 are coefficients that are functions of K, x, and ∆t.

Optimization Algorithm
For the solution of the joint flood control optimal operation model of multi-reservoir, considering the influence of "curse of dimensionality", flood routing, calculation efficiency, initial solution, and convergence, etc., a strategy combining DP, POA, and PSO is proposed. First, DP is used to calculate the reservoir one by one to obtain the initial solution. Second, considering the influence of flood routing, POA is used to optimize the initial solution. Third, from the beginning to the end, PSO is used to optimize the operation process for several consecutive periods in turn. The skeleton of this algorithm is shown in Figure 2. The equations for DP, POA, and PSO methods are given as follows. DP is widely used in reservoir operation of hydropower generation, flood control, water storage, and so on. This paper uses DP to obtain the initial solution of reservoir operation. The recursive function of DP used in the multi-period optimal operation of the reservoir is described as Equation (10); DP is used to optimize reservoir operation from upstream to downstream. In this stage, the influence of river channel flood routing is not considered. The result of DP is used as the initial solution of the model.
where V n,t (Z n,t ) is storage of the n-th reservoir at t-th period; G t (Z n,t ) is the minimum cumulative utility from the G t (Z n,t ) is the minimum cumulative utility from the first period (begin of operation period) to the t-th period (current operation period); g t (Z n,t , Q n,t,out ) is the t-th period utility function, when the DP is used to optimize the operating water level of the n-th reservoir at t-th period, the operating water level in other periods remains unchanged, so only the utility value of the reservoir during the t-th period is calculated, which effectively reduces the amount of calculation [23]; Q h,m,t is the discharge of the h-th reservoir without considering the flood routing at the m-th point at the t-th period, Q h,m,t in Equation (7) that uses the Muskingum method for flood routing is simplified into Q h,m,t .
The detailed steps at this stage are as follows: Step 1: According to the operation requirements, the water level of the reservoirs is discretized into a water level sequence, and the discrete water level is set as a state variable. In this paper, the discrete step length of the reservoir is 0.01 meters.
Step 2: From upstream to downstream, the DP is used to optimize the reservoir operation process in turn.

Progressive Optimality Algorithm (POA)
POA can effectively avoid the "dimensional disaster" suffered by DP, and it can also consider the flood routing. On the basis of the initial solution, the multi-stage optimal operation problem is decomposed into a two-stage optimal operation sub-problem. The decision vector of sub-problem between the t-th and t+1-th period is described in Equation (11): where C n,t is the range that Z n,t can take at the t-th period and is determined by constraint 3, Q n,t−1,out and Q n,t,out are calculated by s.t.(1) using Z n,t−1 , Z n,t , and Z n,t+1 . F(Z n,t−1 , Z n,t , Z n,t+1 ) is the optimal objective at the t-th and t+1-th period. The detailed steps at this stage are as follows [45,46]: Step 1: Select the n-th reservoir from upstream to downstream for optimization calculation, and the operating water levels of other reservoirs remain unchanged.
Step 2: The water level of the n-th reservoir from the second period to the T −1 -th period is optimized in turn. When the water level of the t-th period is optimized, the water level of the n-th reservoir remains unchanged in other periods, the optimization objective is shown in Equation (11).
Step 3: After optimizing the operation water level of the n-th reservoir is completed, the operating water level of the n + 1-th reservoir is optimized using the method in step 2.
Step 4: After optimizing the operation water levels of all reservoirs, repeat steps 1, 2, and 3 until the optimization results remain stable.

Particle Swarm Optimization (PSO)
PSO is derived from the research on the predation behavior of bird swarms. Each particle in the particle swarm represents a possible solution to the problem. Through the simple behavior of individual particles, the information interaction within the particle swarm realizes the intelligence of problem solving. Due to its simple operation and fast convergence speed, PSO has been widely used in many fields such as function optimization, image processing, and geodesy [31]. The paper uses PSO to solve the problem that POA is easy to fall into local optimal solutions. However, if the optimization calculation is performed during the entire operation period, the length of the particles will be too long, thereby reducing the convergence speed. Therefore, this paper combines POA, and PSO, and selects several consecutive operation periods as the "operation window". PSO is used to optimize the operation process in the "operation window", and keeps the operation process in other periods unchanged. When the operation process optimization in the "operation window" is completed, the window moves backward and optimizes the operation process. Each particle in PSO evolves according to the set strategy. The formula for updating the velocity and position of the particles is as follows: where X id (r) is the position of the i-th particle in the r round of optimization, the value of X id (r) in each dimension represents the operating water level of a reservoir in a period; V id (r) is the velocity of the i-th particle with d dimensions in the r round of optimization, the value of V id (r) in each dimension represents the change value of the operating water level of a reservoir in a period; P id is the best position where the i-th particle has appeared so far; P gd is the best position where all particles have appeared so far; C 1 and C 2 are the learning constants; w is the weight of inertia; random(0, 1) is a random number between 0 and 1. The detailed steps at this stage are as follows [47][48][49]: Step 1: Based on experience, the "operation window" is set to contain k periods. If k is set to be relatively large, the convergence speed will be slower. On the contrary, if k is set to be relatively small, the search ability will be reduced. The recommended interval is [10,30], and k is set to 20 in this paper; Step 2: Set the number of particles, and initialize the velocity and position of particles. Each particle contains the water level of k periods of N reservoirs in the "operation window", so the particle contains k × N dimensions. The result of POA in the "operation window" is used as one of the particles.
Step 3: According to the fitness function shown in Equation (14), calculate the fitness value of each particle.
where j is the starting period of the "operation window". When optimizing the reservoir operating water level process in the "operation window" from the j-th period to the j+k-th period, the water levels in other periods remain unchanged, so only the periods in the window need to be considered when calculating the fitness. If j and k are set to 1 and T respectively, it can be discovered that this fitness function is the objective function.
Step 4: According to the fitness value of the particles, update the historical best position of each particle and the best position of the particle swarm [31].
Step 5: Update the velocity and position of each particle according to Equations (12) and (13).
Step 6: Repeat steps 3, 4, and 5 until the best position of the particle swarm remains stable.
Step 7: The "operation window" moves backward p periods, and the water level from the j+p-th period to the j+p+k-th period in the "operation window" is optimized. p is recommended not to be more than one-third of k, which is taken as 5 in this paper.

Introduction of Study Area
The middle and upper Yangtze River is selected as a case study. As the largest river in China, the Yangtze River has a total length of about 6300 km. The area, population, GDP, and grain output of the Yangtze River Basin account for 18.8%, 36%, 40%, and 40% of the country respectively [50]. However, floods in the Yangtze River Basin often occur, killing thousands of people and causing millions of dollars in damage. For example, the flood of 1954 directly caused 30,000 deaths, and more people died of starvation and plague caused by the flood [51]. Therefore, flood control in the Yangtze River Basin has always been the focus of general attention, and flood prevention and disaster reduction in the Yangtze River Basin has been given top priority. The middle and upper Yangtze River are the most severely flooded areas and the key to flood control in the entire basin.
After decades of flood control construction, the Yangtze River Basin has basically formed a comprehensive flood control system with the Three Gorges Reservoir (TG) as the backbone. However, the contradiction between the river discharge capacity and the large flood peaks of the middle and lower reaches of the Yangtze River is still prominent.  Figure 3b, control point is the object that needs to be protected by the flood control operat reservoirs, and the gauging station is the station for observing the water level and f the river.

Data of Reservoir
XLD is located in the lower reaches of the Jinsha River and is the third-l hydropower station in the world. The height, total storage capacity, and flood c storage capacity of the dam are 285.5 meters, 12.67 billion m 3 , and 4.65 billio respectively. XJB is located at the downstream of XLD, which is only 157 km far fr The flood control storage capacity of XJB is 0.90 billion m 3 . The two reservoirs a reservoirs of the Jinsha River, controlling 96% of the Jinsha River basin area.  Table 1.

Data of Flood Control Point
The Chengkun Railway connects Sichuan Province and Yunnan Province with a total length of 1096 km. It is the most important transportation line in southwest China. According to the flood control standard stipulated by law, when the flow is not higher than 11,200 m 3 /s, which is the peak flow of the design flood with a return period of 100 years, PBG should be used to intercept the flow below 8780 m 3 /s, which is the maximum flow that the Chengkun Railway can withstand.
Yibin, Luzhou, and Langzhong are three cities located in Sichuan Province. The populations of the three cities are 5.5 million, 5.8 million, and 0.83 million respectively, and GDP are 37.7 billion, 31.1 billion, 3.7 billion dollars respectively. The current flood control capacity of the three cities is to ensure safety in the event of a flood with a return period of 20 years, the corresponding peak flows of the flood are 51,000 m 3 /s, 52,600 m 3 /s, and 25,100 m 3 /s, respectively. According to the flood control standard stipulated by law, upstream reservoirs should protect the three cities from the design flood with a return period of 50 years, the corresponding peak flow of the flood is 57,800 m 3 /s, 58,100 m 3 /s, and 32,880 m 3 /s.
Chongqing is the largest city in southwestern China and the center of economy, finance, science and technology, shipping, and logistics in the upper reaches of the Yangtze River. It has a population of 31.2 million, a GDP of 352.3 billion dollars, and an area of 82,400 km 2 . The discharge capacity of the river channel in Chongqing is 83,100 m 3 /s which is the peak flow of the design flood with a return period of 50 years, and the upstream reservoirs should protect Chongqing from the design flood with a return period of 100 years, the corresponding flood peak flow is 88,700 m 3 /s, according to the flood control standard stipulated by law.
Sinan is a county located in the northeast of Guizhou Province, with an area of 2230.5 km 2 , a population of 684,500, and a GDP of 2 billion dollars. According to the flood control standard stipulated by law, when the flow is not higher than 11,500 m 3 /s which is the peak flow of the design flood with return period of 5 years, GPT should be used to intercept the flow below 9320 m 3 /s which is the peak flow of the design flood with a return period of 2 years.
Jingjiang is an important production base for grain, cotton, oil, fruits, and aquatic products in China. It has a population of 9.6 million, 16.6 million acres of arable land, many important railways and highways. Once floods occur, they will directly threaten the safety of Wuhan which is the most important city in central China. On the other hand, the flood discharge capacity of the Jingjiang is only 60,600 m 3 /s, which is seriously insufficient compared to the huge amount of inflow. Therefore, the Jingjiang River has always been the difficulty and focus of the Yangtze River flood control. Chenglingji is located in Yueyang City, Hunan Province. It is the exit of Dongting Lake into the Yangtze River. Chenglingji is the area where flood disasters occur most frequently in the Yangtze River Basin.
In this flood control system, the flood control objectives of the Chengkun Railway, Yibin, Luzhou, Chongqing, Langzhong, and Sinan are transformed into constraints to ensure that they meet the flood control standards. Then, this paper focus on solving the flood control problems of Jingjiang and Chenglingji.

Muskingum Parameters of the Basin
The Muskingum parameters of each reach in the basin are shown in Table 2. The flood routing time of each river segment is 6 h.

Design Flood Hydrographs
The flood control operation of the reservoirs is to ensure the safety when a large or extreme flood occurs [52]. The common method is to determine the flood control standard based on the importance of the flood control point, and then the design flood corresponding to the standard is used as the input of the operation model [28,53]. There are two methods commonly used in design flood calculation: homogeneous frequency enlargement and homogeneous multiple enlargement. The homogeneous frequency enlargement is widely used in the case of a large number of reservoirs and flood control points because of its ease of use [54][55][56]. The steps of this method to calculate the design flood are as follows: Step 1: The design flood parameters of the control section are obtained by calculating the frequency of the annual maximum flood peak and n-day volumes data. The control section of this study is the Yichang gauging station, and its design flood parameters are shown in Table 3. In the table, E x , C v , and C s are the mean value, coefficient of variation, and coefficient of skewness of design flood respectively [54].
Step 2: The typical flood hydrographs were determined. An observed flood hydrograph with a large peak, large flood volume, or an operation not conducive to flood control was usually selected [57]. This study selected floods in 1954,1968,1980,1981,1982,1988,1996,1998,1999 Step 3: The design flood hydrographs were determined. According to the ratio of the peak or n-day volumes between the observed typical flood and the design flood at the control section, the typical flood hydrographs of reservoirs and flood control points in the basin are enlarged and become the basin design flood hydrographs. Usually, the flood volume ratio close to the peak ratio is selected for enlarging to ensure that the flood peak and flood volume of the design flood match. In this study, the 1954 and 1998 floods used the ratio of 30-day volume, and other floods used the ratio of 15-day volume.

Optimal Reservoir Operation Results
In this study, a joint flood control optimal operation model of the flood control system composed of six reservoirs and eight points in the upper and middle reaches of the Yangtze River was established, and 1954, 1968, 1980, 1981, 1982, 1988, 1996, 1998, 1999, and 2010 typical year design floods with different return periods were used as the inputs of the model. Then, the DP-POA-PSO proposed in this paper is used to solve the optimal model. The system has three main flood control goals: first, the flood diversion volumes of Jingjiang is reduced to the lowest; second, the flood diversion volumes of Chenglingji is reduced to the lowest; third, the highest operation water level of TG is the lowest. The importance of the three goals decreases in order. The programming tool, programming language, and database use eclipse, java, and MySQL respectively. On a computer with a CPU frequency of 4.0 GHz, the calculation time of using this algorithm to solve the model is about 180 s.

Influence on the Inflow Flood of TG
The operation of the upstream reservoir can reduce the inflow flood of the TG. The operation results of this part are shown in Table 4. In the current method, each reservoir is operated separately and no joint optimal operation is realized. In the table, the design value is the inflow peak of TG without reservoir operation; the reduction is the difference between the design value and the inflow peak after reservoir operation; the comparison is the difference between the peak flows of the current method and the optimal operation. The operation results show that the joint optimal operation of the reservoirs can make full use of the flood control capacity of the upstream reservoirs to reduce the inflow of TG and reduce its flood control pressure. Therefore, the flood control effect of TG on Jingjiang and Chenglingji has been greatly improved. For the designed floods with return periods of 1000 years, 100 years, and 50 years, compared with the current method, the optimal operation increases the inflow peak reduction by 68%, 37% and 36%, respectively. Figure 6 shows the inflow of TG when different methods are used to operate the 1980 and 2010 typical year design flood with a return period of 1000 years. The results show that the joint optimal operation can make the inflow flatter and the peak smaller, which will increase the flexibility of TG for flood control operation of Jingjiang and Chenglingji.

Operation Results of Jingjiang and Chenglingji
Considering that the Chengkun Railway, Yibin, Luzhou, Langzhong, and Sinan can all meet the flood control standards after reservoir operation, many studies [58][59][60][61] have been focusing on the operation effect of the reservoir group on Jingjiang and Chenglingji, especially the flood diversion of these two areas. In order to prove the feasibility of the optimal flood control operation model and DP-POA-PSO algorithm proposed in this paper, the results of the optimal operation model are compared with the current methods. The flood diversions in Jingjiang and Chenglingji are shown in Tables 5 and 6, respectively. The Table 5 show that the flood control optimal operation model of reservoirs can greatly improve the flood control safety of the Jingjiang. When floods with return period of 50 or 100 years occur, the flow after the reservoir operation does not exceed the discharge capacity of the Jingjiang, and there is no need to use the flood diversion areas. In the event of the flood with return period of 1000 years, the optimized model reduced Jingjiang's flood diversion volumes by 25 On the other hand, Table 6 shows that the flood control safety of Chenglingji area has also been greatly improved. In the event of floods with return period of 50 years, there is no need to use the flood diversion areas in Chenglingji, this will greatly reduce the pressure of flood diversion in Chenglingji. Compared with the current method, the Chenglingji's flood diversion volumes has been reduced by 0.16, 0.25, and 1.96 billion m 3 respectively on average when floods with return period of 50, 100, and 1000 years occur.
From the operation results of Jingjiang and Chenglingji, it can be seen that the flood control optimal operation model of reservoirs and DP-POA-PSO algorithm proposed in this paper can significantly reduce the flood diversion of these two flood control points compared with the current method. When the flood is not large, the goal of no flood diversion can be achieved, which will greatly reduce the flood disaster loss.

Operation Results of TG
The highest operating water level of TG is shown in Table 7. The lower water level can reduce the flood control risk of the reservoir itself and increase the flexibility of flood control operation. Compared with the current method, the water level has been reduced by 1.72, 1.58, and 0.38 m respectively on average when floods with return period of 50,100, and 1000 years occur. In the event of a flood with a return period of 1000 years, the maximum highest water level is 175 m, and the minimum highest water level is 171.51 m. These results show that the flood control capacity of the reservoir is fully used in order to reduce the flood diversion of Jingjiang and Chenglingji when a large flood occurs. Conversely, when the flood with return period of 50 years occur, the maximum and minimum highest water level are 166.24, 158.57 m, respectively. These results show that when the flood is not large, the highest water level of the reservoir is reduced as much as possible under the premise of ensuring the minimum flood diversion volumes of Jingjiang and Chenglingji to reduce the flood control risk of the reservoir itself.  Figure 7 shows the results of the TG under different operation methods, using as the designed flood the 1954 typical year with return periods of 50 and 1000 years. The figure shows that the results obtained by the optimal operation model are usually better than that of the current method, because inflow and outflow of TG are much smaller and much flatter. In particular, the operation water level of TG is also lower.
It can be seen from Tables 4-7 that for 10 design floods with different encounter characteristics, the model and algorithm proposed in this paper can effectively reduce the peak inflow of the reservoir, the highest operating water level of the reservoir, and the flood diversion of flood control points. On the other hand, when the flood is large, the goal of reducing the amount of flood diversion is achieved. When the flood is small, the operating water level of the reservoir is further reduced, and the flood control risk of the reservoir itself is reduced. These show that the method has good adaptability.

Conclusions and Discussion
In this study, a joint flood control optimal operation model with multiple reservoirs and flood control points was established, which gave full play to the flood control potential of the reservoirs, and balanced the flood control demand between reservoirs and flood control points. The impact of transposition and attenuation of flood cannot be ignored, so the Muskingum routing method is used by the operation model. A DP-POA-PSO optimization algorithm is proposed to solve the optimal operation model. Then, a flood control system consisting of six reservoirs and eight flood control points in the middle and upper Yangtze River is selected as a case study. The study mainly aims at minimizing the flood diversion volumes of Jingjiang and Chenglingji, the maximum inflow of TG and the highest operating water level of TG. The main conclusions are summarized as follows.
(1) This article has two main innovations. The first innovation is that an optimal operation model for multiple reservoirs and flood control points was proposed. The flood control potentiality of the reservoirs is fully utilized, and the flood loss and risk of the system are significantly reduced. In addition, compared with simplified methods such as time delay, the use of the Muskingum method makes the river channel flood routing more accurate, which makes the operation result more in line with the actual project. The second innovation is the proposed DP-POA-PSO optimization algorithm, which combines the advantages of each algorithm in solving the optimal operation of the reservoir, avoiding the "dimensional disaster", local optimal solution, slow convergence speed, and so on. It shows excellent performance in practical applications.
(2) The optimal operation model and algorithm proposed in this paper can significantly reduce the inflow of the Three Gorges Reservoir. For the design flood with return periods of 1000 years, 100 years, and 50 years, the average reduction of inflow peak is 13,311 m 3 /s, 10,258 m 3 /s, and 10,027 m 3 /s, respectively, the corresponding rate of reduction is 14.6%,13.0%, and 13.4% respectively. Compared with the current method, the optimal operation increases the inflow peak reduction by 5375 m 3 /s, 2794 m 3 /s, and 2628 m 3 /s respectively.
(3) The joint flood control optimal operation of reservoirs can greatly improve the flood control safety of Jingjiang and Chenglingji. When floods with return period of 50 or 100 years occur, the flow after the reservoir operation does not exceed the discharge capacity of the Jingjiang, and there is no need to use the flood diversion areas. In the event of a flood with a return period of 1000 years, the optimized model reduced Jingjiang's flood diversion volumes by 25.3 billion m 3 on average. Compared with the current method, the Jingjiang flood diversion volume has been reduced by 1.78 billion m 3 on average. On the other hand, in the event of a flood with a return period of 50 years, there is no need to use the flood diversion areas in Chenglingji, this will greatly reduce the pressure of flood diversion in Chenglingji. Compared with the current method, Chenglingji's flood diversion volumes have been reduced by 0.16, 0.25, and 1.96 billion m 3 respectively on average when floods with return period of 50,100, and 1000 years occur.
(4) The results show that the optimal operation model has good adaptability and flexibility. When a large flood occurs, the flood control capacity of the reservoir is fully used in order to reduce the flood diversion of Jingjiang and Chenglingji. Conversely, when the flood is not large, the highest water level of Three Gorges Reservoir is reduced as much as possible under the premise of ensuring the minimum flood diversion volumes of Jingjiang and Chenglingji to reduce the flood control risk of the reservoir itself.
In real-time flood operation, the inflow at each point in the basin can be obtained through forecasting, and then the model and algorithm proposed in this paper are used to operate the forecasted inflow. The operation results can provide program support for realtime flood control operation. In addition, the operation results of design floods can provide guidance for the formulation of the joint flood control operation rules of reservoirs and flood control points. It can be seen from the operation results that although the operation of the six reservoirs has reduced the flood diversion volumes of Chenglingji, the flood control pressure is still very high, so more reservoirs and more efficient algorithms need to be considered to solve this problem.
In addition, there are many uncertain factors in the operation, which will affect the flood control effect of reservoirs. Research by Hu et al. [62], Krzysztofowicz et al. [63], Chen et al. [64], and Vermuyten et al. [65] showed that among all the uncertain factors in flood control operation, the uncertainty of flood forecast has the greatest impact on operation results, especially in long-term operation. The uncertainty of the forecast mainly includes three aspects: the error of inflow, the error of flood volume, and the error of flood hydrograph shape. At the same time, the water level-storage curve and the water level-discharge capacity of reservoir curve also have an influence on operation [65]. In order to reduce the difficulty of model solving and increase the speed of calculation, this study did not consider these uncertain factors. Uncertainty can be considered in future research to further improve the operation effect.