Comparative analysis of some evolutionary-based models in optimization of dam reservoirs operation

Deriving optimal operation policies for multi-reservoir systems is a complex engineering problem. It is necessary to employ a reliable technique to efficiently solving such complex problems. In this study, five recently-introduced robust evolutionary algorithms (EAs) of Harris hawks optimization algorithm (HHO), seagull optimization algorithm (SOA), sooty tern optimization algorithm (STOA), tunicate swarm algorithm (TSA) and moth swarm algorithm (MSA) were employed, for the first time, to optimal operation of Halilrood multi-reservoir system. This system includes three dams with parallel and series arrangements simultaneously. The results of mentioned algorithms were compared with two well-known methods of genetic algorithm (GA) and particle swarm optimization (PSO) algorithm. The objective function of the optimization model was defined as the minimization of total deficit over 223 months of reservoirs operation. Four performance criteria of reliability, resilience, vulnerability and sustainability were used to compare the algorithms’ efficiency in optimization of this multi-reservoir operation. It was observed that the MSA algorithm with the best value of objective function (6.96), the shortest CPU run-time (6738 s) and the fastest convergence rate (< 2000 iterations) was the superior algorithm, and the HHO algorithm placed in the next rank. The GA, and the PSO were placed in the middle ranks and the SOA, and the STOA placed in the lowest ranks. Furthermore, the comparison of utilized algorithms in terms of sustainability index indicated the higher performance of the MSA in generating the best operation scenarios for the Halilrood multi-reservoir system. The application of robust EAs, notably the MSA algorithm, to improve the operation policies of multi-reservoir systems is strongly recommended to water resources managers and decision-makers.


List of symbols F(Re)
Objective function Re i,t Release from the ith reservoir during period t (MCM) De i,t Downstream demands of the ith reservoir during period t (MCM) Maximum downstream demand of the ith reservoir in the whole operation period (MCM) i Reservoir index T Total number of operation periods (month) Penalty Penalty function Sp i,t Overflow from the ith reservoir during period t (MCM) S i,t Storages of the ith reservoir at the beginning of period t (MCM) S i, t+1 Storages of the ith reservoir at the end of period t (MCM) S max i Maximum storage in the ith reservoir during period t S min i Minimum storage of the ith reservoir at the beginning of period t Q i,t Inflow into the ith reservoir during period t (MCM) Loss i,t Evaporation loss from the ith reservoir surface during period t (MCM) Ev i,t Net evaporation (evaporation minus precipitation) from the ith reservoir during the period t (m) Water scarcity in many parts of the world makes the careful management of water resources necessary. Reservoirs are one of the most important elements of water management systems designed for storing and regularly releasing water to meet the downstream demands based on the operation policy. A reservoir operation policy includes a set of rules that assign water release in various operation conditions. Developing optimal operation policy for reservoirs is a complex and high challenging engineering problem, particularly for multi-reservoir systems. The main reasons for such complexity are the stochastic nature of the system inputs, the involvement of many decision variables and constraints, and the interfering operations between successive dams (in multi-reservoir systems). Thus, the implementation of an appropriate reservoir operation policy needs a powerful optimization technique. Over the last few decades, several optimization methods have been proposed to solve such complex engineering problems. These methods can be classified into two main categories of classical methods, and evolutionary algorithms (EAs).
Classical methods such as linear programming (LP), non-linear programming (NLP), dynamic programming (DP), and stochastic dynamic programming (SDP), are old, imprecise and time-consuming techniques that often fail to provide reasonable solutions because they have several limitations, particularly for large scale multi-reservoir systems 1 . For example, the LP only can solve optimization problems with linear objective function and constraints, the NLP may get trap in a local optimum, especially in non-convex optimization problems, and the DP and the SDP suffer from the curse of dimensionality and state-space discretization 2,3 .
EAs are iterative search strategies in which a simulation process is repeated and a set of decision values is updated until the values satisfy a given stopping criteria, and, the optimal solutions are obtained by a decisionmaking process 4 . Over the last years, artificial intelligence-based EAs have been successfully employed to solve complex engineering problems [5][6][7][8] . Regarding their high capability in finding the optimal solutions, they have attracted the attention of researchers and water resources managers.
More recently, Qaderi et al. 9 used a water cycle algorithm (WCA) for optimization of the monthly operation of Golestan and Voshmgir consecutive dams in Iran. They compared the results of the WCA with the genetic algorithm (GA), the harmony search algorithm (HS), the particle swarm optimization (PSO), and the imperialist competitive algorithm (ICA). It was found that the WCA was superior to the other algorithms in achieving the global optimum solutions. Mao et al. 10 successfully employed the shuffled complex evolution (SCE) coupled with the stochastic ranking for reservoir scheduling problems. Ahmadebrahimpour 11 used the wolf search algorithm (WSA) to optimize the operation of a four-reservoir system and a single hydropower system in Iran. He reported the supremacy of the WSA to the GA. Mohammadi et al. 12 combined the whale optimization algorithm (WOA) with the GA to produce a hybrid whale-genetic algorithm with higher precision and convergence rate. They successfully applied the developed model to optimize the operation of two multi-reservoir benchmark systems. Zarei et al. 13 integrated the bat algorithm (BA) and the PSO with the game theory for optimal operation of Shahid Dam Reservoir in Iran and reported the better convergence rate and higher reliability of the developed algorithm. Myo Lin et al. 14 proposed a multi-objective model predictive control (MOMPC) for real-time operation of a multi-reservoir system in Myanmar. Chen et al. 15 used the PSO with an adaptive random inertia weight strategy to optimize Panjiakou Reservoir system operation. The results indicated the higher efficiency of the developed model compared to the GA, the conventional PSO, and the improved versions of PSO. Shaikh 16 utilized five different artificial neural network (ANN) models for optimal operation of Damodar Valley multi-purpose multireservoir system in India and documented the high capability of ANN models in this complex problem. Sharifi et al. 17 developed a fitness-distance-balance moth swarm algorithm for optimization of cascade hydropower reservoirs operation. They found that the developed model could successfully increase the hydropower generation by 59.5% compared to the actual generation of energy over a 180-months operation period.
Recently, five robust EAs of tunicate swarm algorithm, TSA 18 , Harris hawks optimization algorithm, HHO 19 , seagull optimization algorithm, SOA 20 , sooty tern optimization algorithm, STOA 21 , and moth swarm algorithm, MSA 22 , were introduced and successfully employed in solving several real-world complex engineering problems. These algorithms showed very promising and competitive results compared to several well-known metaheuristics. They demonstrated many advantages such as the high convergence rate, the high capability in reaching the optimal values for large-scale problems, the minimal time consumption, the minimal computational cost in terms of CPU usage, and several other superiorities.
In the present study, we developed these recently introduced EAs to optimize the operation of Halilrood large-scale multi-reservoir system. To the best of the authors' knowledge, this is the first application of these algorithms in optimization of water resources systems. To evaluate the efficiency of these algorithms, their results were compared with two popular algorithms of GA and PSO. In addition, four statistical metrics were employed to assess the algorithms' performance.

Methodology
Optimization algorithms. Here a brief introduction to the utilized algorithms is presented. More mathematical details of each algorithm are referred to the original works. All the studied EAs were coded and run in the MATLAB R(2014)a platform.
TSA algorithm. Tunicate Swarm Algorithm (TSA), recently proposed by Kaur et al. 18 , is a new bio-inspiredbased EA that is inspired by the jet propulsion and swarm behaviors of tunicates during the navigation and foraging process in the depth of ocean. Tunicate has the ability to find the location of food sources in the sea. Two behaviors of tunicate, jet propulsion and swarm intelligence, are employed for finding the food source. To mathematically model the jet propulsion behavior, a tunicate should satisfy three conditions, including avoiding the conflicts between search agents, the movement towards the position of the best search agent, and remaining close to the best search agent. The swarm behavior will update the positions of other search agents about the best optimal solution. The details of mathematical modeling of these behaviors were documented by Kaur et al. 18 .
In this study, we developed this algorithm to optimize the operation of Halilrood multi-reservoirs system. Figure 1 illustrates the phases of developed TSA in the optimization process.
HHO algorithm. The Harris Hawks Optimization (HHO) algorithm, proposed by Heidari et al. 19 , was inspired by the collaborative behavior and chasing style of Harris' hawks in nature. The main strategy of Harris' hawks to capture a prey is called "surprise pounce" or "seven kills" tactic. Harris' hawks can conduct various chasing styles depending on the circumstances and escaping patterns of a prey. A switching strategy occurs when the best hawk stoops at the prey and gets lost, and the chase will be continued by one of the party members. The main advantage of these collaborative tactics is that the Harris' hawks can pursue the detected prey to exhaustion, which increases its vulnerability. Besides, by perplexing the escaping prey, it loses its defensive capabilities and eventually, it gets caught by the most powerful and experienced hawks who share the captured prey with other party members. More details about this algorithm can be found in Heidari et al. 19 .
The HHO is a population-based, gradient-free optimization technique; hence, it can be applied to any optimization problem subject to a proper formulation. Accordingly, we developed this algorithm for optimization of Halilrood multi-reservoir systems' operation. Figure 2 shows the flowchart of HHO algorithm in the optimal operation of this problem.
SOA algorithm. Seagull optimization algorithm (SOA), proposed by Dhiman and Kumar 20 , stimulates the migration and attacking behaviors of a seagull in nature. Seagulls, with the scientific name of Laridae, are sea birds that are very intelligent. They use bread crumbs to attract fish and produce rain-like sounds with their feet to attract earthworms hidden under the ground. Generally, seagulls live in colonies. The most important behavior of the seagulls is their migrating and attacking behaviors. During migration, seagulls travel in a group. In the group, seagulls can travel towards the direction of the best survival fittest seagull. Therefore, according to the fittest seagull, other seagulls can update their initial positions. More details of the SOA can be found in Dhiman and Kumar 20 . Figure 3 demonstrates the flowchart of the SOA algorithm, including the migration and attacking phases, for the optimal operation of Halilrood multi-reservoirs system.    21 , is a new evolutionary algorithm inspired by the behaviors of sooty tern in nature. Sooty terns, scientifically called Onychoprion fuscatus, are omnivorous birds that eat reptiles, insects, earthworms, fish, amphibians, etc. They live in colonies with a unique migrating and attacking behavior. Their migration is in group with a given distance between every two sooty terns to avoid collisions. In a group, sooty terns can travel towards the direction of best survival fittest sooty tern. Based on the fittest sooty tern, other sooty terns can update their initial positions. The migration and attacking steps are considered as exploitation and exploration phases in a given search space. In this study, the STOA was developed to optimize the operation of Halilrood multi-reservoirs system, as demonstrated in Fig. 4.
MSA algorithm. Moth swarm algorithm (MSA), proposed by Mohamed et al. 22 , was inspired by the behavior of moths in the nature. The moths try to hide from predators during the day, while looking for food resources at night with a celestial navigation technique. They fly in a straight line over a long distance by steering their locomotion at a steady angle relative to moonlight as the celestial far-distant point light. In the MSA, the possible solution is represented by the position of a light source, and the quality of this solution is considered the luminescence intensity of the light source. Three groups of moths (pathfinder, prospectors, and onlookers) are considered in the MSA. Pathfinders are capable of finding the best position over the optimization space with First-In, Last-Out principle to guide the movement of the main swarm. Prospectors tend to wander into a ran-   www.nature.com/scientificreports/ dom spiral path nearby the light sources, which have been marked by the pathfinders. Onlookers drift directly toward the best global solution (moonlight), which has been achieved by prospectors' moths. In each iteration at MSA, each moth enters the problem to find the corresponding luminescence intensity of the light source. The best fitness in the population is considered as the position of pathfinder guiding for the next iteration. Thus, the second and third best groups are called prospectors and onlookers, respectively. The MSA algorithm is performed through three phases of initialization, reconnaissance, and transverse orientation 22 . Figure 5 shows phases of the MSA algorithm in the optimal operation of Halilrood multi-reservoirs system. PSO algorithm. The particle swarm optimization (PSO), proposed by Kennedy and Eberhart 23 , is a metaheuristic algorithm inspired by the swarm behavior of flocks of birds or schools of fish in nature. These swarms follow a cooperative way to find food, and each member in the swarms keeps changing the search pattern according to the learning experiences of its own and other members. While searching for food, the birds are either scattered or go together before they locate the place where they can find the food. While the birds are searching for food from one place to another, there is always a bird that can smell the food very well, that is, the bird is perceptible of the place where the food can be found, having the better food resource information. Because they are continuously exchanging information about the food place, the birds will eventually flock to the place where better food can be found.
In the PSO algorithm, solution swam is equal to the bird swarm, the birds' moving from one place to another is equal to the development of the solution swarm, good information is equal to the most optimist solution, and the food resource is equal to the most optimist solution during the whole course. The position of each particle in the swarm is affected both by the most optimist position during its movement and the position of the most optimist particle in its surrounding. In other word, the movement of each particle is identified in two phases of exploration (global search) and exploitation (local search). In the exploration phase, particle fly across the whole search space to find a limited region containing the global optimum. After the region containing the global optimum has been found, the exploitation phase is begun. The position of each particle is corrected by making small movements in the neighborhood of the global optimum. By adopting the correct sequence of these two phases, it is possible to lead particles towards the global optimum. Figure 6 shows the flowchart of the PSO algorithm in the optimal operation of Halilrood multi-reservoir system. GA algorithm. The Genetic Algorithm (GA), proposed by Holland 24 , is one of the most popular EAs that was inspired by Charles Darwin's theory of natural evolution. This algorithm reflects the process of natural selection where the fittest individuals are selected for reproduction in order to produce offspring of the next generation. GA algorithm starts with an initial set of random solutions, called population. Each solution in the population is called a chromosome, which represents a point in the search space. The chromosomes evolve through successive iterations, called generations. During each generation, the chromosomes are evaluated using some measures of fitness. The fitter the chromosomes, the higher the probabilities of being selected to perform the genetic operations, including crossover and mutation. In the crossover phase, the GA attempts to exchange portions of two parents, that is, two chromosomes in the population to generate an offspring. The crossover operation speeds up the process to reach better solutions. In the mutation phase, the mutation operation maintains the diversity in the population to avoid being trapped in a local optimum. A new generation is formed by selecting some parents and some offspring according to their fitness values, and by rejecting others to keep the population size constant.  www.nature.com/scientificreports/ After the predetermined number of generations is performed, the algorithm converges to the best chromosome, which hopefully represents the optimal solution or maybe a near-optimal solution of the problem. Figure 7 shows the flowchart of the GA algorithm in the optimal operation of Halilrood multi-reservoir system problem.
Simulation-optimization model. The conceptual scheme of the simulation-optimization model was demonstrated by Fig. 8. As seen, this is an integrated and comprehensive model because all the influencing parameters on the problem have been accounted. This optimization-simulation model was developed in a way that the best operational policies for Halilrood multi-reservoir system are achieved. The objective function of the optimization model was defined as the minimization of the total deficit over the operational period. The deficit means the difference between the amount of release from reservoirs and the amount of the downstream demands (drinking, agriculture, environment and industry demands). The greater the amount of deficit, then the greater the failure of the system. Release from each reservoir in each period was considered as a decision variable, and the storage capacity and the inflow to the reservoir were the state variables. The planning horizon was 223 months from October, 2000 to April, 2019. Therefore, the utilized algorithms have 669 variables (223 × 3) in the model of operation of Halilrood multi-reservoir system, making it a large-scale   (I) The constraint on the volume of overflow and evaporation losses: (II) mass balance between inflow and outflow regarding the reservoir storage capacity: (III) Constraints of decision variables S 3,t+1 = S 3,t + Q 3,t + Sp 1,t + Sp 2,t + Dez 1,t + Dez 2,t − Re 3,t − Sp 3,t − Loss 3,t where, Sp i,t is the overflow from the ith reservoir during period t (MCM), S i,t and S i,t+1 are the storages of the ith reservoir at the beginning and end of period t (MCM), S maxi is the maximum storage in the ith reservoir during period t, S mini is the minimum storage of the ith reservoir at the beginning of period t, Q i,t is the inflow into the ith reservoir during period t (MCM), Loss i,t is the evaporation loss from the ith reservoir surface during period t (MCM), A i,t is the area of the ith reservoir during period t (KM 2 ), Ev i,t is the net evaporation (evaporation minus precipitation) from the ith reservoir during the period t (m), a i , b i and c i are the coefficients of storage-area relation for reservoir i. Dez i,t is the environmental demand of ith reservoir during period t (MCM), Re maxi,t and Re mini,t are the maximum and minimum allowable release from the ith reservoir during period t, respectively.

Performance indicators.
Evaluating operation policies is the most important step in utilizing the optimization-simulation model. This is done by four indices of reliability, vulnerability, resiliency, and sustainability, which are formulated as follow: (I) Reliability Reliability index refers to the percentage of periods, during which the system has completely supplied the existing demands. Equation (11) shows the formulation of this index 25 .
where NDe f refers to the total number of failures occurred during the operation period, De t refers to the demand in period t, Re t is the release value of period t, Rel indicates the reliability of the system during the operation period. The higher the value of this index, the greater the reliability of the system.

(II)
Vulnerability This index shows the magnitude of system failure, which is calculated from Eq. (12) 25 .
where Val refers to system vulnerability, De t refers to demand in period t, Re t is the release value in period t, and T indicates the total number of operating periods.

(III)
Resiliency This index shows the system's capability to change the existing conditions. The likelihood that a system will return to normal after failure is called resiliency 25 . The value of this parameter can be calculated by the following equation: where () refers to the number of periods, in which the condition in parentheses has occurred, and Def t refers to water deficit in period t.

(IV)
Sustainability index This index is the summary of system performance criteria in an overall index to facilitate the comparison and decision-making between various scenarios based on reservoir performance indicators that are defined as follows 26 : Case study. Halilrood is one of the main sub-basins of the large Hamoun-Jazmourian basin, which is located in the southeast of Iran between longitude 56°-51′ to 61° 30′ east and latitude 26° 18′ to 29° 30′ north. Its total area is 7224 km 2 (Fig. 9). The range of elevation varies from more than 4000 m at the upstream to less than 500 m at the downstream. Halilrood is the major river in this basin which discharges to the Jazmorian Wetland at the end. Three large multi-purpose dams exist along this river tributaries, including (I) Baft earth dam with the purpose of water supply, (II) Jiroft concrete arch dam for water supply, electricity generation and flood control purposes, and (III) Safarood dam (under construction) for water supply purposes, which together act as a multi-reservoir system. Figure 10 demonstrates the schematic view of the Halilrood multi reservoir system.
The increasing water extraction from various uses and the lack of appropriate water allocation have imposed a significant impact on the downstream areas, especially on the Jazmorian Wetland, and caused serious conflicts in the region. If proper measures are not taken in the operation of existing reservoirs, the situation will get worse. One of the best strategies to overcome these issues is the optimization of the operation of this multi-reservoir if S i,t ≥ S mini and S i,t ≤ S maxi  www.nature.com/scientificreports/ system. In this study, five new robust evolutionary algorithms compared with two optimization algorithms of GA and PSO, were employed to derive the operation policies for Halilrood multi-reservoir system. The term "robust" refers to the algorithm with a set of features such as rapid convergence, time saving in long runs and uncertainty removal about the optimum solution.

Results and discussion
As previously described, to derive the long-term operation policies of Halilrood multi-purpose multi-reservoir dams (including Baft, Safaroud, and Jiroft dams), a series of recently introduced robust evolutionary algorithms in addition to the well-known GA and PSO algorithms were developed in this study. Here, the results of these algorithms were compared together, and the best one was selected based on the values of four performance metrics. A sensitivity analysis on the parameters of the utilized algorithms was done to obtain the best values for each parameter. Table 1 indicates the best values of parameters for Halilrood multi-reservoir system. Figure 11 demonstrates the convergence rate of applied algorithms in reaching the optimum value for the operation of Halilrood multi-reservoir system. It indicates the rapid convergence rate of the MSA in comparison with the HHO, SOA, STOA and TSA.  www.nature.com/scientificreports/ As seen, the highest convergence rate belongs to the MSA and HHO algorithms, respectively, and the slowest convergence rate is seen at the STOA and SOA. If the rate of convergence is higher, then typically fewer iterations are needed to yield an optimum solution, and therefore, the computational costs decreases. Among the utilized www.nature.com/scientificreports/ algorithms, the MSA has performed better in this regard; it could approach to the optimal value before 2000 iterations. The more popular algorithms of GA and PSO were placed in the middle ranks among all. The values of the objective function, the average CPU time, and the algorithms capability to meet the reservoirs downstream demands were presented in Table 2, based on PC with i7 CPU 1.8 GHz/16 GB RAM/2TB HDD. As shown, MSA could yield the best objective function (= 6.96) in the shortest CPU time (= 6738.56 s), while the corresponding values for the other utilized algorithms were not desirable. The MSA could supply 88.31%, 73.23% and 99.75% of downstream demands of Baft dam, Safarood dam and Jiroft dam, respectively, indicating the high performance of MSA in optimizing the studied multi-reservoir system. The next rank belongs to the HHO algorithm with the objective function of 14.68 and the CPU time of 10,223.4 s, which supplied 81.75%, 64.63%, and 98.26% of downstream demands of this multi-reservoir system with the same order. Table 3 shows the values of performance metrics of the reservoirs (reliability, resilience, vulnerability, and sustainability), resulting from the running of MSA, HHO, SOA, STOA, TSA, GA, and PSO algorithms in the long-term operation of the Halilrood multi-reservoir system.
As previously said, the sustainability criterion is an index for determining the performance of water resources systems. From Table 3 Figure 12a-c show the optimized average monthly water release pattern from the Baft, Safarood and Jiroft reservoirs simulated by the studied algorithms during the study period (2000-2019). As can be seen, the optimized values of release by the MSA algorithm have a better estimate of the downstream demands for all the studied reservoirs, and it resulted in the minimum deficit in the water supply. The MSA-based release policy in Baft, Safarood, and Jiroft reservoirs was able to supply the downstream demands, during the 19-years of operation, by the reliability of 74%, 49% and, 99%, respectively. Figure 13a-c show the total annual deficit in the Halilrood multi-reservoir system optimized by the studied algorithms. As can be seen, the MSA algorithm was well capable of minimizing the deficiencies of Baft, Safarood, and Jiroft reservoirs by 38, 151, and 3 million cubic meters, respectively, over 123 months of operation. The MSA-based deficit values had more appropriate estimates of the deficits for all the studied reservoirs, i.e., it demonstrated a high ability to supply demand and reduce the system's deficit. The results of Fig. 13a-c are consistent with those of Fig. 12a-c, in which, the MSA is the best and SOA is the worst algorithm in solving large-scale problems.

Conclusion
Developing optimal operation policies for single or multi-reservoir systems is a complex engineering problem. Solving this complex problem requires a robust, and high-performance method. In this study, four new-introduced evolutionary algorithms of HHO 19 , SOA 20 , STOA 21 , TSA 18 and MSA 22 were developed for the optimal operation of a multi-reservoir system. The performance of these algorithms in the optimization of Halilrood multi-reservoir operation was compared with each other, by four criteria of reliability, resilience, vulnerability, and sustainability. The results showed that the MSA with the objective function of 6.96 has a better performance compared to the HHO (14.68), SOA (31.28), STOA (29.55), TSA (25.32), GA (21.26), and PSO (16.04). Another aspect of the superiority of the MSA algorithm over the other four algorithms was the computational speed (CPU run time); it was 6739 s for MSA, 9347 s for PSO, 10,223 s for HHO, 11,782 s for GA, 12,971 s for SOA, 13,763 s for STOA, and 14,416 s for TSA. The demand deficit for the Baft, Safarood, and Jiroft reservoirs were equal to 37.66, 151.68, and 2.99 million cubic meters by the MSA algorithm, while, the corresponding values for the HHO algorithm, which had the second rank among the utilized algorithms, was 58.78, 200.46, and 20.52 million cubic meters, respectively. Finally, it was shown that the MSA-based optimized operating policy for the Halilrood multi-reservoir system, with sustainability criteria of 57.87 (Baft reservoir), 37.55 (Safarood reservoir), and 86.09 (Jiroft reservoir), demonstrated the highest performance among the utilized algorithms. Therefore, the use of MSA algorithm to extracting the operation policies of complex multi-reservoir systems is strongly recommended. In addition, the dynamic simulation of multi-reservoir systems by other robust optimization algorithms is recommended for future works.

Data availability
All data, models, or code generated or used during the study are available from the corresponding author by request.