Abstract

Dynamic economic dispatch (DED) that considers valve-point effects is a complex nonconvex and nonsmooth optimization problem in power systems. Over the past few decades, multiple approaches have been developed to solve this problem. In this paper, an adaptive differential evolution based on simulated annealing algorithm is proposed to solve the DED problem with valve-point effects. Simulated annealing (SA) algorithm is employed to carry out an adaptive selection mechanism in which the mutation operators of differential evolution (DE) are selected adaptively based on their historical performance. A mutation operator pool consisting of five operators is built to make each operator show its strength at different stages of the evolutionary process. Moreover, a heuristic strategy is introduced to transform infeasible solutions towards feasible ones to enhance the convergence rate of the proposed algorithm. The effectiveness of the proposed methods is demonstrated first on 10 popular benchmark functions with 100 dimensions, in comparison with the classic DE and five variants. Then, it is used to solve four DED problems with 10, 15, 30, and 54 units, which consider the valve-point effects, transmission loss, and prohibited operating zones. The simulation results are compared with those of state-of-the-art algorithms to clarify the significance of the proposed method and verify its performance. Three systems with 100-500 generators are also tested to confirm the advantages of the proposed method on large-scale DED problem.

1. Introduction

Dynamic economic dispatch (DED) is one of the most fundamental problems in modern power system operation and has been the topic of considerable research. The primary objective of DED is to minimize the total fuel cost and meet the load demand in a power system while satisfying system constraints. Solving a DED problem involves formulating a mathematical model of the problem and then selecting an appropriate optimization technique. Conventional mathematical programming based on several deterministic optimization techniques has been employed to solve such problems. These methods include the λ-iteration method, gradient method, and dynamic programming [1].

Solving a DED problem involves formulating a mathematical model of the problem and then selecting an appropriate optimization technique. The methods developed to date can be categorized into conventional mathematical programming methods and heuristic-based optimization algorithms. Previous efforts to solve such problems have employed conventional mathematical programming based on several deterministic optimization techniques, such as the λ-iteration method [2], the gradient method [3], the Lagrange relaxation [4], and dynamic programming [5]. Considering the high nonlinearities and discontinuities in the real input-output characteristics of DED, these conventional methods often fail to adequately solve this problem. Considering the high nonlinearities and discontinuities in the real input-output characteristics of ED, these conventional methods often fail to solve this problem well.

Unlike conventional algorithms, modern heuristics optimization algorithms have the advantages of global search capability and the lack of requirement for certain mathematical properties of the objective function to be satisfied. Thus, these algorithms have proven their effectiveness in solving complex and nonconvex problems such as DED [2]. Representative algorithms include the genetic algorithm (GA) [6], particle swarm optimization (PSO) [7], simulated annealing (SA) [8], artificial immune system (AIS) [9], artificial bee colony (ABC) [10], evolutionary programming (EP) [11], differential evolution (DE) [12], harmony search (HS) [13], seeker optimization algorithm (SOA) [14], imperialist competitive algorithm (ICA) [15], teaching-learning algorithm (TLA) [16], bacterial foraging optimization (BFO) [17], and others. In additional, to further enhance the global search ability of these original algorithms, representative modified variants, such as modified PSOs [1821], modified DE [12, 2224], modified ABCs [2527], and others, have been developed over the past years. In addition to these innovations, hybrid methods that combine two or more approaches have been applied to solve this problem, such as EP-PSO [11], the hybrid genetic algorithm based on differential evolution (HGA-DE) [28], the hybrid PSO with sequential quadratic programming (PSO-SQP) [29], the hybrid GA-pattern search-SQP (GA-PS-SQP) [12], BCO-SQP [30], a modified EP-PSO (MHEP-SQP) [19], and SOA-SQP [14]. By combining various approaches, the hybrid method overcomes a single component’s limitations by applying other components’ advantages to solve the DED problem. One combined approach has been applied to implement coarse searching, and another has been used for fine tuning in a hybrid method [31]. Although many of the aforementioned methods have provided high-quality solutions, they have their own limitations. Heuristic-based approaches are often trapped by the sensitivity of various parameter settings. With hybrid methods, finding an appropriate integration point for the combined approaches is always a challenging problem for designers. Some representative methods are listed in Table 1. The pros and cons of the studied works can be summarized by Table 2.

Differential evolution algorithm (DE) [32], which is proposed by Storn and Price, is a population based algorithm that simulates the biological evolution. Compared with other swarm intelligence algorithms, DE is a simple but powerful method for global numerical optimization with nonsmooth and nonconvex characteristics. Due to its simplicity and efficiency, DE has been successfully applied to many complex engineering fields, such as data mining [33], scheduling [34], energy distribution [35], molecular engineering [36], pattern classification [37], and so on.

However, DE still has some drawbacks. The successful application of DE mainly depends on its mutation operators and parameter settings. Moreover, DE often suffers the problem of premature convergence, which may lead the searching progress into local optimal. In addition, DE lacks a mechanism to deal with the complicated constraints of the DED problem [23].

In this paper, an adaptive differential evolution based on simulated annealing (ADE-SA) algorithm is proposed to solve the large-scale DED problem with valve-point effects. In ADE-SA, DE is employed to search the feasible solution of DED. Five widely used mutation operators are collected to form a mutation operator pool. SA is employed to carry out an adaptive selection mechanism in which the mutation operators are selected adaptively according to their performance in the search history. Thus, individual can adaptively select appropriate mutation operators in the different stages of the search process. Furthermore, a heuristic strategy is introduced to transform infeasible solutions towards feasible ones to enhance the convergence rate of the proposed algorithm. The performance of ADE-SA is evaluated on a series of DED problems, and compared with other eight state-of-the-art algorithms. Comparison results indicate the effectiveness of the proposed algorithm. Three systems with 100-500 generators are also tested to confirm the advantages of the proposed method on large-scale DED problem.

The rest of this paper is organized as follows: Section 2 provides the formulation of the DED problem. The basic SA and DE scheme is briefly described in Section 3, while the framework of the proposed method is also presented in this Section. To illustrate its effectiveness, the proposed approach was tested on a series of systems in Section 4. Section 5 offers the conclusions. Parts of this paper, such as Sections 1 and 3.3.4, are the extended version of our completed research [38].

2. Problem Formulation

DED is a multiconstrained and nonlinear optimization problem. The primary objective of DED is to minimize the total fuel cost and meet the load demand in a power system while satisfying system constraints. Consider a power system with NG generating units. In general, the DED problem can be formulated as the following optimal function [13]:where , , and are the fuel cost coefficients of generator . is the time horizon and is the total number of generators of the system. is the power output of generator at time interval . denotes the total fuel costs over the dispatch period. When the valve-point effect (VPE) is considered, the total production costs can be written as follows [13]:where and are the coefficients of generator reflecting the valve-point effects.

The minimized DED problem should be subjected to the following equality and inequality constraints:

(1)  Real power balance constraintswhere and are the total load demand and the system transmission loss of time interval , which is expressed as a quadratic function using B coefficients [31].

(2)  Power operating constraintswhere is the maximum power outputs of generator .

(3)  Generator ramp rate constraintswhere and are the ramp-up and ramp-down rate constraints of generator , respectively.

3. Methods

3.1. Brief Introduction of Simulated Annealing Algorithm

Simulated annealing (SA) [49] is a metaheuristic that can escape from local optimal with its use of mountain-climbing move. SA is a popular technique over the past 20 years due to its ease of implementation and convergence characteristics. It is usually used for solving discrete and continuous optimization problems. SA comes from the imitation of the physical annealing process. In this process, the crystalline solid is heated and then cooled slowly until it reaches its most stable molecular structure.

At each iteration of SA, the current solution and the new solution are compared. The improved solution is always accepted, while a part of the nonimproved (inferior) solution is also accepted to make the individual escape local optimum. The probability of accepting a nonimproved solution depends on a temperature parameter which is presented in where dif is the difference between the new solution and the original solution. T is the current temperature.

SA reaches steady state through a slow drop in temperature, that is, the temperature determines when the iteration terminates. Therefore, a temperature cooling function is needed to reduce the temperature from high temperature to high temperature. The cooling function is used to lower the temperature which is defined aswhere is a positive number less than 1.

3.2. Brief Introduction of Differential Evolution

The basic DE consist of three phases: mutation, crossover, and selection. In the mutation phase DE randomly selects three distinct individuals , , and from the population such that . G represents the current generation. These three random individuals are employed to generate a mutated individual according towhere is a scaling factor.

After the mutation phase, DE introduces crossover operation to increase the diversity of the population. The offspring individual is generated from its parents and the mutated individual using the following crossover rule:

where Cr is a crossover constant and j denotes the jth element of individual. is a uniformly distributed random number and is a randomly chosen index which guarantees that gets at least one element from the mutated individual .

In the selection phase the fitness value at the individual, , is compared with ; the winner will be the survivor. The selection operation can be expressed as the following rule:

3.3. Adaptive Differential Evolution Based on Simulated Annealing (ADE-SA)
3.3.1. Mutation Operator Pool

It is well-known that DE needs relatively less parameter tuning than many other swarm intelligence algorithms. In the mutation phase, the new individuals are generated by selecting three individuals randomly, which may cause over diversity in population and slow convergence. According to [50], DE can make use of its adaptive mutation strategy to improve its performance. Therefore, a mutation operator pool composed of five different mutation operators is built in this section. These operators will be implemented under the guidance of an adaptive selection mechanism at different stages of search. Through the analysis of the widely used evolutionary strategies in swarm intelligence algorithms, a total of 5 mutation operators are selected to make up the mutation pool (the first operator M1 is presented in (9)), which are as follows:

M2. DE/rand/2

M3. DE/best/1

M4. DE/current-to-best/1 derived from PSO

M5. Employed bees phase in ABCwhere r1, r2, r3, r4, and r5 are mutually integers of randomly selected food sources except for i. is the difference vector, and is a random number in the range of in ABC.

3.3.2. Adaptive Selection Mechanism Based on Simulated Annealing

In fact, each mutation operator may shows different performance levels at different stages of search process. Operators that show effectiveness in some stages of the iterative process may not perform well in other stages. In this section, a learning cycle is utilized to guide the adaptive selection of these operators. It records the historical performance of each mutation operators and calculates the chosen probability of each operator for the next cycle. The learning cycle should be much smaller than the whole iteration process to reflect the different stages of the search process. The learning cycle, which is proposed in [51], is employed in this paper. The details of the learning cycle are presented as follows.

Step 1. Update the total calls of all mutation operators .

Step 2. Evaluate the new individual and the original individual , and calculate the difference of their fitness functions, .

Step 3. Use the Metropolis acceptance criterion to decide whether to accept a new individual and update the parameter and according to the result.

Step 4. Check if the reheating phase () is triggered, update the parameters and calculate the chosen probability of each operator. If not, only calculate the chosen probability.

Step 5. If the iterative process completes a learning cycle, terminate the loop. If not, return to Step 1.
Where is the total number of calls of an operator during the current learning cycle, is the number of individuals that passed the Metropolis acceptance criterion. counts the total number of accepted operator calls in a learning cycle. is the stopping acceptance ratio which is used to trigger the reheating phase of SA. The temperature will be increased to a former temperature. This process continues until a better solution is found. The temperature is increase using

3.3.3. Framework of ADE-SA

The framework of ADE-SA can be summarized as follows.

Step 1. Set initial nonimproving acceptance ratio and stopping nonimproving acceptance ratio . Estimate the starting temperature , stopping temperature , total iterations , number of iterations for DE, and the length of single learning cycle . Set, , , and .

Step 2. Define the mutation strategy pool; initialize the associate parameters , , minimal weight , and an initial weight .

Step 3. Generate the initial population and calculate the fitness value based on the fitness function.

Step 4. Select an operator based on and generate a candidate solution using operator .

Step 5. Judge whether the candidate solution is acceptable or not based on the Metropolis criterion in SA.

Step 6. Check the reheating condition; if necessary, increase the temperature according to the function .

Step 7. Check if is satisfied, cool the temperature based on .

Step 8. Record the history information of current low-level heuristics and check if the learning period is meet, if so, update the associate parameters.

Step 9. If the stopping criteria satisfied, terminate the main loop. If not, return to Step 3.

3.3.4. Heuristic Constraint Handling

DED is a nonlinear constrained optimization problem with many equality and inequality constraints. Some solutions will not satisfy all constraints during the iterative process, such as the real power balance and ramp rate constraints. Penalty method is one of the most frequent-used approaches for complex constraints of DED, but the method is not suitable for improving the computational efficiency. In this section, a heuristic procedure has been employed to accommodate constraints [38]. Taking a solution as an example, the following mechanism is applied to transform an infeasible solution into feasible solution.

Step 1. Randomly select an hour as the initial time interval.

Step 2. The power operating constraints are addressed based on the ramp rate constraints as follows:Note that even if a feasible solution is obtained at current hour, it may make it impossible for to meet (3) and (6) at the next hour. Equations (17) and (18) can make sure that the generation range of the adjacent time interval satisfies the load demand. This prevention strategy helps to avoid generating a potential which affects the equality constraints adjustment at the next hour.

Step 3. Adjust the output of each generator at current time interval to satisfy (5) and (6) using

Step 4. Tackle the demand-supply constraint by the following steps.
Step 4.1. Calculate the system transmission loss at the current time interval using (4).
Step 4.2. Calculate the difference in the power demand constraints at the current time interval asCheck if , go to Step 2 and set ; otherwise, go to the next step. Here, is a tolerance variable which has a large value at the early stages of evolutionary process and is decreased gradually to a small value [31].
Step 4.3. Obtain a random generator from generators as the slack generator, and reset its output at current interval asStep 4.4. Adjust the output to the feasible horizon using (18). If the solution of is still infeasible, recalculate the using (19) and select another generator at random from generators as the slack generator. Repeat this process until a feasible solution is found. Then update the operating range using (20) and (21), and set .
Step 4.5. Repeat Step 4.1 to Step 4.4 until . Here we will obtain a feasible solution at to time intervals.

Step 5. Obtain the solution at the rest of the time intervals by repeating Step 4.1 to Step 4.4. The strategy is similar to the preceding process, except the solution is calculated from the initial time interval to the 1th time interval.

Step 6. If , terminate the above process and output the updated feasible solution .

4. Numerical Experiments and Results

To assess the efficiency of the proposed approach, 10 widely used benchmark functions are applied in this section. ADE-SA is compared with several recent popular modified DE algorithms, including IDE [12], CDE [24], HDE [39], and IMMSADE [52], respectively.

4.1. Benchmark Functions

The details of 10 commonly used benchmark functions are given in Table 3. Table 3 lists the 10 benchmark functions that are selected: F1~F5 are unimodal functions. F6 has one minimum and is discontinuous. F7 is a noisy quadratic function. F8~F10 are multimodal functions. represents the dimension of the benchmark function, S represents the search ranges of the solutions, and is the value of the global optimum. 30 independent runs with a commonly used 300000 max number of function evaluations (FES) are conducted. The simulation parameters chosen in this paper are as follows: , , , , population size , , and .

4.2. Simulation Results and Discussions

The mean and standard variances obtained by each algorithm on each test functions with 100 dimensions are summarized in Table 4. It can been seen in Table 4 that ADE-SA outperformed all other modified DE algorithms on seven functions, including f1, f2, f3, f5, f6, and f7. ADE-SA have the same best results as well as HDE and IMMSADE on functions f9 and f10, while they are both much better than other algorithms on these two test functions. For f4 and f8, HED and IMMSADE achieved the best solutions, respectively. But they are just slightly better than ADE-SA. So we can draw a basic conclusion from Table 4 that ADE-SA is good at solving unimodal and multimodal functions, and it has the advantage of solving high dimension problems.

Figure 1 shows the convergence curve of the compared 7 algorithms on 4 representative benchmark functions with different types. We can easily find that the search accuracy and the convergence rate of ADE-SA are better than the classic DE and all the modified DE variants on those of different types of test functions. The effectiveness of the proposed algorithm is fully demonstrated according to the convergence characteristics, and ADE-SA also shows a stable robustness. Although in some case ADE-SA converged slightly slower than other DE variants, it could also find high-quality result in the later stage of iteration. Overall, the proposed method presents a competitive search capability and rapid convergence on most of the benchmark functions.

5. Experimental Results on DED Problems and Discussions

In this section, to assess the efficiency of the proposed approach, we consider several test problems that involve up to 500 units in DED problems.

Case 1. 10-unit system with and without power loss [31].

Case 2. 15-unit system with prohibited operation zones and power loss [53].

Case 3. 30-unit system generated by combining three 10-unit systems of Case 1.

Case 4. 54-unit system without power loss [46].

Case 5. 100-unit system generated by combining ten 10-unit systems of Case 1 [43].

Case 6. 200-unit system generated by combining twenty 10-unit systems of Case 1.

Case 7. 500-unit system generated by combining fifty 10-unit systems of Case 1.

In all cases, the valve-point effect and the ramp rate constraints are considered. The complete dispatch horizon is selected as one day with 24 hours.

5.1. Experimental Settings

The following values are used to evaluate the quality of the solutions: best/minimum total fuel cost , average total fuel cost , worst/maximum total fuel cost , and standard deviation . The simulation parameters chosen in this paper are as follows: , , , , population size , , and . The hourly load demands for case 1 are shown in Table 5 and the unit data are adapted from [15] and listed in Table 6. For all the cases, the simulation contains 50 independent runs to ensure the robustness of the proposed method. Simulations were executed using the MATLAB 2016a computational environment on a Pentium-core 3.2 GHz personal computer with 8 GB of RAM.

5.2. Result and Comparison

Case 1. 10-generator problem with and without considering the transmission loss.

The first power test system consists of ten thermal generators without considering the transmission loss. The generating system data were taken from [31] and given in Table 6. The hourly system load demand for this case is listed in Table 5. The optimal generation dispatch values of the case without and with considering the transmission loss are shown in Tables 7 and 9, respectively. Tables 8 and 10 give the corresponding statistic comparison among different heuristic methods; 20 previously developed algorithms were used for comparison. It is obvious that ADE-SA can obtain solutions with much lower total generation costs than other reported methods. For the case without transmission loss, although MIQP [31] can obtain lower standard deviation; ADE-SA outperforms it in terms of the total generation costs. Even the maximum costs obtained using these three methods are better than the solution obtained by MIQP.

Case 2. 15-generator problem with considering the prohibited operation zones and the transmission loss.

In this case, the power demand, operating unit characteristics, and B coefficients for calculation of transmission losses were given in [46]. The optimal generation dispatch values of the case with considering the transmission loss are shown in Table 11. From Table 12, it is obvious that ADE-SA gave the best performance on minimum and average fuel cost compared with HS [46], HIS [46], GHS [46], SGHS [46], and NPAHS [46].

Case 3. 30-generator system generated by combining three 10- generator systems of Case 1.

The second case is a lossless electric power system consisting of 30 thermal generators which is tripled from the 10-generator problem of Case 1. The minimum, maximum, and average values of total operating cost for one day are given in Table 13, along with the operating costs from other published results from the literature for the sake of comparison. It is apparently seen that ADE-SA outperforms any other method regardless of minimum cost, average cost, and maximum cost.

Case 4. 54-generator system.

In this case, a challenging problem with 54 generators is employed. The data and information of this system can be obtained from [35]. The simulation results of the proposed methods are compared in Table 14 with GA [15], PSO [15], ICA [15], HS [35], IHS [46], GHS [46] and NPAHS [46]. Note that even the maximum cost obtained by these four methods, which is 1,633,187 $/day, is better than the minimum cost of the other methods.

Case 5 to 7. 100-generator, 200-generator, and 500- generator systems.

To test the performance of ADE-SA to deal with large-scale DED problems, we apply ADE-SA to three larger-scale DED problems in case 5 to case 7, which contain 100 generators, 200 generators, and 500 generators derived from the 10- generator system of case 1. Compared with previous cases, the increase of the number of process variables and constraints in Case 5 to Case 7 makes the DED problem much more complex. For Case 5, the simulation results obtained by ADE-SA are compared with GA [45], PSO [45], FA [45], FAPSO [47], SAFA [47], and CSO [43]. It can be seen from Table 15 that ADE-SA outperforms any other method regardless of minimum cost, average cost, and maximum cost. For Case 6 and Case 7 with 200 and 500 generators, only CSO reported results so far.

5.3. Learning Cycle Study

In the learning cycle, the performance of each mutation operators is evaluated and used to determine the chosen probability in the latter stage of search process. Figure 2 shows that five mutation operators are adapting to the whole evolution process. The chosen probability of each operator can be seen from the figure. In the initial of the evolution process, M1 and M2 showed high chosen probability whereas the other operators showed low probability. As the iteration proceeds, the chosen probability of M4 and M5 increased continuously.

5.4. Adaptive Selection Mechanism Study

To demonstrate the effect of the proposed adaptive selection mechanism, the average of best objective function values over 50 independent runs of each mutation operator are recorded in Table 16 (take Case 1, for example). From this table, it can be seen that each operator shows poor performance when runs depend on themselves. When the adaptive selection mechanism is applied, the operators are able to obtain high-quality solutions.

5.5. Heuristic for Constraint Handling Study

The proposed heuristic procedure can help the infeasible solutions to transform into high-quality feasible solutions. To test the effect of this heuristic procedure, we record the average cost over 20 independent runs and the failure times of constraint handling with and without the heuristic procedure in the same experiment situation. The penalty method, which is one of the most frequent methods for dealing with DED problem, is used for comparison. Table 17 shows that each operator with the heuristic procedure can find better solution. In the constraint handling process, the infeasible solution is adjusted hour by hour. When a feasible solution is obtained at current hour, it may make it impossible for to meet the power balance and ramp rate constraints at the next hour. Failure times counts the times of the situation. A lower failure times means an effective constraint handling strategy. All the experiments have the same 2000 generations.

6. Conclusions

This paper proposes and successfully applies a novel method to solve the DED problem considering value-point effects. DE is employed to search the feasible solution of DED. Five widely used mutation operators are gathered to form a mutation operator pool. Simulated annealing algorithm is used to carry out an adaptive selection mechanism in which the mutation operators of DE are selected adaptively. The effectiveness of the proposed methods is demonstrated first on 10 popular benchmark functions with 100 dimensions, in comparison with the classic DE and five variants. Then, seven test systems having 10, 15, 30, 54, 100, 200, and 500 generators are analyzed to illustrate the application of the proposed method. The computational results show that the proposed method has the capacity to provide stable and high-quality solutions efficiently for DED problems. Comparison of the results with algorithms recently reported in the literature confirms the superiority of the proposed method and its higher probability of finding accurate and optimal solutions for DED problems. Finally, the performance of different components is discussed.

In future work, we plan to extend the present work in the following directions. First, adaptive techniques for parameter setting of the mutation operator will be investigated. Second, we will apply the proposed method to more complex practical DED problems that consider the uncertainty of load demand, wind speed, and solar irradiation.

Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

This work was supported in part by the National Natural Science Foundation of China under Grants nos. 61773105, 61374147, 61703085, 61004083, and 61533007 and in part by the Fundamental Research Funds for the Central Universities under Grants nos. N150402001 and N120404014.