A cuckoo bird-inspired meta-heuristic algorithm for optimal short-term hydrothermal generation cooperation

Optimal short-term hydrothermal generation cooperation (OSTHTGC) is used to determine the optimal operation strategy for cascaded hydropower plants and thermal plants. In this problem, equality constraints consisting of water time delay, reservoir volume, continuity water and power balance, as well as inequality constraints related to the limits of thermal and hydro generations are taken into account. In this paper, a Cuckoo Bird-Inspired Meta-Heuristic Algorithm (CBIA), which is advanced by using only few control parameters and high success rate, has been proposed to solve the OSTHTGC problem. The proposed method is applied for four test cases where the valve point loading effects of thermal units and the combination operation of a four-cascaded reservoir system are considered. The comparison of obtained results has revealed that that CBIA method is an efficient method with high quality solution and fast convergence speed compared to other relevant approaches. Subjects: Algorithms & Complexity; Artificial Intelligence; Power & Energy


Introduction
Optimal short-term hydrothermal generation cooperation (OSTHTGC) problem can be classified into two different problems, fixed head and variable head short-term hydrothermal scheduling. The former considers the water head of the reservoir as a constant while the latter treats the water head as Additional information is available at the end of the article ABOUT THE AUTHOR Thang Trung Nguyen received his BE and ME degrees in Electrical Engineering from University of Technical Education Ho Chi Minh City, Ho Chi Minh City, Vietnam in 2008 and 2010, respectively. Now, he is teaching at department of electrical and electronics engineering, Ton Duc Thang University, Ho Chi Minh city, Vietnam. His research interests include optimization of power system, power system operation and control and renewable energy. He has co-authored research papers and articles in national journals, international journals, conference proceedings, conference paper presentation and book chapters.

PUBLIC INTEREST STATEMENT
Cuckoo Bird-Inspired Meta-Heuristic Algorithm (CBIA) is an effective and robust technique for solving optimization problems where complicated objectives and complex constraints are considered. It comprised two phases of generating new solutions at each iteration, Lévy flights random walk to explore prospective domains of the search space and selective random walk to rapidly converge to a global optimal solution. Optimal short-term hydrothermal generation cooperation (OSTHTGC) minimizes the electricity generation fuel cost of thermal units while using the natural water of a river effectively with complicated hydro reservoir constraints. Applying CBIA approach to solve the OSTHTGC problem is essential and very promising. a variable. The head of the reservoir is fixed if the reservoir volume within the entire scheduling horizon is constant. This assumption is reasonable based on the fact that the hydropower plant has a large reservoir and the volume difference between inflow and discharge (via turbine) is very small. On the contrary, when the reservoir capacity is quite small or the difference between the inflow and discharge is high enough, the head of reservoir will be regarded as a variable. Obviously, the variable head short-term scheduling is more complex than the fixed head short-term scheduling because the power generation of hydro plants involves a function with respect to water discharge and reservoir volume, which are varying during the optimal operation process (Wood & Wollenberg, 1984). Furthermore, in this paper the hydro plants called cascaded reservoirs are considered where the discharge water of the upper reservoirs is the inflow water of the lower reservoirs.
Among these methods, Decomposition and Coordination techniques (SC) are the earliest methods employed to deal with the complex OSTHTGC problem. The methods have used Lagrange optimization function and divided a large problem into two sub-problems such as thermal sub-problem and thermal sub-problem (Soares et al., 1980;Wang & Shahidehpour, 1993). A principal difficulty which is considered over the optimal interval in the studies is the effect of stochastic load demand. Therefore, the stochastic sub-problem must first be solved to determine the fixed input data for the two other sub-problems. Then, based on the Lagrange function, the thermal sub-problem is solved for lamda value which is used as input data in the hydro sub-problem. Finally, the solution for the hydrothermal system scheduling including thermal and hydro generations is obtained. In spite of obtaining accurate solutions and insignificant constraint violations, the mentioned methods still suffer from the basic drawback of the Lagrange optimal approach, which is unable to deal with the complex power systems with the non-convex fuel cost function of thermal units. In Yang et al. (1996), the conventional evolutionary (CEP) has been applied to solve an OSTHTGC problem with the non-convex objective function. Other methods based evolutionary programming, such as the fast EP (FEP) and improved fast EP (IFEP) have been developed to solve the more complex OSTHTGC problem considering the non-convex objective function and prohibited operating zone of hydro units (Sinha et al., 2003). The CEP was stated more efficient and robust than simulated annealing via comparisons of result obtained from two different systems. However, there was no comparison among these modified EP methods with other research to verify the ability of the methods to deal with large scale and complex hydrothermal systems. The first classical GA (CGA) and its improved version, real-coded genetic algorithm (RCGA) applied for the OSTHTGC problem are presented in Orero and Irving (1998), Kumar and Naresh (2007), Basu (2010), respectively. The study did not aim to demonstrate any advantages of the CGA over other methods but testing the ability of CGA to deal with constraint violations on a power system of four cascaded hydropower plants and one thermal plant with quadratic fuel cost function of thermal units and without transmission losses. Dissimilarly, the RCGA has been tested on a larger power system of four cascaded hydropower plants and three thermal plants with non-convex fuel cost function. However, if only paying attention on the computational time, the GA based methods are relatively weak because they have taken much time to reach optimal solutions. Another method based TPNN taking the scheduled water discharge of hydro reservoirs as the states of analogue neurons (Naresh & Sharma, 1999) was developed to deal with the OSTHTGC problem and compared to the standard augmented Lagrange method (ALM). Although TPNN could obtain better solutions than ALM, but the method have the same disadvantage as applied to the non-differential problem. The modified DE (Lakshminarasimman & Subramanian, 2006) has proposed several modifications to deal with equality constraints like the load balance, especially the reservoir's volume at the final subinterval. Moreover, a new approach modified from hybrid DE (Lakshminarasimman & Subramanian, 2008) was developed to deal with equality constraints and to reduce computational time. The hybrid DE was built as two extra operations, acceleration and migration, where the former has improved the fitness quality and the later has exploited the search space. As a result, this method could converge to global optimal solutions very fast compared to conventional DE, MDE, HDE and several other methods. An adaptive chaotic differential evolution (ACDE) integrating an adaptive dynamic control mechanism for the crossover factor (Youlin et al., 2010) has been proposed to control the recombination and chaotic local search operation to avoid the effect of premature convergence (local trap). It is one of the best modified DE versions with the high quality solution and fast computational time. Similarly, the PSO method has been applied for a large scale hydrothermal system consisting of four hydro plants and three thermal plants with non-convex fuel cost function (Mandal et al., 2008). Several improved versions of PSO have also been developed by recommending new various factors, such as the inertia weight and constriction factor, and a selection strategy of best solutions from either all particles or just several specific particles. However, despite many attempts, the modified versions of PSO could not get the better solution than the improved versions of DE in the same test systems.
New modifications of DE, such as Clonal Selection Algorithm (Swain et al., 2011) and Hybrid Differential Evolution and Sequential Quadratic Programming (HDE-SQP) (Sivasubramani & Shanti Swarup, 2011), have been developed to improve the solution quality and convergence speed for hydrothermal systems with fixed head and variable head. Several test cases have been performed to verify the efficiency of the proposed methods in the complex power systems considering the nonconvex objective function and the prohibited zone of hydro units. An adaptive chaotic artificial bee colony (ACABC) (Xiang et al., 2013) has been implemented for searching the optimal solution of the short-term hydrothermal scheduling problem considering nonlinear constraints and nonlinear objective functions. The method could avoid the effect of premature convergence thank to the chaotic search and adaptive coordinating mechanism. A novel teaching learning based optimization (TLBO) has been applied to the problem with the non-convex fuel cost function and the prohibited zone of hydro units (Roy, 2013). The TLBO is mainly based on teaching phase and learning phase, and does not need any algorithm to determine the control parameters. Similarly, some recent state-of-the-art algorithms have been presented in Roy et al. (in press), Das and Bhattacharya (in press), Basu (2016), Dubey et al. (2016) to solve the OSTHTGC problem with promising results as the lower fuel costs and acceptable simulation time. However, the studies have verified their performance by the convergence time only, but there was no any comparison of computer processors which used in the performances. This issue will solved and presented in detail in our paper.
In this paper, a Cuckoo Bird-Inspired Meta-Heuristic Algorithm (CBIA), which is advanced by using only few control parameters and high success rates, has been proposed to solve the OSTHTGC problem. The CBIA, inspired from the intelligent reproduction behaviour of cuckoo birds (Yang & Deb, 2009), has several advantages over PSO and GA for benchmark functions such as obtaining better solution quality, higher success rate, and fewer easily selected control parameters. Consequently, CBIA has been widely and successfully applied for various engineering optimization problems such as the economic load dispatch (Dieu, Schegner, & Ongsakul, 2013), the hydrothermal scheduling Thang, Dieu, & Anh, 2015) and the distribution network reconfiguration. The OSTHTGC problem considered the non-convex fuel cost function of thermal units and the cascaded reservoirs of hydro plants. The four test systems have been applied to verify the performance of the proposed CBIA and the simulation results in terms of total costs as well as run time are used to compare to those of other reported methods. The comparisons have shown that the proposed CBIA is a very strong method for solving the short-term cascaded hydrothermal scheduling problem.

Fuel cost objective
The main objective of the OSTHTGC problem is to minimize total generation fuel cost while satisfying the constraints of hydraulic, load power balance, and generation capacity limits. The OSTHTGC problem having N 1 thermal units and N 2 hydro units scheduled in M sub-intervals is formulated as follows: where F im is the fuel cost of the ith thermal unit for one hour at the mth subinterval. Traditionally, the fuel cost of thermal units is approximately represented as a quadratic function, such as Recently, the fuel cost of thermal units with valve-point loading effects has been widely used in the power system optimization problems. This curve contains higher order nonlinearity and discontinuity due to the valve-point loading effects as follows:

Operational constraints
The proposed system has to satisfy following constraints: • Equality constraint of power balance The total power generation of all thermal and hydro units must satisfy the load demand and power losses in transmission lines as shown by where P L,m and P D,m are load demand and transmission loss at the mth subinterval. P hj,m is the power output of hydro plant j at the mth subinterval which is related to the water discharge and the reservoir volume as described by the function below: where C 1hj , C 2hj , C 3hj , C 4hj , C 5hj , C 6hj are the six coefficients of the jth hydropower plant.
• Hydraulic continuity equation The water flow of hydro plants must satisfy the hydraulic constraints as (1) where V j,m, I j,m and S j,m are reservoir volume, water inflow and spillage discharge rate of jth hydropower plant in mth interval. τ i,j is the water delay time between reservoir j and its up-stream i at the mth subinterval and N u is the number of up-stream units above hydro-plant j.
• Initial and final reservoir volumes • Reservoir storage and water discharge Limits where V j,max and V j,min are the maximum and minimum reservoir storage of the hydro plant j, respectively; Q j,max and Q j,min are the maximum and minimum water discharge of the hydro plant j.
• Generator operating limits where P si,max , P si,min and P hj,max , P hj,min are maximum, minimum power output of thermal plant i and hydro plant j, respectively.

The slack water discharge and the slack power output
In the optimization methods based meta-heuristic, to keep the balance constraints in power systems, most operational variables are first determined except the slack values. The slack water discharge of jth reservoir at the mth subinterval, Q j,M,d , and the slack power output of thermal unit 1 at the mth subinterval, P s1,m , the mth used to exactly keep power balance constraint (4) and volume constraint (7), respectively, are obtained by

Cuckoo Bird-Inspired Algorithm
Cuckoo bird is one of brood parasite species so it does not build its own nest and female cuckoo will lay her own eggs to other host bird nests. The cuckoos are very intelligent to choose the host bird, whose eggs with the same color as Cuckoos eggs. The strategy allows the Cuckoo egg to trick the host bird since the host bird cannot identify any alien eggs in their own nests. The fact demonstrates why there are more than 120 species of other birds can be cheated and continue to incubate the Cuckoo eggs until they are hatched.
Not every host bird is totally tricked, however, about 20% of Cuckoo eggs will be recognized as alien eggs and thrown away out of the nests or the host bird forsakes them and the host nest. Each female Cuckoo can lay between 12 and 22 eggs per season and lays each one in each nest. On the other hand, before laying Cuckoo eggs into other nests, the Cuckoos carefully observe the routine and the behavior of the other species to select the specie which has longer timing of hatching than them. Thanks to this strategy, Cuckoo chicks are very aggressive toward the host chicks; therefore, the first instinct action that Cuckoo chicks will do is to propel the host eggs out of the nest, increasing the food host bird provide the Cuckoo chicks.
Therefore, the Cuckoo optimisation algorithm (Yang & Deb, 2009), inspired from the cuckoo's behavior in the real life, consists of three idealized rules below: Rule 1: Each cuckoo lays eggs and put each egg in a nest of other species.
Rule 2: The best nest with the highest quality of cuckoo egg will be carried over to the next generation.
Rule 3: A fraction of the initial cuckoo eggs may be discovered as alien eggs by the host bird. The probability of the discovery is in range of [0, 1]. In this case, the host bird either propels the alien egg out of its nest or forsakes its nest to build a new one elsewhere. To simplify it is supposed that a fraction P a of the number of nests is replaced by new nests in this rule.
Therefore, the method based CBIA is developed corresponding to the three important stages below: • Initialization: A population of N p host nests is randomly initialized by using Rule 1.
• The first phase of generating new solutions: The new solutions will be randomly generated via Lévy Flights corresponding to Rule 2.
• The second phase of generating new solutions: The new solutions will be generated via the discovery of alien eggs corresponding to Rule 3.

Implementation of Cuckoo Bird-Inspired Algorithm
Based on the three rules in Section 3.1, the method based CBIA for solving OSTHTGC problems is described as follows.

Initialization
Similar to other meta-heuristic algorithms, each cuckoo nest (solution) in N p nests is represented by a vector X d = [P si,m,d, Q j,m,d ] (d = 1, …, N p ). Certainly, the valid solution has to satisfy the upper and lower limits shown in Equations (8)÷(11). Consequently, each solution X d (cuckoo nest) should be randomly initialized within the limits P si,min ≤ P si,m,d ≤ P si,max (i = 2, … , N 1 ; m = 1, …, M) and Q j,min ≤ Q j,m,d ≤ Q j,max (j = 2, …, N 2 ; m = 1, …, M − 1).
From Equation (6), the reservoir volume at the mth subinterval is obtained by: The values of Q j,M,d is obtained by (12) and hydro generations can be then calculated using Equation (5). The slack thermal unit is obtained using Equation (13).
Based on the initial population of nests, the fitness function to verify the quality of optimization solutions for the considered problem is calculated by where K s and K h are respectively penalty factors for the slack thermal unit 1 and all hydro units; K V and K Q are respectively penalty factors for reservoir volume over M − 1 subintervals and water discharge at the subinterval M. Futhermore, the limits of variables in Equation (15) can be calculated as below.

The first phase of generating new solutions
In this section, the new solutions generated by Lévy flights  as shown by where X best and X d are the best egg (the current best solution) and the dth egg in the population (N p eggs); α > 0 is an updated step size.
The value of α has a significant influence on the searching process because different new solutions would be produced according to specific values of α. If this parameter is set to a high value, there is a huge difference between the old and new solutions and the searching process is either obtained a chance of fast convergence or left outside of the feasible zone. On the contrary, if the value is set too small the location for the new solution is very close to the previous and the optimal search strategy is also not effective due to long computational time. Similar to other meta-heusistic methods, there is no criterion to appropriately select the parameter α, but it should be adjusted via experiment work. The new solutions generated from (20) must satisfy the required limits and need to be redefined if violated according to the equations below

The second phase of generating new solutions
In this section, the second phase of solution generation is used to improve quality of the previously obtained solution. This mechanism differs from other meta-heuristic methods, leading to a better optimization solution and faster computational time. In the cuckoo bird's behavior, there is a possibility that an alien egg may be identified by the host bird. Thus, it could either throw the exposed egg out of its nest or forsake its own nest with all eggs and build a new one. Like the Lévy flights, the discovery action of alien eggs in the nests is randomly implemented with a probability of P a to generate a new solution for an optimization problem. This new solution is created by: The newly obtained solutions also need to be redefined using Equation (21) above in case they violate the upper or lower values.

Stopping criteria
The searching process mentioned above will stop when the maximum number of iterations is reached.

Overall procedure of solving the OSTHTGC problem using the proposed CBIA
Step 1: Select CBIA parameters including the number of host nests, N p , the probability of a host bird discovering an alien egg in its nest, P a , and the maximum number of iterations, G max .
Step 2: Initialize the population of N p host nests as mentioned in Section 3.3.1.
• Calculate the slack water discharge and the slack thermal unit 1 using Equations (12) and (13).
Step 3: Evaluate the fitness functions using Equation (15) to choose the best nest (current best solution) corresponding to the lowest fitness function value, X best .
• Set the initial iteration G = 1.
Step 4: Generate a new set of solutions via Lévy flights as described in Section 3.3.2 and repair any violated solutions by using Equation (21).
Step 5: Calculate the slack water discharge and the slack thermal unit 1 corresponding to the new set of solutions using Equations (12) and (13), respectively.
• Calculate new values of all hydro generations using e.g. (5).
Step 6: Calculate the functions of the newly obtained solutions using Equation (15).
• Compare the fitness of the new solution and that of the old solution (at the same nest) to retain a better one at each nest.
Step 7: Generate a new set of solutions based on the discovery action of alien egg as described in Section 3.3.3 and repair any violated solution using Equation (21).
Step 8: Calculate the slack water discharge and the slack thermal unit 1 corresponding to the new set of solutions using (12) and (13), respectively.
• Calculate new values of all hydro generations using (5).
Step 9: Calculate the fitness functions for the newer set of solutions using (15).
• Compare the fitness of the new solution and that of the current solution (at the same nest) to retain a better one at each nest.
Step 10: Evaluate all updated current solutions to choose the best current solution, G best .
Step 11: Check if G < G max , G = G + 1 and return to Step 4. Otherwise, stop.

Numerical results
In this paper, the performance of the proposed CBIA has been verified by using four test cases of the OSTHTGC problem. The proposed CBIA is coded in MATLAB platform and implemented 50 independent trials for each specific value of P a which run on the PC configuration of 1.8 GHz processor and 4 GB RAM.

Two test cases with the quadratic fuel cost function of thermal plants
In this section, two test cases comprising four cascaded hydropower plants and one thermal plant with the quadratic fuel cost function scheduled in 24 one-hour sub-intervals are considered. The transportation delay times in hour considered in test case 1 are τ 13 = 2, τ 23 = 3, τ 34 = 4, and in test case 2 are τ 13 = 1, τ 23 = 2, τ 34 = 2. The full data for the test case 1 and test case 2 are taken from Lakshminarasimman and Subramanian (2008), Youlin et al. (2010), respectively.
The simulation has been performed via various values of P a in a range of [0.1, 0.9] with a step of 0.1 to analyze the effect of P a for the CBIA method. The number of nests and the maximum number of iterations are respectively set to 100 and 15,000. Tables 1 and 2 present the statistic data of obtained results for the test case 1 and 2, respectively. As observed from Table 1, the best minimum and average costs are respectively obtained at P a = 0.3 and 0.1 whereas the best standard deviation costs are gotten at P a = 0.4. Table 2 shows that the best minimum cost and average cost are respectively $154590.1657 and $154595.566 obtained at the same value P a = 0.6 meanwhile the standard deviation costs is reached at P a = 0.3 for system 2. Moreover, Tables 3 and 4 show the comparison results of CBIA and other methods for test cases 1 and 2, respectively. As shown in these tables, CBIA can obtain the optimal solution better than all other relevant methods (Kumar & Naresh, 2007) and the computation time faster than most of mentioned methods except BCGA and and RCGA for only test case 1. To analyze in more details about the CBIA performance, the fitness convergence characteristic for the test system 1 is depicted in Figure 1 as well as the optimal solutions in terms of water discharge, thermal and hydro generations for the two test systems are shown in Figures 2-5. All of these results demonstrate the advanced efficiency of CBIA compared to other relevant methods. Furthermore, to fairly compare the computational time among the considered methods for solving the same problem but using different programming languages and/or different computer processors, we convert the execution time corresponding to CPU speed (GHz) of the mentioned methods into a common base unit which is defined by CPU time of CBIA with 1.8 GHz processor. Thus, the adjusted CPU time in PU is determined by: It is noted that the value of 1.8 (GHZ) is the processor of the CPU chip used to run CBIA and the CPU time obtained by CBIA is used as the common base time. Therefore, the adjusted CPU time determined for other methods is a number of relative times compared to CPU time of the proposed CBIA. Of course, according to Equation (23), the adjusted CPU time of the CBIA is absolutely 1.0 pu.
Adjusted CPU time = Given CPU speed (GHz) 1.8 (GHz) × Given CPU time (second) CPU time from CBIA (second)     Consequently, if the adjusted CPU time of any method is 2.0 pu, it means that this method executes two times slower than the CBIA in the same problem. On the contrary, if that value is only 0.5 pu, it indicates that the method runs two times faster than CBIA. Table 5 presents the adjusted CPU times in pu for the test case 1 corresponding to the CPU chip used to run these methods. However, it is unable to do the same for the test system 2 because the relevant papers of mentioned methods did not report their execution time. As observed from Table 5, the adjusted CPU time for test system 1, the proposed CBIA method is the relatively fastest, which is from 1.17 to 13.77 times faster than other mentioned methods.

Two test cases of complex power systems with the non-convex fuel cost function of thermal plants
In this section, the test case 3 has applied for a power system comprising four cascaded hydropower plants and one thermal plant with the non-convex fuel cost function and test case 4 has comprised four cascaded hydropower plants and three thermal plants with the non-convex fuel cost function. The optimization period is scheduled in 24 one-hour sub-intervals. The full data for the test case 3 and test case 4 are taken from Kumar & Naresh (2007) and Basu (2004), respectively. Similar to the investigation of the CBIA method for the test systems with the quadratic fuel cost function of thermal plants, as desribed in Section 4.1, this simulation has been performed via 9 values of P a in a range of [0.1, 0.9] with a step of 0.1. The number of nests and the maximum number of iterations are respectively set to 100 and 15,000. Tables 6 and 7 present the statistic data of obtained results for the test case 3 and 4, respectively. As observed from these tables, the best minimum cost of test systems 3 is $947056.12 at P a = 0.4 and that of test system 4 is $41064.897 at P a = 0.5. To analyze in more details about the CBIA performance, the fitness convergence characteristic for test system 4 is depicted in Figure 6 as well as the optimal solutions in terms of water discharge, thermal and hydro generations for the two non-convex test systems are shown in Figures 7-10. Moreover , Tables 8 and  9 show the comparison results of CBIA and other relevant methods for test case 3 and 4, respectively. As shown in Table 8, the CBIA method can obtain the optimal solution better than other methods but taking more time for convergence. As seen from the Table 9, the fuel cost of CBIA is the  smallest one while the second best fuel cost from MHDE (Lakshminarasimman & Subramanian, 2008) is $41,856.50 and the poorest from SA (Mandal et al., 2008) is $47,306.00. It implies that the saving cost, thanks to using CBIA, is approximately 2 and 15% compared to using SA and MHDE, respectively. Table 10 reports the adjusted CPU time in pu of relevant methods to fairly compare among their execution time. As seen from test system 3, the CBIA and two versions of GA including BCGA and RCGA (Kumar & Naresh, 2007) have nearly the same values of the adjusted CPU time. On the contrary, at test system 4 the adjusted CPU time of CBIA is higher than that from most other Figure 7. Optimal water discharge for test case 3 (non-convex system).

Figure 8. Hydro and thermal generations for test case 3 (non-convex system).
methods except DE and MDE (Lakshminarasimman & Subramanian, 2008) and Clonal selection algorithm (Swain et al., 2011). However, the mentioned methods with faster compuation time have implemented in the C++ and JAVA programming which are from two to three times faster than Matlab.