A Phase-Based Adaptive Differential Evolution Algorithm for the Economic Load Dispatch Considering Valve-Point Effects and Transmission Losses

A phase-based adaptive differential evolution (PADE) algorithm is proposed to solve the economic load dispatch (ELD) considering valve-point effects (VPE) and transmission losses. To a great extent, PADEmakes up for the drawbacks of the traditional differential evolution (DE) through three improvements. First, we establish an archive of storing successful individuals to improve the quality of offspring. Second, to balance the exploring and exploiting ability of the algorithm, a phase-basedmutation operation is carried out. Third, two control parameters are adaptively adjusted, which is helpful for enhancing the robustness of the algorithm. In addition, two types of repair methods of constraint handling are employed for the ELD without or with transmission losses to help PADE find feasible solutions more efficiently. A performance comparison between PADE and other DE approaches from the literaturewas carried out on six ELD test cases which consider a set of operating constraints including the VPE and transmission losses. Results show a competitive PADE performance in all test cases regarding the compared DE approaches. Compared to methods from the literature, the costs obtained by PADE are lower in most cases while the corresponding constraint violations reach a lower level.


Introduction
In power systems, economic load dispatch (ELD) [1] is an important, complex optimization problem that allocates the required generation among the available generators. The main purpose of ELD is to find the optimal distribution of the generating units to make the fuel cost minimize while meeting all inequality and equality constraints in practice. Therefore, an accurate mathematical model is needed for adapting the ELD problems with different characteristics. So far, the quadratic function has been widely used as the classical mathematical model for solving the cost of the ELD problems. Moreover, addition conditions are often considered in the actual system, which constitute various versions for the ELD problems such as the valve-point effects (VPE) [2], transmission losses [3], and prohibited operating zones (POZs) [4].
Many classical methods, such as linear programming [5], quadratic programming [6], lambda iteration method [7], and gradient method [1], have been applied to solve the conventional ELD problems with the quadratic function. The above classical methods require that the generators exhibit convex and smooth characteristics. However, the inputoutput characteristics of generators are generally nonconvex and non-smooth in practical thermal generation plants. The reasons are caused by the VPE, POZs and so on, which may make the cost function produce a number of local optima. Therefore, these classical methods can hardly find the global optima for such types of problems. Dynamic programming [8] is different from the above classical methods. It can handle nonconvex and nonsmooth ELD problems, but it is restricted by the curse of dimensionality that inevitably leads to extra computational cost. Thus it is not also a suitable one for practical ELD problems.
In order to seek the optimal solution of ELD, more and more scholars tend to use modern intelligent algorithms. These modern intelligent algorithms have huge advantage over classical methods, especially for the ELD problems with various operation constraints. They are essentially approximation methods, which do not need the derivative 2 Mathematical Problems in Engineering information of the ELD problems. The practical ELD problems are quite complex, so most of scholars are committed to continuously improving or hybrid modern intelligent algorithms and developing new constraint handling mechanisms at present. Qin and Cheng et al. [9] adopted an orthogonal designed method and proposed auxiliary vector generation based on multiple strategies to enhance the effectiveness of orthogonal designed operation. Based on the adaptive adjustment of acceleration coefficients, the algorithm's robustness and global search capability can be improved by employing a tent chaotic map. Furthermore, a repair method is designed to handle the practical constraints. Coelho and Mariani [10] put forward an improved harmony search (IHS) algorithm based on exponential distribution for the ELD problems. In IHS, a repair process is adopted to find feasible solutions. Al-Betar and Awadallah et al. [11] developed a modified version of harmony search algorithm that is called tournament-based harmony search (THS) algorithm, where the tournament selection process replaces the random selection process in the memory consideration operator. The tournament selection process is beneficial to improve the convergence performance of the algorithm. The introduced repair strategy can completely guarantee the feasibility of each generated solution before calculating the fuel cost. Therefore, to some extent, the THS algorithm can achieve better solutions compared to many state-of-the-art optimization approaches. Mandal and Roy et al. [12] proposed a new krill herd algorithm (KHA) to solve the ELD problems with various constraints. The crossover and mutation operations of DE are embedded into the KHA algorithm. The four versions of KHA are presented to solve different scale ELD problems. Zou and Li et al. [13] proposed a new global particle swarm optimization (NGPSO) algorithm for the economic emission dispatch with or without transmission losses. First, NGPSO devises a new position updating equation that utilizes the global best particle to guide the convergence of all particles. Second, the randomization based on the uniform distribution is used to avoid the premature convergence in the late optimization process. Simultaneously, different treatment methods are adopted for the equality and inequality constraints, respectively. For the inequality constraints, a common and simple penalty function is used. More importantly, the equality constraints are handled by three important steps, where the equality constraints are completely desirable. Different from the reported variable selection [14,15], a variable is randomly selected to treat all variables equally in the repair process. As a result, a large number of high-quality feasible solutions can be obtained. Kumar and Chaturvedi [16] solved the ELD problems using a fuzzy particle swarm optimization (FPSO). The quadratic cost function with valve point loading and prohibited operating zone is studied in the models. The suitable value of inertia weight is yielded by a fuzzy logic controller (FLC). The FLC is composed of a knowledge base, a fuzzification interface, an inference system, and the defuzzification interface. As a result, FPSO efficiently avoids premature convergence. Mahdad and Srairi [17] adopted a hybrid method based on genetic algorithm (GA), differential evolution (DE), and pattern search (PS) [18], namely, GA-DE-PS. First, GA and DE are jointly used to find the optimal solutions. Second, the optimal solutions achieved by GA and DE are adjusted by PS. Third, the new solutions obtained by PS are employed to GA and DE for the new initial solutions of GA and DE. The interactive mechanism allows individuals to communicate with each other, which is very helpful to balance the exploiting and exploring capability of GA-DE-PS. Jeddi and Vahidinasab [19] proposed a modified harmony search algorithm (MHSA). In MHSA, a new improvising method based on wavelet mutation combined with a new memory consideration scheme based on the roulette wheel mechanism is adopted, which improves the accuracy, accelerates convergence speed, and enhances robustness of the classical HSA. For the constraint handling strategy, a common penalty function method is used to find feasible solutions. Singh and Dhillon [20] proposed an oppositionbased greedy heuristic search (OGHS) to solve practical optimization problems. The initial population is composed of several randomly generated candidate individuals based on the uniform distribution, and its opposite population is formed by the opposite-based learning strategy [21]. Then the better half of the two populations is adopted to carry out subsequent optimization process, which improves the convergence speed of the algorithm. The opposition-based learning is also applied to its migration to maintain the diversity. The mutation strategy applied is based on the greediness and randomness to find the global best solution.
According to the review above, most works are focused on enhancing or hybridizing intelligent algorithms, and some of them studied the constraint handling improvement for practical ELD problems. The inequality constraints of ELD problems are relatively easy to meet. In contrast, the equality constraints of the ELD problems are hard to handle. Thus they mainly study the treatment method of the equality constraints for the ELD problems. However, their methods of treating the equality constraint with transmission losses are usually the same as those without transmission losses. In fact, the complexity of both equality constraints without transmission losses and equality constraints with transmission losses has a great difference. Therefore, the method of handling the equality constraint with transmission losses should be different from that of handling the equality constraint without transmission losses, which can differently repair the infeasible solutions. For our studies, the equality constraint refers to the power balance constraint. Meanwhile, in order to further improve the performance of the original DE algorithm and repair the infeasible solutions, the contributions of this paper contain two aspects. First, a phasebased adaptive differential evolution algorithm (PADE) is proposed to improve the performance of the original DE algorithm. Second, two repair methods are used to properly handle the equality constraints of two types of ELD problems, respectively.
The remainder of the paper is organized as follows: Section 2 describes the mathematical model of the ELD problems; Section 3 introduces the traditional differential evolution algorithm and the proposed DE variant; simulation results and the corresponding comprehensive analysis are presented in Section 4; Section 5 concludes the paper.

Economic Load Dispatch (ELD)
2.1. Objective Function. The objective of ELD is to minimize the total fuel cost of the entire power system under the condition of meeting various constraints. The cost of ELD can be described as follows: where represents the total fuel cost of the power generating (in $/h). $/h represents the total fuel cost of the power generating per hour. represents the fuel cost of the generating unit (in $/h). The variable represents the real output power of the generating unit (in MW). MW means megawatt.
is the number of units in the power system. In general, the traditional cost function of ELD is represented as the quadratic function, which can be given by Coelho and Mariani [10]: where , , and denote the cost coefficients of unit . Under most circumstances, the valve-point effects (VPE) are contained in the practical power system. In view of this condition, the sinusoidal component is combined with the quadratic function, and the practical cost function can be amended as follows [10]: where and denote the cost coefficients of unit considering VPE.

Constraints.
The ELD problems are solved subject to some inequality and equality constraints, including power output limits, power balance, and prohibited operating zones (POZs).

Power Output Limits.
The power output of each unit should be between its upper and lower bounds [10] in the power system.
where min (in MW) and max (in MW) stand for the minimum and maximum power output of unit , respectively.

Power Balance Constraints.
The total generated power should be equal to the sum of the total load demand and transmission losses [10].
where denotes the total load demand of power system (in MW); denotes the total transmission losses (in MW), and it can be commonly obtained by Kron's loss formula as follows [13]: where ℎ denotes the transmission loss coefficient. Additionally, the simplified version of power balance constraint has been commonly studied by scholars. This version does not include the transmission loss , namely, = 0. Furthermore, it can be expressed by

Prohibited Operating Zones (POZs).
Under certain circumstances, a generating unit may not operate within certain ranges known as prohibited operating zones due to physical operation limitations. As a result, the feasible operation zones of the generating unit may be discontinuous. At this time, the feasible power output ranges of the generating unit can be properly described as follows [4]: where , (in MW) and , (in MW) denote the lower and upper bounds of prohibited operating zone of unit , respectively. is the number of prohibited operating zones of unit .

Constraint Handling Mechanism
2.3.1. Inequality Constraint Handling. The inequality constraints (4) are easy to meet, because the optimization variables are directly set within their ranges in the initialization process, and any variable beyond their bound can be limited to the bound which it exceeds for the proposed algorithm. For the inequality constraints (8), a common penalty function approach is used to handle it. If the output power of the generating unit lies in the POZs, the new cost function of the ELD problems with the POZs is expressed by where is a penalty factor. min( − , −1 , , −1 − ) denotes constraint violations when the unit operates in the POZs.

Equality Constraint Handling for the ELD Problems without Transmission
Losses. Contrary to inequality constraints, the equality constraints are hard to satisfy. The penalty function method can be also used to help intelligent algorithms to find the optimization solutions, but the quality of solutions cannot be ensured. Aim at the equality constraints, scholars proposed diverse repair methods for the ELD problems with or without transmission losses [9,10,13]. Due to the difference between the ELD problems without transmission losses and the ELD problems with transmission losses, two proposed repair methods will be, respectively, introduced. A repair method of this paper for the ELD problems without transmission losses is first introduced in this section.
To obtain more promising feasible solutions, a repair method and penalty function method are used to together cope with the equality constraint, and the equality constraint handling process is clearly shown in Algorithm 1.
Δ represents the adjustment limit of the component . | | represents the equality constraint violation. sign is called sign function. When is positive, the sign( ) value is equal to 1. On the contrary, when is negative, the sign( ) value is equal to -1. The repair method will select a random component from the individual vector ( = [ 1 , 2 , ⋅ ⋅ ⋅ , ], = 1, 2, ⋅ ⋅ ⋅ ), and let = {1, 2, ⋅ ⋅ ⋅ , }. If the selected component is equal to the upper bound max (or lower bound min ) when the total generating power output is smaller (larger) than demand, this component will be removed from the set . Then a new component will be randomly selected from the set that has excluded the previous component (s). This process will be repeated unless a component that is not equal to its upper (or lower) bound is found. The way of selecting a component can avoid some ineffective operations which do not repair those infeasible solutions. The smaller one in | | and Δ is selected to adjust the output power of the unit, which is to avoid the violation of inequality constraints after executing the repair process. Similar to the reported literature [10], a penalty method is still used to further handle equality constraint. Furthermore, a sign flag is set to judge whether the equality constraint has been completely met after the repair operation. If the equality constraint violations are surely equal to zero, the penalty method will not be used to further handle the equality constraint. On the contrary, if the equality constraint violations are not equal to zero, the penalty method will be used to further handle the remaining constraint violations. The inequality constraints related to the POZs are also handled by the penalty method. Because the environment of all experiments is based on MATLAB, MATLAB precision has a certain effect on the optimization results. The advantage of the above setting flag is that it can decrease the effects of MATLAB precision on the optimization results as much as possible. If another programming language is used to implement the proposed mechanism, the flag variable can be also used to eliminate the error due to the other operating environments precision. Thus, the flag variable cannot be ignored when other programming languages are used.

Equality Constraint Handling for the ELD Problems
with Transmission Losses. Qin and Cheng et al. [9] presented a repair method that can reduce the equality constraint violations until these violations are within the error tolerance ( = 0.001), which means that the repaired solutions are considered to be the feasible ones. A component is randomly selected from the set . , , , | |, and sign represent the concept same as the aforementioned repair method. rand represents a randomly generated number in the interval [0, 1]. The method also compares the | | value with the Δ value to select a smaller one that is used to update solutions. If the | | value is larger (smaller) than the Δ value, will be subtracted to sign(w) × Δ (or ). Then the transmission losses and equality constraint violations are recalculated. This above procedure is repeated until the | | value is smaller than the error tolerance . The proposed repair method is similar to the repair method of Qin and Cheng et al. [9]. The only difference is that Δ is calculated based on the adjustment limit on the component , and the repair method for the ELD problems with transmission losses is shown in Algorithm 2. This modification contributes to reduce more constraint violations and obtains the better solutions satisfying the tolerance error compared to the original mechanism. It will be indicated from Section 4.8 that the modification is effective.

Differential Evolution Algorithm and Its Improved Version
3.1. Differential Evolution Algorithm. Differential evolution (DE) algorithm, presented by Storn and Price [22], is a very competitive intelligent algorithm. It has simple structure and powerful capacity of searching solution space. Therefore, DE has been frequently applied to many fields for solving those real-word problems, such as optimal reactive power dispatch [23], parameter estimation problems in computational systems biology [24], iris recognition [25], and constraint optimization [26]. In this section, the traditional DE algorithm and its improved version will be introduced in detail.
The operation steps of DE are initialization, mutation, crossover and selection, respectively. In addition, a few parameters also play a key role on the DE's performance. These steps and parameters are explained as follows.
3.1.1. Initialization. The goal of this step is to generate candidate solutions consisting of a population. All candidate solutions are randomly generated within the lower and upper bounds of searching space. For example, the th ( = 1, 2, ⋅ ⋅ ⋅ , ) component of the th ( = 1, 2, ⋅ ⋅ ⋅ , ) candidate solution is initialized by where and denote the lower and upper bounds of the th component, respectively. (0) denotes the th variable of the th individual at the initial generation. (2) If < 0 (3) Randomly select a component from the set ; Remove from , and let the new set be . Then randomly select a new component from , and ← , ← ; End While (7) Calculate the adjustment limit of the component , namely Δ = max − ; Randomly select a component from the set ; (10) While = min (11) Remove from , and let the new set be . Then randomly select a new component from , and ← , ← ; End While Randomly select a component from the set ; End If (14) Recalculate the transmission losses and equality constraint violations | |; (15) End While; Algorithm 2: The pseudo code of the equality constraint handling process for the ELD problems with transmission losses. 6 Mathematical Problems in Engineering formula. More than common and popular, and it is the standard or canonical mutation operator in DE.
where 1, 2, and 3 are three randomly selected integers in [1, ], and 1 ̸ = 2 ̸ = 3 ̸ = . The control parameter represents the scale factor that is used to scale the difference vector ( 1 ( ) − 2 ( )). The scale factor is a user-defined parameter ( > 0) [27] and usually set in the interval (0, 1]. Except for this mutation strategy DE/rand/1, another five popular mutation strategies are listed below [28]: D DE/current-to-best/1 E DE/rand-to-best/1 where ( ) represents the best individual vector at generation .

Crossover Operation.
Crossover has two types of strategies including binary crossover and exponential crossover. In this work, the binary crossover is adopted to obtain the trail vector ( ) ( ( ) = ( 1 ( ), 2 ( ), ⋅ ⋅ ⋅ , ( ))) whose components are got from those of both parent vector ( ) and mutant vector V ( ) based on the crossover rate. This operation is expressed by where another important control parameter stands for the crossover rate which determines the probability of selecting components of mutant vector V ( ) ( = 1, 2, ⋅ ⋅ ⋅ , ). is a user-defined parameter set in the interval [0, 1].
is a randomly generated integer between 1 and .

Selection
Operation. The behavior of this operation is similar to a greedy strategy, and it only reserves better individuals into the next generation. For DE, if the updated trail vector ( ) is better than parent vector ( ), the former will replace latter into the next generation. Otherwise, the corresponding offspring individual vector is the same as the parent individual vector, which means that the updating of the individual is not successful. For a minimization problem, the offspring individual ( + 1) ( = 1, 2, ⋅ ⋅ ⋅ , ) is generated by where (⋅) stands for the fitness function, which is to calculate the fitness value of an individual. In this problem, the smaller the fitness value of individual is, the better the performance of individual is.

A Phase-Based Adaptive Differential Evolution Algorithm.
A phase-based adaptive differential evolution algorithm (PADE) is proposed to enhance the DE's ability in finding the global best solution. PADE and DE are different in three aspects.

An Archive of Storing Successful Individuals.
In the DE's mutation strategy, the parent individual is usually randomly selected from the current population. This parent-selecting way is beneficial for maintaining the population diversity, but it may not jump out of local optima when stagnation is happening to DE. Based on this issue, Guo and Yang et al. [29] established an archive of recently updated individuals to improve the parent-selecting way by utilizing the successful individuals. If the current parent individual fails to update beyond the stagnation tolerance on the number of consecutive unsuccessfully updates, the parent individuals of mutation strategy will be randomly selected from the archive instead of the current population. When the archive exceeds its size, the oldest individual is removed. On the contrary, Tanabe and Fukunaga [30] adopted an external archive that stores the failed solutions at each generation to maintain the population diversity. Inspired by the study, we also set up an archive that stores the recently updated individuals. Different from the archive of [29,30], the randomly selected individuals from the archive are used for the base vector of mutation strategy in the current generation. Then the recently updated individuals are stored to that has been set to empty after each generation. In other words, only stores the recently updated individuals. As long as the recently updated individual exists, the base vector of each generation is randomly chosen from the archive in the early phase. In addition, if there are not updated individuals in a generation, the base vector of mutation strategy will be randomly selected from the current population. The new parent-selecting way has a greater chance to generate successful individuals, and maintains the population diversity. For the initialization of , it is first initiated to be a copy of the initial population. Additionally, in case the initial population contains infeasible solutions, they will be also used for the mutation operation in the early evolutionary phase, because those infeasible solutions are gradually evolved to feasible ones as the evolution progresses.

A Phase-Based Mutation Operation.
In DE, mutation is an important step for obtaining the mutant individuals in the area around the base individual. In this area, the base individual locates in the center, the difference vector acts as the movement direction and step, and the scale factor scales the movement step. In fact, there are many different types of the base individuals and difference vectors, which consist of several mutation strategies of different characteristics, such as DE/rand/1, DE/rand/2, DE/best/1, DE/best/2, and DE/current-to-best/1. According to [31], DE/rand /1 and DE/rand/2 are very helpful to maintain the diversity of population but they perform poor in convergence rate and are easy to gradually approach the local optima. On the contrary, it is known from [31] that DE/best/1 and DE/best/2 are beneficial to accelerate population convergence, because the base individual is the global best individual with good information, which can guide other individuals to approach it. However, it is obvious that they are very easy to fall into the local optima. Specially, it can lead to premature convergence when they are used in the early phase of evolution. DE/current-to-best/1 [32] employs the current individual as the origin individual and the best individual as the guiding individual, which constructs a mixed base individual. DE/current-to-best/1 is called composite mutation strategies, which has better diversity than DE/best/1, but its global exploration is poor. In the entire evolutionary process, global exploration and local exploitation are two vital aspects, which should cooperate well with each other to find the promising solutions. First, the algorithm should widely search in the global area and find some local areas with promising individuals in the search space. Then, the algorithm should focus on local exploitation in these local areas. It is to enhance the ability of finding the optimal solution as fast as possible. In DE, the population individuals spread over the searching space in the early phase of evolution. They can update themselves based on other individuals through global search but it can also cause premature convergence. In the late phase of evolution, most individuals gather in a relatively small area, and they implement local exploitation instead of global exploration.
Based on the above analysis, Miranda-Varela and Mezura-Montes [33] relied on the number of feasible solutions to select one between DE/rand/1 and DE/best/1 to balance the exploration and exploitation ability. Qin and Huang et al. [31] constructed a mutation strategy candidate pool that contains several mutation strategies with diverse characteristics. In the evolution process, one strategy will be selected from the candidate pool according to a probability learned from its previous experience of generating promising solutions for each parent individual. The more successfully one strategy behaved in previous generations, the more probably it will be selected in the current generation. Different from the above studies, to match different phases of evolution, DE/srand/1 and DE/ best/1 are respectively used to execute the mutation operation in the early and late phase of evolution in PADE. For DE/srand/1, its base vector is randomly selected from the archive , and its difference vector is the same as DE/rand/1. Therefore, DE/srand/1 is more easy to find the successful individuals than DE/rand/1. Now, the key task is how to select a generation index threshold between 1 and to divide the entire evolution process into the early and late phase. Many scholars tend to set a probability rule based on the number of iterations [34,35] or the fitness value of random individual [36] to determine whether the current evolution lies in the early phase or late phase of iteration. Different from their methods, Tang and Dong et al. [37] introduced the success ratio and an auxiliary generation index threshold to help select the generation index threshold . The threshold is set based on the recently updated individuals, which is more heuristic for enhancing the performance of DE than the above two methods. In view of the advantage of the success ratio, this paper employs it to divide the evolution process into the early and late phase. The success ratio of each generation is defined by: where is the number of the updated individuals at generation .
stands for the success ratio at generation . In general, is a relative large value in the early phase of evolution. As the generation advances, is prone to decrease continually, so it can be used to determine whether the late phase begins. When is smaller than a threshold of success ratio for a given number of consecutive generations max , the late phase of the algorithm begins. If is larger than at generation , the current number of consecutive generations is denoted as = 0. Otherwise, = + 1. If reaches max , keeps the fixed value max until the evolution ends. That means that PADE only uses the mutation strategy DE/best/1 in the late phase. For most problems, can decrease to for consecutive generations max in the early phase of evolution. However, it should be noted that cannot decrease to for consecutive generations max in the early or entire evolution phase for a few problems. Therefore, under these circumstances, a threshold is forced to set to divide the process into two phases. The threshold determined by success ratio is given by In order to comprehensively describe the phase-based mutation operation producing the mutant vector, the following formula is given by 8 Mathematical Problems in Engineering where the index of successful individual is randomly integer between 1 and the size of the archive , and To visualize this phase-based mutation operation from an overall perspective, its schematic diagram is presented as follows (see Figure 1).
In Figure 1, the entire evolution phase is divided into the early and late phase. There are two situations on determining the threshold, which are represented as S1 and S2, respectively. Specifically, above the 0 to axis, the way of dividing the entire evolution phase is represented as S1. Below the 0 to axis, the way of dividing the entire evolution phase is represented as S2. The threshold can lie in [ , ], or it does not exist in the [0, ], where the is used as the threshold at this time just as S2.

Adaptive Control Parameters.
The two control parameters, scale factor and crossover rate , play a significant role on the performance of DE. Their values are traditionally fixed as constants to participate in evolution. A suitable and are related to a specific problem, and all individuals share their fixed control parameters in the evolution process for the original DE. The result is that the robustness of the algorithm is poor, and a lot of tedious optimization trials must be done to find a suitable and for practical problems. To enhance the robustness of the algorithm and avoid tedious optimization trials, many researchers adopted different selfadaptive strategies to make control parameters self-adaptively change with individual evolution [36,38]. Some of them randomly take a value within the range of parameters owned by the individual when the individual doesn't survive into the next generation. The other do not have the guidance rule at all, and the control parameters are changed based on a probability in the search process. In addition, there are some works [39,40] only for this type of problems, where selfadaptive approaches have proven to be competent in different problems, for example, constraint optimization problem and biobjective tooth profile spur gear optimization problem. Qin and Huang et al. [31] proposed that the scale factor values are randomly sampled from a normal distribution and applied to each parent individual in the current population, which is helpful for the global and local of population. An adaptation scheme is used to generate the crossover rate and is based on the previous experiences of generating promising individuals. However, these self-adaptive strategies cannot well guide the population evolution for this type of ELD problems. It is clearly different that the two adaptive parameters with the learning ability are used to evolve with individuals in this paper. Because the best individuals usually own best control parameters in the improved population, control parameters of other individuals are required to learn from those of the best individuals, which is beneficial to find the global best individual. Same as other self-adaptive strategies, two control parameters are first encoded as shown in Figure 2.
Every individual vector ( ) owns two control parameters ( ) and ( ). In the initial process, ( ) and ( ) are randomly generated in a given range [ min , max ] and [ min , max ], respectively. Since the population is not updated during the first generation, the two control parameters of the failed to update individuals are regenerated, which is the same as their initialization process. In the next evolution, if the current population has been improved in the th generation (marked as fail = 0), two control parameters of the failed to update individuals will learn from those of the best individual. Otherwise (marked as fail = 1), their control parameters are randomly generated in their ranges. If the current individual is successfully improved at generation , its corresponding parameters keep the previous value. That means that control parameters of the current individual are helpful to improve individual. Based on the above explanation, adaptive parameters with learning ability are described by , If ( ( )) > ( ( )) and > 1 and fail = 0 , If ( ( )) > ( ( )) and > 1 and fail = 0 where ( ) and ( ) are the scale factor and crossover rate of the current th individual at generation . ( + 1) and ( + 1) are the scale factor and crossover rate of the th individual at generation + 1.
( ) and ( ) are the best control parameters owned by the best individual in the current population, min and max denote the minimum and maximum scale factor, respectively. min and max denote the minimum and maximum crossover rate, respectively. When control parameters of the current individual learn from the best control parameters, its learning step is a Mathematical Problems in Engineering The late phase The entire evolution phase The early phase The early phase The late phase t  t  t force Figure 1: Schematic diagram of dividing the evolution phase.

Experimental Cases and Parameter Settings.
In this section, there are six ELD cases selected to verify the superiority of PADE compared to the other five DE variants and methods from the literature, and the parameters of experiment are set in advance. To be more specific, these cases include (1) the 6-unit system ( = 1263MW) [41]; (2) the 6-unit system with transmission losses ( = 1200MW) [3]; (3) the 10-unit system with VPE and transmission losses ( = 2000MW) [14]; (4) the 14-unit system with transmission losses ( = 950MW) [42]; (5) the 40-unit system with VPE ( = 10500MW) [2]; (6) the 40-unit system with VPE and POZs ( = 10500MW) [4]. Cases 2-4 consider the transmission losses of power system, and case 3 also considers VPE. Meanwhile, the cases with VPE also include the last two cases, and the last case takes the POZs into consideration. The detailed information of these cases can be obtained from the corresponding literature mentioned above [2-4, 14, 41, 42], and the prohibited operation zones are taken from Table 12 of [4].
The PADE algorithm is compared with other five state-ofthe-art DE variants including differential evolution algorithm with self-adapting control parameters (jDE) [43], differential evolution algorithm with controlled search direction and selfadapting control parameters (jdDE) [44], modified differential evolution algorithm (mDE) [45], differential evolution with composite trial vector generation strategies and control parameters (CODE) [46], and improving differential evolution with a successful-parent-selecting framework (SPS-DE) [29]. The parameters of the six DE variants are set as follows: for jDE, = 0.1, = 0.9, and 1 = 2 = 0.1. The scale factor and crossover rate are, respectively, set to 0.6 and 0.3 in the initial generation, which cannot be obtained from [43]. For jdDE, the values of , , 1 , 2 , and the initial control parameters are the same as the above corresponding values in jDE.
= 0.1 and = 1.0. For mDE, the fixed crossover rate is taken as = 0.8. For CODE, except for the maximum function evaluation number MAX FE [47,48], the other parameters are obtained by [46]. For SPS-DE, the stagnation tolerance is set to 32, and the SPS framework is incorporated with the standard DE with the strategy DE/rand/1/bin, scale factor = 0.6 and crossover rate = 0.3. For the proposed algorithm, min = 0.1, max = 0.9, = 0.1, = 90% × , and max = ⌈3 ⌉. min = 0.1 and max = 0.3 are set for the last two cases, and the crossover rate ranges between 0.3 and 0.6 for the remaining cases. It should be noted that these ELD cases can find the most promising solutions by our parameter settings. For all simulation experiments, the population size is set to 40, and they have the same MAX FE when all six DE variants optimize the same ELD case, which make the results comparison fair. In addition, the maximum generation number of jDE, jdDE, mDE, SPS-DE, and PADE are same for a ELD case. The maximum generation number is equal to 1000, 300, 300, 300, 4000, and 4000 for cases 1-6, respectively. The MAX FE is equal to 40000, 12000, 12000, 12000, 160000, and 160000 for cases 1-6, respectively. When the penalty method is used to handle constraints, the penalty factor is set to 10 20 . Additionally, Matlab8.3 simulation software is used to execute all experiments under the environment of an Intel (R) Core (TM) i5-2450M CPU @ Algorithm: phase-based adaptive differential evolution (1) Initialize a population pop 0 ; (2) Set population size and the dimension of the individual ; the maximum generation number ; the maximum function evaluation number MAX FE; the initial individual control parameters and their bounds; max ; ; an archive and = pop 0 ; = 0; = 0; the cases including the 6-unit system ( = 1263MW) [41], the 6-unit system with transmission losses ( = 1200MW) [3], the 10-unit system with VPE and transmission losses ( = 2000MW) [14], the 14-unit system with transmission losses ( = 950MW) [42], the 40-unit system with VPE ( = 10500MW) [2], the 40-unit system with VPE and POZs ( = 10500MW) [4]; Obtain the global best solution ( ); If is empty (10) Randomly generate three integers 1, 2, and 3 in the range [1, ], Randomly generate two integers 1 and 2 in the range [1, ], and randomly generate a integer between 1 and the size of the archive , and 1 ̸ = 2 ̸ = ̸ = ; (14) V ( ) = ( ) + ( ) × ( 1 ( ) − 2 ( )); (15) End If (16) Else (17) Randomly generate two integers 1 and 2 in the range [1, ], and End If End If (31) End For (32) If ( ( )) > ( ( )) (33) ( + 1) = ( ); 2.50GHZ. The simulation results of 30 independent runs can be obtained in Table 1.

Comparison among Six DE Variants for Six ELD Problems.
In order to improve the quality of solutions, an improved DE and one of the proposed repair methods are combined to find the best feasible solution for a specific case. In this section, all six DE variants are combined with the corresponding one of the proposed two repair methods that are different for the ELD problems without transmission losses and the ELD problems with transmission losses. Thus the performance of the six DEs is only compared in this section while the quality of repair methods will be discussed in Section 4.8.
In Table 1, several criteria are adopted to compare the performance of six improved DEs, including " min ", " max ", " mean ", " median " and " std ", which stand for the minimal, maximal, mean, median and standard deviation values of 30 test results, respectively. For case 1, all six DE variants are able to exactly find the global best solution, and all six DE variants have low std values, indicating that their stability is strong. For case 2, the minimal cost obtained by the first three comparison algorithms are the same as that of PADE, but the performance of PADE is better than that of the three algorithms according to the max and std . The min values of the remaining two comparison algorithms are slightly larger than that of PADE. In terms of the mean , the results of PADE and jdDE are both equal to the minimal cost function value, and they are smaller than those of the other four algorithms. The median values of PADE, jDE, and jdDE are also exactly equal to the minimal cost function value, which are slightly smaller than those of mDE, CODE, and SPS-DE. Overall, the max , mean , and median values of PADE are equal to its min value, and its std value is close to 0. It means that PADE is easy to find the smallest objective function value in each run and has strong stability for this case. For case 3, both jdDE and PADE can find the same lowest cost compared to the other four DEs. Furthermore, the max , mean , median , and std values of PADE are smaller than those of jdDE. Thus PADE is superior to jdDE in stability and success rate. All the results obtained by jDE, mDE, CODE, and SPS-DE are worse than those of the other two competitors, and mDE performs the worst among all methods in terms of the std , explaining that the stability of mDE is poor. With regard to case 4, PADE can obtain the best min , mean and median compared to the other approaches. The min value of mDE is slightly larger than that of PADE, and it is acceptable. However, unfortunately, the max , mean , and std values of mDE are far larger than those of another five approaches. Therefore, although mDE has the strong global search capability close to PADE, it performs poor in stability. The max and std values of PADE are slightly larger than those of jDE, and they are smaller than those of the other four DEs. Therefore, the results of PADE are acceptable. jDE has the advantage over another four comparison algorithms in the average and median cost function values. On the contrary, jDE has poor performance according to the min . In addition, CODE has the worst global search ability among all six algorithms, but its std value is acceptable. It means that CODE has a great chance to get trapped into a local optimum over 30 runs. For case 5, PADE has absolute advantage over another five algorithms, because it is able to find the best min , max , mean , median , and std values, and these values are much smaller than those of another five algorithms. It suggests that PADE performs the best among all six DEs on solving the case 5. For case 6, although system considers POZs, PADE still has great advantage over five competitors. The min , mean , median , and std values of PADE are better than those of another five competitors, which exhibits that the performance of PADE is very competitive. In fact, the ELD problems with POZs are more complex than those without POZs. However, it is clearly seen that the min , max , mean , median , and std obtained by PADE for case 6 are slightly smaller than, slightly larger than, or equal to, those of PADE for case 5. Therefore, the proposed algorithm still has strong search ability of solution space when the ELD problems become more complex. Additionally, the min , mean , and median obtained by jDE for case 6 are even smaller than those of it for case 5. It explains that jDE also has strong search ability for the complex ELD problems. Contrary to PADE and jDE, the results of jdDE and mDE on solving the case 6 are far larger than those of them on solving the case 5 except the std value of mDE. Thus the POZs of the ELD problems make it more difficult to find the feasible solutions for jdDE and mDE. In addition, the results of CODE and SPS-DE are worse than those of PADE except the max value of SPS-DE for case 6. In short, jDE has good performance, but PADE performs the best among them.
In order to obtain the statistical significant difference when the proposed algorithm is compared to those comparison algorithms, Wilcoxon rank-sum test at the 0.05 significance level [51] is used in this section. The marks "+," "=," and "-" represent that the performance of PADE is better than, equal to, and worse than that of the corresponding comparison algorithms, respectively. It should be emphasized that Wilcoxon rank-sum test can only obtain whether the two algorithms have the difference. In other words, the better one cannot be still determined when the two algorithms have the Mathematical Problems in Engineering 13 significant difference. If PADE and one of its competitors do not have significant difference, they are marked as "=." Based on the Wilcoxon rank-sum test characteristic, when PADE and any of the other five DEs have the significant difference, the criterion min is used to determine whether PADE is statistically better than the other algorithm. The better one is marked as "+," and the worse one is marked as "-." The results of Wilcoxon rank-sum test statistical analysis are listed in the last column of Table 1 In sum, PADE is statistically no worse than the other five algorithms in all cases. In addition, PADE significantly beats all competitors for cases 5-6. Therefore, PADE can efficiently solve the six ELD problems.

Convergence of Six DEs for The Fuel Cost Minimization.
Under the condition that the proposed repair methods have been combined with six DEs to easily meet all the requirements of different constraints for the ELD problems, we now devote to study convergence of six DEs for the fuel cost minimization. In order to visualize the convergence characteristics of six DEs, in Figures 3(a)-3(f), we show the average convergence curves of six DEs for the fuel cost minimization of the selected six ELD problems over 30 runs, and these curves are based on the new objective function where the fitness function and penalty function are added for all ELD problems. The abscissa represents the function evaluations, and the ordinate represents the average fuel cost.
From Figures 3(a)-3(f), we can find that the convergence speed of PADE isn't the fastest among all six DEs for all ELD cases, but the average curves of PADE can gradually become better than, or as good as those of the other five DEs in all cases as the evolution progresses. For the first case, the convergence speed of PADE ranks third among all six DEs. In Figure 3(b), mDE, CODE and SPS-DE do not converge to the feasible zone while the other three DEs can find the feasible solutions, and the average curve of PADE is the lowest among three visible curves. In Figure 3(c), the average curve of CODE is the highest among all six DEs, because it does not lie in the feasible range. On the contrary, PADE performs the best according to the average convergence curve. In Figure 3(d), it can be clearly seen that the average curve of PADE reaches a lower level compared to those of the other five DEs. For the first four figures, the feasible solutions are only contained, because their constraint violations descend rapidly. For cases 5-6, the average curves of all six DEs descend very fast at the beginning of evolution. To be more specific, some infeasible solutions contain in the figures for the last two cases. Their constraint violations cannot descend rapidly. Thus, the new objective function values are very large at the beginning of evolution. The declining speed of constraint violations can be seen from the upper bound of ordinate. If constraint violations exist, the penalty function value combining the penalty factor (10 20 ) multiplies and the constraint violations will be very large. However, they can gradually find the feasible solutions. Furthermore, the solutions obtained by PADE are better than those obtained by the other five DEs in complex ELD cases.

The Influence of the Introduced Archive.
In this paper, the introduced archive that stores the recently updated individuals is used to generate the promising feasible solutions with a higher probability in the early phase of evolution. To show the effect of the archive , PADE without archive is compared with PADE with archive, where the archive component of the algorithm is eliminated to test its contributions in improving the quality of solutions. The two PADE algorithms adopt the same proposed repair methods. Furthermore, three criteria are used to compare the performance of the two PADE versions, and they are mean and std of 30 test results and the total constraint violation of the best solution in 30 test results, respectively. The value contains the violations of all the equality and inequality constraints for a solution. The experimental parameters are the same as those of Section 4.2. The comparison results are listed in Table 2.
According to Table 2, the mean obtained by PADE with archive is the same as that of PADE without archive for cases 1-2, and their values are equal to 0 for case 1. The std value of PADE with archive is slightly smaller than that of PADE without archive for the first two cases, meaning that the archive slightly enhances the stability in PADE. For cases 3-4, the mean and std values of PADE with archive are slightly smaller than those of PADE without archive, which indicates that PADE with archive performs more stable than PADE without archive. With regard to the last two cases, the archive improves the quality of solutions more obvious compared to the first four cases. To be more specific, the mean and std values of PADE with archive are far smaller than those of PADE without archive for cases 5-6. For cases 1-4, the performance of the two PADE versions is similar in the constraint violation. PADE with archive achieves a lower level than PADE without archive for cases 5-6 in the constraint violation, which suggests that the archive motivates PADE to find a more satisfactory solution. In all, PADE with archive performs more outstanding than PADE without archive.

The Influence of the PADE Parameters on the Solution
Quality and Stability. In order to investigate the influence of success ratio on the solution quality and stability, the value is set to 0.1, 0.3, 0.5, 0.7, and 0.9, respectively. The other experimental parameters are the same as those of Section 4.2.   To investigate the influence of their adaptive mechanism for two control parameters on the solution quality and stability, PADE with adaptive control parameters is compared with PADE with fixed control parameters by the criteria mean and std . The fixed and values are set to 0.6 and 0.3, respectively. The other experimental parameters are the same as those of Section 4.2. The comparison results are listed in Table 4.
According to Table 4, it is clearly seen from the terms mean and std that PADE with adaptive control parameters is superior to PADE with fixed control parameters for all cases except case 1. For case 1, the solution quality and stability of PADE with adaptive control parameters and PADE with fixed control parameters are same. In all, the proposed adaptive control parameters have a significant effect in improving the performance of PADE for most cases.

Comparison between PADE with Phase-Based and Five Popular Mutation Operators.
To test the effect of the proposed phase-based mutation compared to five popular mutation operators, PADE with phase-based mutation are compared with PADE with different mutation operators on the solution quality and stability. The other experimental parameters are the same as those of Section 4.2. The comparison results are presented in Table 5 In Table 5, the mean value of PADE with phase-based mutation is better than, or as good as that of PADE with other mutation operators for the first two cases. For cases 3-4, the mean value of PADE with phase-based mutation is slightly larger than that of PADE with DE/best/1 and is lower than that of PADE with other mutation operators. It means that the solutions obtained by PADE with phase-based mutation are slightly worse than those of PADE with DE/best/1 but are better than those of PADE with other mutation operators. PADE with phase-based mutation outperforms PADE with another five mutation operators on the solution quality and stability for the last two complex ELD problems.

Comparison between PADE and the Other Methods
Presented in the Literature. In this section, the best results of PADE are compared with those presented in the literature to evaluate the performance of different methods. Different from most literatures, we not only give the best cost function value ( ), but also calculate and list the constraint violations of their best solutions listed in the paper. Here, the constraint violations are caused by power output limits, power balance constraints with or without transmission losses, or POZs for six studied ELD problems. The constraint violations of power output limits, power balance constraints with or without transmission losses and POZs are denoted as , and POZ . For the simplest ELD case, the best solution of case 1 can be obtained according to [41], and its corresponding cost function value is 15275.93039, which is the same as the best cost function value of PADE. Nevertheless, the values of and calculated are, respectively, equal to 0 and 6e-08, which makes the best solution incompletely feasible. The and values of PADE are exactly equal to 0 according to its best solution (446.707273, 171.257992, 264.105658, 125.216766, 172.118858, and 83.593453). Therefore, PADE performs better than that method from [41]. For the remaining five cases, most of the best solutions and corresponding costs have been presented in the literature, and all constraint violations , and POZ are calculated according to those best solutions. The best solutions, , and corresponding constraint violations of inequality and equality constraints are listed in Tables 6-10. Besides, the minority of reported minimal costs exist calculation errors, and those costs are amended according to their best solutions, which is to make the comparisons fair and meaningful. To be more specific, for a reported method, its best cost is verified by bring its best distribution of the generating units in Tables 6-10 into formulas (1)-(3) to recalculate its objective function value. AC denotes the amended cost, and the term "NA" represents "not available." As we all know, the inequality constraints are easy to completely meet, so its error tolerance ( ) is set to 0. However, because the actual ELD problems usually perform nonlinear characteristics especially for the ELD problems with transmission losses, the equality constraint, that is power balance constraint, can be slightly relaxed this restriction of the error tolerance. For the ELD problems with transmission losses, is set to 0.001 according to our repair method. Accordingly, if ≤ , corresponding constraint violations are determine to be acceptable and vice versa for cases 2-4. In addition, since the proposed repair method is easy to completely repair the constraint violations for the ELD problems without transmission losses, the smaller the constraint violations are, the finer corresponding solution is. If constraint violations are equal to 0, the feasible solution can be obtained.
In Table 6, the best results produced by PADE and the other methods from the literature for case 2 are presented. The best cost obtained by PADE is lower than those of QOTLBO [3], TLBO [3], and BSA [49], and it is higher than that of OGHS [20]. OGHS [20] obtains the lowest cost among all methods, but its value is equal to 2.11, which is far larger than the error tolerance ( = 0.001) for this case. Thus the solution obtained by OGHS [20] is infeasible. All values obtained by PADE and the other methods are equal to 0, indicating that their inequality constraints are strictly met. Furthermore, the values of PADE, QOTLBO [3], TLBO [3], and BSA [49] meet the error tolerance for this case, so their solutions are considered as the feasible ones. The best cost obtained by TLBO [3] is the largest on all the best costs obtained by the other methods producing the feasible solutions. In all, PADE outperforms the other methods in achieving the minimal fuel cost while meeting the restriction of constraints.
Regarding the optimization of case 3, Table 7 presents the comparison between PADE and the other methods reported   in the literature. For OGHS [20], its best cost is reported to be 111490. However, its actual cost is equal to 111497.612522 after careful verification, which is slightly higher than that of PADE. The best cost obtained by PADE is lower than those of the other methods, and it is equal to 111497.562722. Besides, the constraint violations of all these methods reach our error tolerance set, suggesting that all solutions are the feasible ones for case 3. Therefore, PADE can achieve the most satisfactory generation dispatch of case 3. Table 8 lists the best results obtained by the three methods for case 4. Among the three methods, the cost obtained by PADE is the lowest, and it is equal to 4303.507984. For the constraint violations of the three methods, all of their values are equal to 0, and all of their values are smaller than 0.001, revealing that all three methods have strong ability in finding the feasible solution. Overall, the performance of PADE is the best compared to the other two methods.
The best results yielded by PADE and the other methods from the literature are shown in Table 9 for case 5. Most methods from the literature report the corresponding best solution, but unfortunately the best solution, and values of IFEP [2], is not available. We get their best cost from the IFEP [2], and it is equal to 122624.35, which is larger than that of PADE. The best cost obtained by -MBA [50] is equal to 121412.5355, which seems similar to that of PADE. However, the value of -MBA [50] is incorrect after verification based on the corresponding solution, and the precise cost is equal to 121578.483653, which is much larger than that of PADE. Except for the -MBA [50] amended, the best cost obtained by FPSO [16] is also amended to 123860.094144. The best costs obtained by GA-DE-PS [17] and NAPSO [4] are slightly higher than that of PADE, and the best costs obtained by the remaining reported methods are far higher than that of PADE. The value of FPSO [16] is nonzero and equal to 461.85; the value of any of the other methods is exactly equal to 0. If the accuracy required is equal to 1e-05, the solutions obtained by -MBA [50], THS [11], and PADE are acceptable. However, if the accuracy required is equal to 1e-06, the solutions obtained by -MBA [50] and THS [11] are infeasible, and the solutions obtained by PADE are still feasible. Therefore, it depends on the accuracy required whether the solutions obtained by other methods are feasible. The values of FPSO [16], NAPSO [4], GA-DE-PS [17], THS [11], and -MBA [50] are larger than that of PADE, indicating that PADE is better than these methods in achieving the feasible solution.
The proposed algorithm and several approaches from the literature are used to optimize case 6, and their best results are given in Table 10. The best costs obtained by NAPSO [4] and KHA-IV [12] are slightly higher than that of PADE, and the best costs obtained by the remaining reported methods are far higher than that of PADE, which suggests that PADE has stronger search ability compared to the other methods. PSO [4], FAPSO [4], and KHA-IV [12] do not report their best solutions and constraint violations. However, they are inferior to PADE in the best cost. On the obtained value, the value of PADE equal to 1.82e-12 is the lowest, and the value of -MBA [50] equal to 72 is the largest, which cannot be neglected. All methods meet the power output limits and POZs completely. Based on the best cost and constraint violation, PADE is more efficient for case 6 compared to the other methods.

Comparison among Different Constraint
Handling Methods. Finally, the selected six ELD cases are employed to test the efficiency of the proposed repair methods and the other reported repair methods on handling constraints, which mainly refers to the equality constraint handling methods. Extensive experiment is carried out to compare the performance of PADE with different constraint handling strategies. Due to difference of the selected problems, we separately propose two types of repair methods for the ELD problems without transmission losses and the ELD problems with transmission losses, and the two repair methods are called PRM1 and PRM2, respectively. In addition, the penalty method (PM) is combined with PRM1 to further handle constraints for the ELD problems without transmission losses, which is represented as PM+PRM1. For the ELD problems with transmission losses, the obtained solution is feasible when ≤ , so PRM2 does not need to combine PM. The penalty method (PM), the penalty method based on the first repair method (PM+RM1) [10], and the penalty method based on the second repair method (PM+RM2) [52] are used to compare with PM+PRM1 on solving the ELD problems without transmission losses. It should be explained that PM does not need to combine with RM2 for case 1 and case 5, and PM is necessary to combine with RM2 to handle POZs for case 6. The penalty method (PM), the penalty method based on the third repair method (PM+RM3) [53], and the fourth repair method (RM4) [9] are used to compare with PRM2 on solving the ELD problems with transmission losses. RM4 also uses an approximate method to repair infeasible solutions, where its error tolerance is the same as that of PRM2, so PM is not combined with RM4. All constraint handling methods are separately incorporated into PADE, and the population size and maximum generation number are the same as the above experiments. The results of 30 independent experiments are separately listed in Tables 11-12, where the terms AOFV, OFVSTD, and TACV represent the average objective function value, standard deviation of the objective function values, and total average constraint violation, respectively.
It can be clearly seen that PADE based on PM+PRM1 performs better than PADE based on the other constrainthandling methods in AOFV and OFVSTD for the ELD cases without transmission losses, and PADE based on PRM2 also performs better than PADE based on the other constrainthandling methods in AOFV and OFVSTD for the ELD cases with transmission losses. It means that the proposed two repair methods make PADE stronger search ability and stable compared to the other constraint-handling methods. In contrast, the AOFVs and OFVSTDs of PM are the largest on all cases, indicating that the search ability and stable of PADE based on PM are poor. The TACVs of PM+PRM1are lower than that of PM for case 1, case 5, and case 6, and the TACVs of PM+PRM1 are equal to 0 for case 1. Although the TACVs of PM+PRM1 are larger than those of PM+RM1 and (PM)+RM2 for cases 5-6, its TACVs will also reach desirable accuracy if the error tolerance is set to 0.00001 for cases 5-6. For cases 2-4, the TACVs of all these methods meet the requirement of accuracy, because the error tolerance is set to 0.001 for the three cases. Therefore, to some extent, the proposed two repair methods exhibit good efficiency in finding feasible solutions for the ELD problems. In all, PM+ PRM1 and PRM2 are very competitive compared to the other constraint-handling methods.

Conclusions
We improve the original differential evolution (DE) algorithm to solve the ELD problems. First, an archive of storing successful individuals is set to have a greater chance to obtain the updated individuals. Second, a phase-based mutation operation is executed, where DE/srand/1 and DE/best/1 are adopted in the early and late phase, respectively. To determine the threshold at the end of the early phase, a success ratio and auxiliary generation index threshold are used together to divide the entire evolution phase into the early and late phase. The mutation operation is beneficial to balance the global and local search when the algorithm solves the ELD problems. Third, the two control parameters are adaptively generated for every individual. Different from other selfadaptive control parameters, the control parameters owned by the individuals have the ability to learn from the best control parameters owned by the global best individual, which further enhances PADE's robustness. In addition to the algorithm's improvement, two types of repair methods are used to meet the requirement of the equality constraint for the ELD problems without or with transmission losses, respectively. As a result, the proposed repair methods can completely eliminate the constraint violation on solving some ELD cases or reduce a lower level compared to the ones for the other ELD cases. Due to the efficient repair methods, PADE can obtain a lot of high-quality feasible solutions. Finally, "-" represents that the presented does not need to be amended.

22
Mathematical Problems in Engineering

Data Availability
The data used to support the findings of this study will be provided by the author upon request.

Conflicts of Interest
The authors declare that there are no conflicts of interest.