Differential evolution for economic load dispatch problems

Available online at www.sciencedirect.com

Electric Power Systems Research 78 (2008) 1322–1331

Differential evolution for economic load dispatch problems Nasimul Noman a,∗ , Hitoshi Iba b a

Department of Computer Science and Engineering, University of Dhaka, Dhaka, Bangladesh b Iba Laboratory, Graduate School of Frontier Sciences, University of Tokyo, Japan

Received 15 March 2007; received in revised form 15 October 2007; accepted 22 November 2007 Available online 19 February 2008

Abstract In this work, differential evolution (DE) algorithm was studied for solving economic load dispatch (ELD) problems in power systems. DE has proven to be effective in solving many real world constrained optimization problems in different domains. ELD problems are complex and nonlinear in nature with equality and inequality constraints and here special measures were taken to satisfy those. Five ELD problems of different characteristics were used to investigate the effectiveness of the proposal. Comparing with the other existing techniques, the current proposal was found better than, or at least comparable to, them considering the quality of the solution obtained. © 2007 Elsevier B.V. All rights reserved. Keywords: Economic load dispatch; Differential evolution; Non-smooth cost function; Transmission loss

1. Introduction Economic load dispatch (ELD) has become an essential function in operation and control of modern power system. The ELD problem can be defined as determining the least cost power generation schedule from a set of online generating units to meet the total power demand at a given point of time [1]. Though the core objective of the problem is to minimize the operating cost fulfilling the load demand, various types of physical and operational constraints make ELD a highly nonlinear constrained optimization problem, particularly for larger systems [2]. However, careful and intelligent scheduling of the units can not only reduce the operating cost significantly but also assure higher reliability, improved security and less environmental impact [3]. Traditionally, in ELD the input–output characteristics (or cost function) of a generator is approximately represented by using a single quadratic function. Practically, operating conditions of many generating units need the cost function to be represented as piecewise quadratic function [4]. However, higher-order nonlinearities and discontinuities are observed in real input–output characteristics due to valve-point loading in fossil fuel burning plants [5]. Besides, the units may have prohibited operating ∗

Corresponding author. E-mail addresses: [email protected] (N. Noman), [email protected] (H. Iba). 0378-7796/$ – see front matter © 2007 Elsevier B.V. All rights reserved. doi:10.1016/j.epsr.2007.11.007

zones due to faults in machines or associated auxiliaries, such as boilers, feed pumps, etc. leading to instabilities in certain ranges of the unit loading [6]. Furthermore, the operating range for online units is actually restricted by their ramp-rate limits [7]. These and other constraints transform an ELD problem into a hard nonconvex optimization problem. Because of the highly nonlinear characteristics of the problem with many local optimum solutions and a large number of constraints, the classical calculus-based method and Newton-based algorithms cannot perform very well, respectively, in solving ELD problems. Though dynamic programming is not affected by the nonlinearity and discontinuity of the cost curves, it suffers from the “curse of dimensionality” and local optimality [8]. Among the other artificial intelligence approaches, neural networks have been successful to solve the ELD problems with different constraints [4,9]. However, some of these neural network-based approaches may suffer from excessive numerical iterations, resulting in huge calculations. Evolutionary algorithms (EAs), such as genetic algorithm (GA), evolutionary strategy (ES) and evolutionary programming (EP), are faster than simulated annealing (SA) because of their inherent parallel search technique. Besides, other advantages of EAs, such as global search capability, robust and effective constraints handling capacity, reliable performance and minimum information requirement, make it a potential choice for solving ELD problems. Consequently, EAs have received

N. Noman, H. Iba / Electric Power Systems Research 78 (2008) 1322–1331

much attention in solving ELD problems. However, as the epsitasis among the parameters increases with the dimension of the problem, GAs may converge prematurely in local optima [10]. Another EA, namely particle swarm optimization (PSO), inspired by social behavior of bird flocking or fish schooling, has been found to be more robust in solving such nonlinear optimization problems. Therefore, many researchers have tried PSO and its hybrids in solving ELD problems [2,10–12]. Differential evolution (DE) is one of the most prominent new generation EAs, proposed by Storn and Price [13], to exhibit consistent and reliable performance in nonlinear and multimodal environment [14] and proven effective for constrained optimization problems [15–18]. The advantages of DE over other EAs, like simple and compact structure, few control parameters, high convergence characteristics, have made it a popular stochastic optimizer. Empirical studies have shown that DE can outperform many other well known EAs [13,19] and also PSO [20]. Because of its powerful and reliable search capability, it has got many real work applications in various fields such as pattern recognition, communication, mechanical and chemical engineering and biotechnology [14]. Even some variants of DE have been attempted to solve ELDs [21–23] though without any generalized conclusion about the behavior of the algorithm in solving them. In this work, we have studied the prospect of classic DE in solving ELD problems. Specialized constraint handling mechanisms have been used for handling different equality and inequality constraints. Performance comparison with other algorithms demonstrates the effectiveness of DE in solving ELD problems of different characteristics. The paper is organized as follows. The next section of the paper presents the formulation of ELD problem as a constrained optimization problem. Section 3 contains a brief overview of DE with the proposed constraint handling mechanism. Section 4 reports the different experimental result comparing with the other methodologies. The paper ends in Section 5 with a brief discussion on results.

1323

where the fuel input–power output cost function of is i-th unit is represented by the function Fi . The most simplified fuel cost function Fi (Pi ) for generator i loaded with Pi MW is approximated by a quadratic function as follows: Fi (Pi ) = ai + bi Pi + ci Pi2

(2)

where ai , bi and ci are the fuel cost coefficients of the i-th generatic unit. In reality, the generating units with multi-valve steam turbine have very different input–output curve compared with the smooth cost function. Therefore, the inclusion of the valve-point loading effects makes the representation of the incremental fuel cost function of the generating units more practical. The incremental fuel cost function of a generating unit with valve-point loadings is represented as follows: Fi (Pi ) = ai + bi Pi + ci Pi2 + |ei sin(fi (Pimin − Pi ))|

(3)

where ei and fi are the coefficients of generator i reflecting the valve-point effects. Another case of using non-smooth cost function, for a more realistic representation of generating units, is the multiple fuel problem where the fuel cost function is expressed as the piecewise quadratic functions as follows: ⎧ ai1 + bi1 Pi + ci1 Pi2 ⎪ ⎪ ⎪ ⎪ ⎨ ai2 + bi2 Pi + ci2 Pi2 Fi (Pi ) = . ⎪ . ⎪ ⎪ ⎪. ⎩ aik + bik Pi + cik Pi2

if Pimin ≤ Pi ≤ Pi1 if Pi1 ≤ Pi ≤ Pi2 .. .

(4)

if Pik−1 ≤ Pi ≤ Pimax

where aik , bik and cik are the cost coefficients of generating unit i for fuel type k, where k = 1, 2, . . . , k. The economic dispatch of generating units for the given load demand to be done by satisfying the following constraints. 1. Power balance constraint: N 

2. ELD problem formulation The primary objective of ELD problem is to determine the most economic loading of the generating units such that the load demand in the power system can be met [24]. Additionally, the ELD planning must be performed satisfying different equality and inequality constraints. In general, the problem is formulated as follows. Consider a power system having N generating units, each loaded to Pi MW. The generating units should be loaded in such a way that minimizes the total fuel cost FT while satisfying the power balance and other constraints. Therefore, the classic ELD problem can be formulated as an optimization process with the objective: min FT = min

N  i=1

Fi (Pi )

(1)

Pi − (PD + PL ) = 0

(5)

i=1

where PD is the total load demand and PL is the total transmission loss. Calculation of PL using the B-matrix loss coefficients is expressed as a quadratic function: PL =

N  N 

Pi Bij Pj +

i=1 j=1

N 

B0i Pi + B00

(6)

i=n

2. The generating capacity constraint is given by Pi min ≤ Pi ≤ Pi max

(7)

where Pi min and Pi max are the minimum and the maximum power outputs of the i-th unit. 3. Ramp-rate limit constraint: max(Pimin , Pi0 − DRi ) ≤ Pi ≤ min(Pimax , Pi0 + URi )

(8)

1324

N. Noman, H. Iba / Electric Power Systems Research 78 (2008) 1322–1331

where Pi0 is the previous output power, URi and DRi is the upramp and downramp limit of the i-th generating unit. 4. Prohibited operating zones constraint: ⎧ min l ⎪ ⎨ Pi ≤ Pi ≤ Pi,1 u l Pi ∈ Pi,k−1 ≤ Pi ≤ Pi,k ⎪ ⎩ u Pi,ni ≤ Pi ≤ Pimax

k = 2, . . . , ni

(9)

l and where ni is the number of prohibited zones for unit-i, Pi,k are the lower/upper bounds of the k-th prohibited zone of unit-i. Inclusion of additional constraints such as line flow (not considered in this study), etc. can result in a more accurate modeling of ELD problems. u Pi,k

3. Differential evolution In this section, we review the differential evolution (DE) algorithm that was used for searching the optimum solution of ELD problems. DE is a population-based stochastic search technique that works in the general framework of EAs. The design principles of DE are simplicity, efficiency and use of real coding. It starts to explore the search space by randomly choosing the initial candidate solutions within the boundary. Then the algorithm tries to locate the global optimum solution for the problem by iterated refining of the population through reproduction and selection. For the ELD problem with N generating units, each individual of the DE population consists of a N-dimensional trial vector x = {P1 , P2 , . . . , PN }, where the i-th element of x represents the power output of the i-th generating unit. In DE, unlike many other EAs, each individual of the current generation breeds its own offspring by mating with other randomly selected individuals. Specifically, for each individual xiG , i = 1, . . . , P, where G denotes the current generation, three other random individuals xjG , xkG and xlG are selected from the population such that j, k and l ∈ {1, . . . , P} and i = j = k = l. These three random individuals are employed to generate a mutated individual vG i using a differential mutation operation, according to the following equation: G G G vG i = xj + F (xk − xl )

(10)

where F, commonly known as scaling factor or amplification factor, is a positive real number, typically less than 1.0, that controls the rate at which the population evolves. After the mutation phase, DE uses a binomial crossover operation, in which the G mutated individual vG i is mated with xi and the offspring or G trial individual ui is generated. The genes of uG i are inherited from xiG and vG , determined by a parameter called crossover i probability (Cr ∈ [0, 1]), as follows:  vG i,t if r(t) ≤ Cr or t = rn(i) G ui,t = (11) G xi,t if r(t) > Cr and t = rn(i)

where t(= 1, . . . , N) denotes the t-th element of individual vectors. r(t) ∈ [0, 1] is the t-th evaluation of a uniform random number generator and rn(i) ∈ {1, . . . , N} is a randomly chosen G index which ensures that uG i gets at least one element from vi . From the above description, it is clear that DE does not select individuals for reproduction, rather applies selection pressure exclusively for picking survivors. A knock-out competition is played between each individual xiG and its offspring uG i and the winner, selected deterministically based on fitness values, is promoted to the next generation. The selection operation can be expressed as follows:  G if f (uG uG i i ) ≤ f (xi ) xiG+1 = (12) xiG otherwise The above steps of reproduction and selection are repeated generation after generation until some stopping criteria are satisfied. 3.1. Handling the constraints of ELD problems Michalewicz and Schoenauer [25] grouped all the constrained handling techniques that came out over the last few years into four categories: (1) methods based on preserving feasibility of solution often using some specialized operators to transform infeasible solutions into feasible ones, (2) methods based on penalty functions where fitness of infeasible individuals are penalized in different ways, (3) methods which make a clear separation between feasible and infeasible solutions and often prefers a feasible one with lower objective value over an infeasible one with higher objective value and (4) hybrid methods that combines evolutionary techniques with deterministic procedures. Among these, the most common approach to handle constraints is the methods based on exterior penalty functions. However, penalty method has several limitations such as difficulty in selecting the appropriate penalty factors, poor performance under certain specialized conditions, etc. Similarly, constraint handling techniques belonging to other categories also have advantages and disadvantages. In handling constraints while optimizing using DE, different techniques have been used [15,16,26]. Some of these approaches are generalized and some others are problem specific. In this work we have used a specialized method for handling the constraints of ELD problem. The technique was used with a target of satisfying the equality and inequality constraints. Though this method was used for ELDs, it can be readily applied to other problems with similar equality and inequality constraints. The constraint handling procedure can be articulated as follows. 3.1.1. Handling the generating capacity constraints The initial individuals (xiG=1 , i = 1, 2, . . . , P) of DE algorithms are created randomly, where each element PjG=1 (j = 1, 2, . . . , N) of xiG=1 , representing the power output of the generating unit j, is uniformly sampled from [Pjmin , Pjmax ]. Therefore, the generating capacity constraints of Eq. (7) can be easily satisfied, during the initialization process, by randomly selecting the power outputs within the capacity ranges.

N. Noman, H. Iba / Electric Power Systems Research 78 (2008) 1322–1331

1325

Fig. 1. Power balance constraint handling mechanism.

However, when offspring are generated applying recombination operators of DE, their components may not fulfill the constraint of Eq. (7). To keep the offspring feasible in terms of the aforementioned bound-constraints, offspring variables that violate boundary constraints are reflected back from the violated boundary by the amount of violation using the following rule [27]:

uG i,j

⎧ min 2P − uG ⎪ i,j ⎨ j max = 2Pj − uG i,j ⎪ ⎩ G ui,j

min if uG i,j < Pj

max if uG i,j > Pj

(13)

up to a certain threshold (here 0.01 MW in Fig. 1(b)). However, all individuals must pass through the appropriate adjustment treatment before evaluating their fitness. 3.2. Fitness evaluation criteria To evaluate the suitability of a DE individual for survival we need to define a fitness function (f ()) that will assign the fitness value. The generalized fitness function used in this work is as follows:

otherwise

3.1.2. Satisfying the power balance criteria In order to treat all the constraints in a unified way, often the equality constraints are converted into inequality constraints using a very small value . Here, we have used a mechanism that meets the power balance equation by adjusting the power output of generating units within their capacity range. After generating an individual X by the initialization procedure or by the breeding process and after refining it for satisfying the capacity constraints, it is adjusted for meeting the power balance criteria. The procedure of such adjustments, for ELD problems without and with considering transmission losses, are shown in Fig. 1(a) and (b), respectively. As Fig. 1 suggests, our adjustment procedure is slightly different for ELD problems with  and without PL . For ELD problems without PL , the difference ( N i=1 X[i] − PD ) is adjusted by subtracting it from a randomly chosen unit’s load such that the generating constraint is met. In contrast, the adjustment scheme for ELD problems with PL repeatedly attempts to distribute the  difference ( N X[i] − PL − PD ) evenly among all the generi=1 ating units such that their loads do not exceed the corresponding boundaries. Such even distribution helps to keep the change in PL minimum because of the small changes in Pi s and minimize the effect of adjustment thereby. Moreover, such iterated effort of even distribution helps to satisfy the power balance criteria

f (x) =

N 

Fi (Pi ) + k1 (

i=1

N  i=1

2

Pi − P D − P L ) + k 2 (

N  i=1

ViR ) (14)

where ViR , reckons the violation of the prohibited zone constraints for the individual i, can be defined as  1 if Pi violates the prohibited zones ViR = (15) 0 otherwise In the fitness function of (14), k1 and k2 are penalty factors associated with power balance and prohibited constraints, respectively. For the ELD problems without transmission loss and prohibited zone constraints, the setting k1 = k2 = 0 is most rational. And for the ELD problems with transmission loss and prohibited zone constraints these factors are tuned empirically and their values were set as k1 = 1 and k2 = 5 ∗ N in all cases. 4. Simulation results We verified the feasibility of the proposed DE method by applying it to different ELD problems. The studied problems can be classified into three categories depending on the type of their fuel cost functions and constraints. In the first category, there are two ELD problems with non-smooth cost functions due to considering valve-point loading. The first problem is a 13 generator system with original load demand 1800 MW. We call

1326

N. Noman, H. Iba / Electric Power Systems Research 78 (2008) 1322–1331

it SYS-1. The characteristics description of SYS-1 can be found in [8,28]. We also studied SYS-1 for a different load demand of 2520 MW [2,28]. The other system with valve-point loading had 40 generating units with higher nonlinearities. The description of the system, denoted as SYS-2, is available from [8] with the power demand of 10,500 MW. The second category contains an ELD problem with multiple-fuel effects. We call it SYS-3 for the input data and related constraints available from [4,9]. The third category includes two ELD problems that consider prohibited operating zones, ramp-rate limits as well as transmission losses. The first system consists of 6-thermal units with a load demand of 1263 MW. We denote it as SYS-4 and the characteristics information is available from [7,10]. The other problem in this category, named SYS-5, is a 15-thermal unit system. The load demand of SYS-5 is 2640 MW and the constraint data is available from [10]. The results obtained from the proposed DE algorithm are compared with the reported results by other methods, namely particle swarm optimization (PSO) [10], modified PSO (MPSO) [12], EP with sequential quadratic programming (EP-SQP) [2], PSO-SQP [2], GA [2], evolutionary strategy optimization (ESO) [7], differential evolution (DEvol) [21], self tuning hybrid DE (ST-HDE) [29], chaotic DE with SQP (DECSQP) [22], simulated annealing (SA) [2], GA-SA [2], CEP, FEP, MFEP, IFEP [8], improved EP (IEP) [30], hierarchical numerical method (HM) [31], modified hopfield neural network (MHNN) [9] and adaptive hopfield neural network (AHNN) [4]. Our code was written in Java and executed on a Intel(R)Core(TM)2 CPU 1.67 GHz personal computer with 1GB RAM. 4.1. Selection of DE parameters Though DE works with only a few control parameters, values of these parameters often significantly influence the quality of the solution obtained [16,32]. In general, selecting the most appropriate values for the parameter set is mostly a problem dependent task requiring previous experience. As described in Section 3, the three parameters upon which DEs performance is dependent are F , Cr and P. In classic DE, once selected all the three parameters remain fixed throughout the optimization procedure. The parameter F controls the speed and robustness of the search, i.e., a lower value of F increases the convergence rate but also increases the risk of getting stuck into a local optimum [27]. On the other hand, if F > 1.0 then solutions tend to be more time consuming and less reliable [14]. The parameter Cr which controls the crossover operation, can also be thought

of as a mutation rate, i.e., the probability that a variable will be inherited from the mutated individual. The role of Cr is to provide the means of exploiting decomposability, if it exists. Population size P can be a critical choice for the performance of DE because of its one-to-one reproduction strategy. Storn and Price [13] suggested a larger population size (between 5N and 10N) for DE. However, smaller population size can be useful given a maximum number of function evaluations[19]. Earlier studies on ELD using DE [21–23] did not investigate the effect of parameter settings for the algorithm, rather tuned the parameters for different problems. Therefore, due to lack of such studies with DE on the ELD problems, we made an extensive study for selecting the most suitable parameter set for DE for the problem in hand. In this work, we studied the SYS-1 (13 generating unit) with load demand PD = 1800 MW, by varying the parameters within the suggested ranges. In order to select the most suitable { F , Cr } pair, we fixed the P = 10N = 130 and experimented by varying F ∈ [0.1, 2.0] and Cr ∈ [0.1, 0.9] with a step size of 0.1 for each. To assure convergence we allowed enough generations (MAXGEN = 1000) in every experimental run. In every setting, 100 random trails were repeated and some selected results are presented in Table 1 in the form of Avg(Std). The results presented in Table 1 suggest that for most of the settings of Cr and F, DE can exhibit a performance better than many other methods (as shown in Table 3) and also for other settings the performance is reasonable. However, according to Table 1, the best settings are F = 0.1 and Cr = 0.2. Though some parameter studies on DE [33] have suggested that the setting of F < 0.4 is not useful [14], there are also counter example of success using F ranging from 0.0001 to 0.4, with F = 0.2 often proving effective [34]. On the other hand, optimizing an extensive test beds Storn and Price [13] suggested the most appropriate value for crossover probability to be 0.0 ≤ Cr ≤ 0.2 or 0.9 ≤ Cr ≤ 1.0 which agree with our settings. However, the suggested settings for the DE parameters may appear to be atypical but not completely aberrant. Moreover, successes on other ELD problems with this setting justify its effectiveness, as will be shown in subsequent case studies. The other parameter of DE needs to be tuned is population size. We experimented by setting P = qN where q = {1, 2, . . . , 10}, with a couple of settings for F and Cr , for the same problem of SYS-1. In order to make the comparison fair enough, we used equal number of fitness evaluations FE = 130,000 in each trial. So the maximum number of generations (MAXGEN) in each experiment was set to FE/P. In each experiment we

Table 1 Effect of F and Cr on DE-performance (SYS-1 with PD = 1800 MW) F

0.1 0.2 0.3 0.5 0.7 0.9

Cr 0.1

0.2

0.3

0.5

0.7

0.9

17966.84 (1.421) 17968.87 (0.769) 17970.14 (1.403) 17970.74 (1.537) 17979.05 (3.040) 17975.20 (2.403)

17965.91 (1.016) 17969.50 (1.319) 17974.84 (3.019) 17980.22 (3.988) 18016.42 (15.531) 17992.09 (8.114)

17967.49 (2.374) 17974.85 (3.124) 17993.03 (12.030) 18021.09 (17.156) 18089.39 (20.981) 18035.66 (15.587)

17974.69 (2.782) 18010.11 (22.159) 18080.59 (32.654) 18134.42 (25.009) 18215.88 (30.010) 18135.50 (24.851)

17974.99 (2.949) 18052.90 (50.884) 18180.01 (38.026) 18242.05 (36.479) 18327.14 (39.666) 18238.61 (32.708)

17999.59 (31.581) 18002.96 (34.054) 18033.13 (71.538) 18363.89 (36.946) 18437.52 (38.301) 18364.53 (34.895)

N. Noman, H. Iba / Electric Power Systems Research 78 (2008) 1322–1331

1327

Table 2 Effect of population size on DE-performance (SYS-1 with PD = 1800 MW) FE

130,000 130,000 130,000 130,000 130,000 130,000

MAXGEN

10,000 3,333 2,000 1,667 1,250 1,000

P

13 39 65 78 104 130

CR = 0.1, AF = 0.1

CR = 0.2, AF = 0.1

CR = 0.3, AF = 0.1

Average

S.D.

Average

S.D.

Average

S.D.

17974.95 17965.39 17964.42 17964.73 17965.55 17966.84

15.409 2.296 1.127 1.406 1.447 1.421

17976.58 17966.81 17964.65 17964.26 17964.95 17965.91

16.939 3.706 2.193 1.013 1.153 1.016

17981.88 17968.97 17967.42 17966.00 17966.74 17967.49

20.358 4.879 4.527 3.813 3.336 2.374

performed 100 repeated trials and some selected results are presented in Table 2. According to Table 2, the best result is obtained for the setting P = 6N, AF = 0.1, Cr = 0.2, though some other neighboring settings gave very similar results. Nevertheless, we chose this setting for all other experiments in this work. 4.2. ELD problems with cost functions considering valve-point effects In order to verify the effectiveness of DE in solving ELD problems, we compared the results obtained for SYS-1 with that reported using other algorithms. The setup for the algorithm was as follows: P = 6N = 78, AF = 0.1, Cr = 0.2 and MAXGEN = 1200. The presented results are from 100 repeated trials. We compared the proposed method with different types of EPs, namely CEP, FEP, MFEP and IFEP, from [8], EP-SQP, PSO, PSO-SQP from [2], DEC-SQP [22] and HDE, ST-HDE from [29] in terms of average cost, maximum cost and minimum cost, in Table 3. From Table 3 it is very evident that DE not only has found the highest quality results among the all algorithms compared, but also possesses the highest probability of finding the better solution for the problem under consideration. Moreover, the average execution time for this experiment (1050.8 ms) suggests that DE is capable to reach the solution at a very high speed.

Table 3 Comparison in terms of fuel cost for SYS-1 with PD = 1800 MW Method

Minimum cost

Maximum cost

Average cost

CEP FEP MFEP IFEP EP-SQP PSO PSO-SQP DEC-SQP HDE ST-HDE DE

18048.21 18018.00 18028.09 17994.07 17991.03 18030.72 17969.93 17963.94 17975.73 17963.89 17963.83

18404.04 18453.82 18416.89 18267.42 – – – 17984.81 – – 17975.36

18190.32 18200.79 18192.00 18127.06 18106.93 18205.78 18029.99 17973.13 18134.80 18046.38 17965.48

In order to further investigate the robustness of the proposed method we experimented with SYS-1 with a different load demand of 2520 MW. Experimental setup was exactly the same as in previous case except MAXGEN = 1000 was used. Over 100 repeated trials, the proposed DE algorithm was successful to achieve an average cost 24169.9177 with standard deviation 4.45E −5. The minimum and maximum costs among these 100 trials were 24169.9177 and 24169.9180, respectively. Table 4 compares the best result with those from other algorithms as reported in [2] and [7]. From the comparative results

Table 4 Comparison in terms of fuel cost for SYS-1 with PD = 2520 MW GU

GA

SA

GA-SA

EP-SQP

PSO-SQP

ESO

DE

z1 z2 z3 a1 a2 a3 a4 b1 b2 b3 c1 c2 c3 TP TC

628.32 356.49 359.43 159.73 109.86 159.73 159.63 159.73 159.73 77.31 75.00 60.00 55.00 2519.96 24398.23

668.4 359.78 358.2 104.28 60.36 110.64 162.12 163.03 161.52 117.09 75.00 60.00 119.58 2520.00 24970.91

628.23 299.22 299.17 159.12 159.95 158.85 157.26 159.93 159.86 110.78 75.00 60.00 92.62 2519.99 24275.71

628.3136 299.1715 299.0474 159.6399 159.656 158.4831 159.6749 159.7265 159.6653 114.0334 75.0000 60.0000 87.5884 2520.00 24266.44

628.3205 299.0524 298.9681 159.468 159.1429 159.2724 159.5371 158.8522 159.7845 110.9618 75.0000 60.0000 91.6401 2520.00 24261.05

628.3164 299.1971 299.1940 159.7323 159.7314 159.7309 159.7323 159.7281 159.7316 112.2280 112.6778 55.0000 55.0000 2520.00 24179.59

628.3185 299.1993 299.1993 159.7331 159.7331 159.7331 159.7331 159.7331 159.7331 77.3999 77.3999 92.3999 87.6845 2520.00 24169.9177

GU: generating unit; TP: total power; TC: total cost.

1328

N. Noman, H. Iba / Electric Power Systems Research 78 (2008) 1322–1331

Table 5 Comparison in terms of fuel cost for SYS-2 with PD = 10, 500 MW

Table 6 Estimated fuel cost for SYS-3 with different load demands

Method

Minimum cost

Maximum cost

Average cost

Demand (PD )

Minimum cost

Maximum cost

Average cost

CEP FEP MFEP IFEP EP-SQP PSO PSO-SQP MPSO DEvol DEC-SQP HDE ST-HDE DE

123488.29 122679.71 122647.57 122624.35 122323.97 123930.45 122094.67 122252.27 121412.91 121741.98 121813.26 121698.51 121416.29

126902.89 127245.59 124356.47 125740.63 – – – – 121464.40 122839.29 – – 121431.47

124793.48 124119.37 123489.74 123382.00 122379.63 124154.49 122245.25 – 121430.00 122295.13 122705.66 122304.30 121422.72

2400 2500 2600 2700

481.723 526.239 574.381 623.809

481.723 526.239 574.547 623.809

481.723 526.239 574.385 623.809

of Table 4, it is evident that the worst cost found by DE for this problem was better than the best cost by any other method. Therefore, from the experimental results of SYS-1, it can be stated that the proposed DE algorithm was capable of producing the minimum cost solution for two different power demands. Our second case study, considering the valve-point effects, was a 40-unit system. The DE algorithm with the same parameter settings (except P = 6N = 240 and MAXGEN = 1000) was applied to solve SYS-2. Table 5 presents the minimum cost, maximum cost and average cost produced by the proposed algo-

rithm comparing with the other reported results by different EPs, EP-SQP, PSO, MPSO, PSO-SQP, DEC-SQP, DEvol, HDE and ST-HDE. The average computational time required by DE to obtain the above results was 72.94 s. Comparing the results of Table 5, it can be found that the worst case result (maximum cost) by DE is even better than the best case result (minimum cost) by any other method except that of DEvol. Nevertheless, the average generation cost by DE was much better that that of DEvol, though DEvol produced slightly better minimum generation cost. Moreover, using the proposed algorithm and a very simple tuning of populations size to P = 50 and same number of fitness evaluation, produced the following results: min cost: 121412.6, max cost: 121421.2 and avg cost: 121414.5. In contrast, DEvol generated the above result by tuning all DE parameters even the DE mutation strategy. Consequently, the superiority of the proposed method compared to the existing techniques is undeniable. Moreover, these results from SYS2 once again demonstrate the appropriateness of the chosen parameter values for DE.

Table 7 Comparison of optimization methods for SYS-3 with PD = 2400 MW and 2500 MW S

U

PD = 2400 1

2

3

TP TC PD = 2500 1

2

3

TP TC

HM

MHNN

AHNN

IEP

MPSO

DE

F

GEN

F

GEN

F

GEN

F

GEN

F

GEN

F

GEN

1 2 3 4 5 6 7 8 9 10

1 1 1 3 1 1 1 3 1 1

193.2 204.1 259.1 234.3 249.0 195.5 260.1 234.3 325.3 246.3 2401.2 488.500

1 1 1 2 1 3 1 3 1 1

192.7 203.8 259.1 195.1 248.7 234.2 260.3 234.2 324.7 246.8 2399.8 487.870

1 1 1 3 1 1 1 3 1 1

189.1 202.0 254.0 233.0 241.7 233.0 254.1 232.9 320.0 240.3 2400.0 481.700

1 1 1 3 1 3 1 3 1 1

190.9 202.3 253.9 233.9 243.8 235.0 253.2 232.8 317.2 237.0 2400.0 481.779

1 1 1 3 1 3 1 3 1 1

189.7 202.3 253.9 233.0 241.8 233.0 253.3 233.0 320.4 239.4 2400.0 481.723

1 1 1 3 1 3 1 3 1 1

189.74 202.34 253.89 233.05 241.83 233.05 253.28 233.05 320.38 239.40 2400.0 481.723

1 2 3 4 5 6 7 8 9 10

2 1 1 3 1 3 1 3 1 1

206.6 206.5 265.9 236.0 258.2 236.0 269.0 236.0 331.6 255.2 2501.1 526.700

2 1 1 3 1 3 1 3 1 1

206.1 206.3 265.7 235.7 258.2 235.9 269.1 235.9 331.2 255.7 2499.8 526.130

2 1 1 3 1 3 1 3 1 1

206.0 206.3 265.7 235.9 257.9 235.9 269.6 235.9 331.4 255.4 2500.0 526.230

2 1 1 3 1 1 1 3 1 1

203.1 207.2 266.9 234.6 259.9 236.8 270.8 234.4 331.4 254.9 2500.0 526.304

2 1 1 3 1 3 1 3 1 1

206.5 206.5 265.7 236.0 258.0 236.0 268.9 235.9 331.5 255.1 2500.0 526.239

2 1 1 3 1 3 1 3 1 1

206.52 206.46 265.74 235.95 258.01 235.95 268.87 235.96 331.49 255.06 2500.0 526.239

N. Noman, H. Iba / Electric Power Systems Research 78 (2008) 1322–1331

1329

Table 8 Comparison of optimization methods for SYS-3 with PD = 2600 MW and 2700 MW S

U

PD = 2600 1

2

3

TP TC PD = 2700 1

2

3

TP TC

HM

MHNN

AHNN

IEP

MPSO

DE

F

GEN

F

GEN

F

GEN

F

GEN

F

GEN

F

GEN

1 2 3 4 5 6 7 8 9 10

2 1 1 3 1 3 1 3 1 1

216.4 210.9 278.5 239.1 275.4 239.1 285.6 239.1 343.3 271.9 2600.0 574.030

2 1 1 3 1 3 1 3 1 1

215.3 210.6 278.9 238.9 275.7 239.1 286.2 239.1 343.5 272.6 2599.8 574.260

2 1 1 3 1 3 1 3 1 1

215.8 210.7 279.1 239.1 276.3 239.1 286.0 239.1 342.8 271.9 2600.0 574.370

2 1 1 3 1 1 1 3 1 1

213.0 211.3 283.1 239.2 279.3 239.5 283.1 239.2 340.5 271.9 2600.0 574.473

2 1 1 3 1 3 1 3 1 1

216.5 210.9 278.5 239.1 275.5 239.1 285.7 239.1 343.5 272.0 2600.0 574.381

2 1 1 3 1 3 1 3 1 1

216.55 210.91 278.54 239.10 275.52 239.10 285.71 239.10 343.48 271.99 2600.00 574.381

1 2 3 4 5 6 7 8 9 10

2 1 1 3 1 3 1 3 3 1

218.4 211.8 281.0 239.7 279.0 239.7 289.0 239.7 429.2 275.2 2702.2 625.180

2 1 3 3 1 3 1 3 1 1

224.5 215.0 291.8 242.2 293.3 242.2 303.1 242.2 355.7 289.5 2699.7 626.120

2 1 1 3 1 3 1 3 1 1

225.7 215.2 291.8 242.3 293.7 242.3 302.8 242.3 355.1 288.8 2700.0 626.240

2 1 1 3 1 1 1 3 3 1

219.5 211.4 279.7 240.3 276.5 239.9 289.0 241.3 425.1 277.2 2700.0 623.851

2 1 1 3 1 3 1 3 3 1

218.3 211.7 280.7 239.6 278.5 239.6 288.6 239.6 428.5 274.9 2700.0 623.809

2 1 1 3 1 3 1 3 3 1

218.24 211.66 280.72 239.63 278.50 239.63 288.59 239.63 428.52 274.87 2700.0 623.809

4.3. ELD problems with cost functions considering multiple fuels The proposed DE algorithm was also verified by applying it to a multiple fuel ELD problem SYS-3. The system consists of 10 generating units where the fuel cost functions of the generators are represented as the piecewise quadratic functions. In this problem, generators use three types of fuels: Type 1, 2 and 3. The whole system is organized in a hierarchical structure consisting of three subsystems, where the first subsystem consists of four units and each of the other two subsystems consists of three units. The system characteristics data, used in our study, was the same as in [4,9]. We studied the system by varying the load demand from 2400 MW to 2700 MW with 100 MW increments. We used the same set of parameters for DE as before except P = 6N = 60 and GENMAX = 200 was used. Again, 100 repeated trials were performed for each load demand and the summary of these experiments are presented in Table 6. The average execution time required by our method for this problem was 83.1 ms. From the presented results it seemed that the performance of the algorithm did not vary at all over different random trials which points out the consistent performance of the algorithm for multi-fuel ELD problems. The best results from our algorithm was also compared with those of HM [31], IEP [30], MHNN [9], AHNN [4] and MPSO [12] and presented in Tables 7 and 8. Consulting Tables 7 and 8, it can be found that AHNN, IEP, MPSO and DE always satisfied all of the equality and inequality

constraints whereas HM and MHNN could not satisfy the power balance constraint exactly. From the presented results it can be found that DE has always found better results than IEP. Compared to HM, DE always produced cheaper cost solution except for PD = 2600 MW in which case the total cost by DE was slightly higher. The best solutions given by MHNN could not satisfy the power demand exactly, though, the minimum cost was slightly better than DE in case of PD = 2500 and 2600. AHNN was successful to present better solution than DE in case of PD = 2400, 2500 and 2600 but the difference in total cost was in second or third digit after decimal point. On the other hand, DE found a solution for PD = 2700 case that was much better than that of AHNN. When compared to MPSO, almost exact solution was produced by DE in every case, in terms of fuel type for each unit, generating power of each unit and the total generation cost. 4.4. ELD problems with transmission loss and other constraints In this section we evaluated our algorithm for two ELD problems considering TL , prohibited operating zones and ramp-rate limit constraints. To solve the first system containing six thermal units (SYS-4), P = 6N = 36 and MAXGEN = 100 were used. All other parameters were the same as before. The computational results are shown in Table 9 comparing with other methods. From Table 9, it is observed that the solution found by DE is the best known result for the problem. Though the generation

1330

N. Noman, H. Iba / Electric Power Systems Research 78 (2008) 1322–1331

Table 9 Comparison of different algorithms for SYS-4 Unit power (MW)

GA

PSO

ESO a

DE

P1 P2 P3 P4 P5 P6 Total power PL Total cost Maximum cost Average cost

474.807 178.636 262.209 134.283 151.904 74.181 1276.030 13.022 15459.000 15524.000 15469.000

447.497 173.322 263.475 139.059 165.476 87.128 1275.960 12.958 15450.000 15492.000 15454.000

451.560 173.440 263.990 147.460 164.680 71.320 1272.460 12.824 15407.527 15470.000 15430.000

447.744 173.407 263.411 139.076 165.364 86.944 1275.947 12.957 15449.766 15449.874 15449.777

a

In case of ESO power balance constraint is not satisfied.

cost found by ESO is the minimum, the power balance constraint is violated by −3.374 MW which makes the solution defective. The average CPU time required by our algorithm to find such results was 33.5 ms. Not only the best result, but also the average result found by DE is better than any other method which suggests the reliability of the algorithm. Furthermore, such success of the method in an ELD problem of different characteristics justifies the robustness of the chosen parameters for DE. Our last study involves a 15-unit system (SYS-5) in the third category. For economic dispatch of SYS-5, we employed our algorithm with the same setting as before, i.e., F = 0.1, Cr = 0.2, P = 6N = 90, MAXGEN = 500. The proposed method took 1162.13 ms in an average to find the solution for the problem. Table 10 compares the best solution of the proposed method with that of ESO, PSO and GA. From the results presented in Table 10 it is obvious that the worst result (max cost) obtained by DE is better than the best result (min cost) found by any other method. And the best result found by DE is the best known soluTable 10 Comparison of different algorithms for SYS-5 Unit power (MW)

GA

PSO

ESO

DE

P1 P2 P3 P4 P5 P6 P7 P8 P9 P10 P11 P12 P13 P14 P15 Total power PL  N P − (PD + PL ) i=1 i Total cost Maximum cost Average cost

415.31 359.72 104.42 74.98 380.28 426.79 341.32 124.79 133.14 89.26 60.06 50.00 38.77 41.94 22.64 2663.44 38.39 −4.96 33063.54 33337.00 33228.00

439.12 407.97 119.63 129.99 151.07 459.99 425.56 98.56 113.49 101.11 33.91 79.96 25.00 41.41 35.61 2662.40 32.43 −0.01 32858.00 33331.00 33039.00

455.00 380.00 130.00 130.00 170.00 460.00 430.00 69.21 56.60 159.19 80.00 80.00 25.00 15.00 15.00 2655.01 30.30 −5.30 32640.86 32710.00 32620.00

454.997 419.997 129.997 129.998 269.917 459.990 429.995 60.007 25.001 63.111 79.973 79.983 25.001 15.001 15.000 2657.966 27.975 −0.007 32588.865 32641.419 32609.851

tion for the problem. The robustness of DE in finding the best solutions independent of ELD problem types, dimensions and constraints, suggests the suitability of the tuned parameters once again. 5. Conclusion Differential evolution (DE) has proven to be a high class technique for solving different types of real world problems from communication engineering to biotechnology. In this work we have investigated the potential of the algorithm in solving ELD problems studying different cases. An improved way of satisfying the power balance constraint has been also proposed. Other boundary constraints were satisfied using a reflection mechanism that is commonly used in constrained optimization with DE. First of all, we made a parameter study on DE to choose the best set of parameters which was fixed for the rest of the studies. Though, the chosen set of parameters for DE is atypical, the success of the algorithm justifies its suitability. For two ELD problems with valve-point effects our method found better solutions compared to what was known as best until now. The robustness of the method was verified by the change in load demands of the problems. In case of ELD problems with piecewise quadratic cost functions, DE exhibited an overall better performance compared to conventional numerical methods, conventional Hopfield neural network approach and evolutionary programming approach. Nevertheless, the performance of DE for these problems was very similar to MHNN and exactly the same as MPSO. Moreover, for the problems considering transmission loss, ramp-rate limits and prohibited operating zones, DE found solutions better than so far best known results by any other method in terms of cost, power loss and mean performance. And in all these experiments the average performance of DE was outstanding and the required time was very much suitable for online solving. Considering all these results of the study with ELD problems with different characteristics, dimensions, demands and constraints it can be concluded that DE-performs better, at least matching many of the previously proposed methods. Acknowledgments The authors would like to thank the associate editor and the anonymous reviewers for their constructive comments and criticisms that helped a lot to improve the quality of the paper. References [1] A.J. Wood, B.F. Wollenberg, Power Generation, Operation and Control, Wiley-Interscience, Newyork, 1996. [2] T.A.A. Victoire, A.E. Jeyakumar, Hybrid PSO-SQP for economic dispatch with valve-point effect, Electric Power Syst. Res. 71 (1) (2004) 51–59. [3] M.A. Abido, Multiobjective evolutionary algorithms for electric power dispatch problem, IEEE Trans. Evol. Comput. 10 (3) (2006) 315–329. [4] K.Y. Lee, A. Sode-Yome, J.H. Park, Adaptive hopfield neural networks for economic load dispatch, IEEE Trans. Power Syst. 13 (2) (1998) 519–526. [5] W.-M. Lin, F.-S. Cheng, M.-T. Tsay, An improved tabu search for economic dispatch with multiple minima, IEEE Trans. Power Syst. 17 (1) (2002) 108–112.

N. Noman, H. Iba / Electric Power Systems Research 78 (2008) 1322–1331 [6] S.O. Orero, M.R. Irving, Economic dispatch of generators with prohibited operating zones: a genetic algorithm approach, IEE Proc. Gen. Transm. Distrib. 143 (6) (1996) 529–534. [7] A. Pereira-Neto, C. Unsihuay, O.R. Saavedra, Efficient evolutionary strategy optimisation procedure to solve the nonconvex economic dispatch problem with generator constraints, IEE Proc. Gen. Transm. Distrib. 152 (5) (2005) 653–660. [8] N. Sinha, R. Chakrabarti, P.K. Chattopadhyay, Evolutionary programming techniques for economic load dispatch, IEEE Trans. Evol. Comput. 7 (1) (2003) 83–94. [9] J.H. Park, Y.S. Kim, I.K. Eom, K.Y. Lee, Economic load dispatch for piecewise quadratic cost function using Hopfield neural network, IEEE Trans. Power Syst. 8 (3) (1993) 1030–1038. [10] Z.-L. Gaing, Particle swarm optimization to solving the economic dispatch considering the generator constraints, IEEE Tran. Power Syst. 18 (3) (2003) 1187–1195. [11] T.A.A. Victoire, A.E. Jeyakumar, Deterministically guided PSO for dynamic dispatch considering valve-point effect, Elect. Power Syst. Res. 73 (3) (2005) 313–322. [12] J.-B. Park, K.-S. Lee, J.-R. Shin, K.Y. Lee, A particle swarm optimization for economic dispatch with non-smooth cost functions, IEEE Tran. Power Syst. 20 (1) (2005) 34–42. [13] R. Storn, K.V. Price, Differential evolution a simple and efficient heuristic for global optimization over continuous spaces, J. Global Optim. 11 (4) (1997) 341–359. [14] K.V. Price, R.M. Storn, J.A. Lampinen, Differential Evolution: A Practical Approach to Global Optimization, Springer, Berlin, Heidelberg, 2005. [15] T. Takahama, S. Sakai, Constrained optimization by the epsilon constrained differential evolution with gradient-based mutation and feasible elites, in: Proceedings of the 2006 IEEE Congress on Evolutionary Computation, 2006, pp. 308–315. [16] V.L. Huang, A.K. Qin, P.N. Suganthan, Self-adaptive differential evolution algorithm for constrained real-parameter optimization, in: Proceedings of the 2006 IEEE Congress on Evolutionary Computation, 2006, pp. 324–331. [17] E. Mezura-Montes, J. Velazquez-Reyes, C.A.C. Coello, Modified differential evolution for constrained optimization, in: Proceedings of the 2006 IEEE Congress on Evolutionary Computation, 2006, pp. 332– 339. [18] M.F. Tasgetiren, P.N. Suganthan, A multi-populated differential evolution algorithm for solving constrained optimization problem, in: Proceedings of the 2006 IEEE Congress on Evolutionary Computation, 2006, pp. 340–347. [19] N. Noman, H. Iba, A new generation alternation model for differential evolution, in: Genetic and Evolutionary Computation Conference (GECCO 2006), July, 2006, pp. 1265–1272.

1331

[20] J. Vesterstrom, R. Thomsen, A comparative study of differential evolution, particle swarm optimization, and evolutionary algorithms on numerical benchmark problems, Congr. Evol. Comput. (2004) 1980–1987. [21] R.E. Perez-Guerrero, R.J. Cedenio-Maldonado, Economic power dispatch with non-smooth cost functions using differential evolution, in: Proceedings of the 37th Annual North American Power Symposium, October, 2005, pp. 183–190. [22] L. Dos, S. Coelho, V.C. Mariani, Combining of chaotic differential evolution and quadratic programming for economic dispatch optimization with valve-point effect, IEEE Trans. Power Syst. 21 (2) (2006) 989–996. [23] J.-P. Chiou, Variable scaling hybrid differential evolution for large-scale economic dispatch problems, Elect. Power Syst. Res. 77 (3/4) (2007) 212–218. [24] K.S. Swarup, P.R. Kumar, A new evolutionary computation technique for economic dispatch with security constraints, Int. J. Elect. Power Energy Syst. 28 (4) (2006) 273–283. [25] Z. Michalewicz, M. Schoenauer, Evolutionary algorithms for constrained parameter optimization problems, Evol. Comput. 4 (1) (1996) 1–32. [26] R. Storn, System design by constraint adaptation and differential evolution, IEEE Trans. Evol. Comput. 3 (1) (1999) 22–34. [27] S. Kukkonen, J. Lampinen, Constrained real-parameter optimization with generalized differential evolution, in: Proceedings of the 2006 IEEE Congress on Evolutionary Computation, 2006, pp. 207–214. [28] K.P. Wong, Y.W. Wong, Genetic and genetic/simulated-annealing approaches to economic dispatch, IEE Proc. Gen. Transm. Distrib. 141 (5) (1994) 507–513. [29] S.-K. Wang, J.-P. Chiou, C.-W. Liu, Non-smooth/nonconvex economic dispatch by a novel hybrid differential evolution algorithm, IET Gen. Trans. Distrib. 1 (5) (2007) 793–803. [30] Y.-M. Park, J.R. Won, J.B. Park, A new approach to economic load dispatch based on improved evolutionary programming, Int. J. Eng. Intell. Syst. Elect. Eng. Commun. 6 (2) (1998) 103–110. [31] C.E. Lin, G.L. Viviani, Hierarchical economic dispatch for piecewise quadratic cost functions, IEEE Trans. Power Appl. Syst. PAS-103 (6) (1984) 1170–1175. ˇ [32] J. Brest, S. Greiner, B. Boˇskovi´c, M. Mernik, V. Zumer, Self-adapting control parameters in differential evolution: a comparative study on numerical benchmark problems, IEEE Trans. Evol. Comput. 10 (6) (2006) 646–657. [33] R. G¨amperle, S.D. M¨uller, P. Koumoutsakos, A parameter study for differential evolution, in: WSEAS International Conference on Advances in Intelligent Systems, Fuzzy Systems, Evolutionary Computation, 2002, pp. 293–298. [34] N. Chakraborti, K. Misra, P. Bhatt, N. Barman, R. Prasad, Tight-binding calculations of Si-H clusters using genetic algorithms and related techniques: studies using differential evolution, J. Phase Equilibr. 22 (5) (2001) 525–530.