A new modified artificial bee colony algorithm for the economic dispatch problem

Energy Conversion and Management 89 (2015) 43–62

Contents lists available at ScienceDirect

Energy Conversion and Management journal homepage: www.elsevier.com/locate/enconman

A new modified artificial bee colony algorithm for the economic dispatch problem Dinu Calin Secui ⇑ Department of Energy Engineering, University of Oradea, Universitatii, 1, 410087 Oradea, Romania

a r t i c l e

i n f o

Article history: Received 2 June 2014 Accepted 11 September 2014

Keywords: Economic/emission dispatch Artificial bee colony algorithm Chaotic maps Valve-point effects Transmission losses

a b s t r a c t In this paper a new modified artificial bee colony algorithm (MABC) is proposed to solve the economic dispatch problem by taking into account the valve-point effects, the emission pollutions and various operating constraints of the generating units. The MABC algorithm introduces a new relation to update the solutions within the search space, in order to increase the algorithm ability to avoid premature convergence and to find stable and high quality solutions. Moreover, to strengthen the MABC algorithm performance, it is endowed with a chaotic sequence generated by both a cat map and a logistic map. The MABC algorithm behavior is investigated for several combinations resulting from three generating modalities of the chaotic sequences and two selection schemes of the solutions. The performance of the MABC variants is tested on four systems having six units, thirteen units, forty units and fifty-two thermal generating units. The comparison of the results shows that the MABC variants have a better performance than the classical ABC algorithm and other optimization techniques. Ó 2014 Elsevier Ltd. All rights reserved.

1. Introduction Economic dispatch (EcD) is an important optimizing problem in the electric power system operating. The aim of the EcD problem is to determine the power outputs of all generating units from a system, for a given time interval (one hour), in order to have minimum fuel cost, and to meet the required constrains. Another problem that plays an important role in the operating of the power systems is the emission dispatch (EmD). The EmD problem is similar to the EcD problem, but in this case it aims determining the power outputs of the generating units that minimize the emissions’ quantity – nitrogen oxides (NOx), carbon dioxide (COx), sulfur oxides (SOx) – released in the environment (caused by the burning of fossil fuels within thermal power plants), in certain imposed operating conditions. Searching for an optimal solution in the dispatching of the power outputs of the generators taking into consideration both fuel cost and the environmental emissions level, is known as Combined Economic and Emission Dispatch (CEED) problem. Basically, EcD and EmD problems aim to find a way to reduce operating costs and emissions levels for a power system without making investments in new equipment and/or technologies. The cost and emissions reduction is achieved

⇑ Tel.: +40 751912527; fax: +40 0259414626. E-mail address: [email protected] http://dx.doi.org/10.1016/j.enconman.2014.09.034 0196-8904/Ó 2014 Elsevier Ltd. All rights reserved.

by the optimal planning of the load dispatched among the system’s generating units. Solving EcD/EmD involves formulating a mathematical model of optimization and then selecting and/or developing an appropriate optimization technique. The simplest model for the EcD (EmD) problem is one in which the fuel cost (or the emissions level) of the generating units is defined by a quadratic function, and the constraints are limited to only two: the equality between the powers generated and demanded in the system, respectively the generating units operation between the minimum and maximum limits of power. However, in order to take into account more practical features of the operating units, the mathematical model has been improved by considering also the valve-point effects and by introducing the constraints regarding the ramp rate limits, prohibited zones of the units, the transmission losses. These additions determine a non-linear, non-smooth and non-continuous mathematical model of optimization. To solve the EcD/EmD or CEED problem, several methods, classic or based on artificial intelligence, have been used over time. Among the classical methods the following can be mentioned: quadratic programming [1], non-linear programming [2,3]. The classical methods may have difficulties in finding the optimal solution due to the non-linear and non-continuous nature of the optimization model [4,5]. As a result, several methods based on artificial intelligence were developed and applied for the above mentioned problems. These methods have the ability to identify higher quality solutions

44

D.C. Secui / Energy Conversion and Management 89 (2015) 43–62

[4,6,7] and can be grouped into three categories. The first category consists of methods applied in their original version, the second refers to methods/algorithms derived from the original version by changing the relations to updating the solutions or by adapting some of the algorithm’s parameters, and the third includes the hybrid methods. From the first category several methods to solve the EcD/EmD or CEED problems were applied, such as: Genetic Algorithm (GA) [8,9], Particle Swarm Optimization (PSO) [4,8], Differential Evolution (DE) [5,8], Virus Optimization Algorithm (VOA) [8], Harmony Search (HS) [8], Seeker Optimization Algorithm (SOA) [8], Ant Colony Optimization (ACO) [10], Artificial Bee Colony (ABC) [11], Bacterial Foraging Optimization (BFO) [12], Gravitational Search Algorithm (GSA) [8,13], Firefly Algorithm (FA) [8,14], Imperialist Competitive Algorithm (ICA) [15]. In the second category of methods the following are included: Modified HS algorithm (MHSA) [16], Self-Adaptive Real Coded Genetic Algorithm (SARGA) [17], PSO with Time Varying Acceleration Coefficients (PSO_TVAC) [18], Improved Bacterial Foraging Algorithm (IBFA) [19], Modified Bacterial Foraging Algorithm (MBFA) [20], Improved Harmony Search with Wavelet Mutation (IHSWM) [21], Modified ABC (MABC) [22], Incremental ABC (IABC), ABC with Dynamic Population (ABCDP), ABCDP with Local Search (ABCDP-LC) [23], Modified Shuffled Frog Leaping algorithm (MSFL) [24]. Among the hybrid methods a few are mentioned: Hybrid Differential Evolution algorithm (HDE) [25], Hybrid Genetic Algorithm based on differential evolution (HGA) [26], hybrid PSO with Sequential Quadratic Programming (PSO-SQP), Evolutionary Programming-SQP (EP-SQP) [27], hybrid GA–Pattern Search–SQP (GA–PS–SQP) [28], Fuzzy Adaptive PSO with Variable Differential Evolution (FAPSO-VDE) [29], Hybrid Multi-Agent based PSO (HMAPSO) [30], hybrid Chaotic PSO algorithm and SQP technique (CPSO– SQP) [31], hybrid that combines Shuffled frog leaping algorithm and Differential Evolution (SDE) [32], Hybrid swarm intelligence based Harmony Search (HHS) [33], hybrid Simulated Annealing with PSO (SA-PSO) [34]. Each of the methods above have their advantages and disadvantages regarding the quality and stability of the solution, the convergence of iterative process, the calculus efficiency or the simplicity of implementation. This paper proposes a method from the second category, being a modified ABC method. The ABC algorithm is a biological-inspired optimization technique, based on population which imitates the foraging behavior of the real honey bee, being introduced by Karaboga in 2005 [35]. Initially, the ABC algorithm has been used in solving unimodal and multi-modal numerical optimization problems on a limited set of test functions [35]. Due to its many positive features – simple in concept and easy to implement, flexible, the possibility of using chaotic maps and of developing hybrids from combinations with other techniques – the ABC algorithm has been successfully applied to the optimization of complex mathematical functions with or without constraints [36–43], or solving engineering problems, such as: the optimization of truss structures [44], engineering design optimization [45], automatic voltage regulator system [46], design and economic optimization of shell and tube heat exchangers [47], fault section estimation in power systems [48], optimal reactive power flow [49], economic dispatch [50]. Several studies show that the ABC algorithm performs better or just as good as other biological-inspired algorithms, such as genetic algorithm, differential evolution and particle swarm optimization [36–38]. Besides the already mentioned advantages, ABC algorithm may encounter a number of challenges (sub-optimal solutions, low robustness, slow convergence) in the optimization of composite functions, non-separable functions or problems that require

constraints [51]. Also, in [41] is concluded that ABC algorithm is good at exploration but poor at exploitation. To overcome these problems several variants of ABC algorithm have recently been proposed in order to enhance the performance of the original version. Some of the ways used to achieve this purpose are: the introduction of new relations (inspired by other algorithms such as DE or PSO) to update the solutions in the search space [39–43,51], ABC algorithm hybridization with other classical algorithms (derivative methods [52]) or metaheuristics (PSO, DE) [53,54], using chaotic maps [55], fuzzy theory [56], contraction strategies of the solutions’ search space [57], multi-objective ABC [58], self-adaptation of mutation step size [59]. In this paper, a modified ABC algorithm (called MABC) is proposed, to solve the economic dispatch problem, taking into account systems having various characteristics (different operating constraints, the emission pollutions, the valve-point effects). The MABC algorithm has the structure of the original ABC algorithm, but relies on a new relation for the update of the solutions in the search space, based on the relations presented in [42] and [51]. In [42] an ABC variant it is proposed – Global best artificial bee colony algorithm – that updates the solutions using the information from the neighborhood of the best solution. To improve the convergence rate of the basic ABC algorithm, in [51] two parameters are introduced to control the frequency and magnitude of the perturbation in the update process of the solutions. The MABC algorithm behavior is also investigated under two conditions: (i) incorporation of chaotic components (generated by cat and logistic maps) in the relation used for the solutions update; (ii) the introduction of two schemes for the solutions selection (disruptive selection and classical proportional selection). It must be mentioned that the cat map and disruptive potential selection scheme have not been investigated when solving the economic/emission dispatch problem. The performance of the MABC method was tested on systems with different characteristics and then compared with several other optimization techniques from the solutions stability and quality point of view. 2. Formulation of the problems (EcD, EmD, CEED) 2.1. The economic dispatch (EcD) problem The mathematic optimization model for the EcD problem includes three elements – optimizing variables, the objective function and constraint relations – listed below. The optimizing variables are the real output powers Pj, j = 1,2, . . . n of the generating units, presented as a solution vector P = [P1,P2, . . . Pj, . . . , Pn], where n is the total number of units from the analyzed system. The objective function F[P] aims to minimize the fuel total cost for the entire system:

min F ¼

n X F j ðPj Þ

ð1Þ

j¼1

where Fj(Pj) is the fuel cost, in $/h, for the jth unit, which is modeled by a quadratic polynomial function such as:

F j ðPj Þ ¼ cj P2j þ bj Pj þ aj

ð2Þ

If the valve-point effects are taken into consideration, then the cost function for unit j includes a sine factor, besides the quadratic polynomial function [19,27,60]:

F j ðPj Þ ¼ cj P2j þ bj Pj þ aj þ jej sinðf j ðPj;min  Pj ÞÞj

ð3Þ

45

D.C. Secui / Energy Conversion and Management 89 (2015) 43–62

where aj, bj and cj are the fuel cost coefficients of unit j, and ej and fj are the coefficients of unit j reflecting the valve-point effects. The EcD problem constrains can be expressed by the relations (4)–(9) [4,5]: (i) Minimum and maximum operating power limits are:

Pj;min  Pj  Pj;max ; j ¼ 1; 2; . . . ; n

min E ¼

Pj 

P0j

 DRj ;

if the output power increases

ð5Þ

if the output power decreases

ð6Þ

where P0j is the output power of unit j in the previous hour. DRj and URj are the down-ramp and up-ramp limits of the j unit (in MW/time-period). Combining the relations (4)–(6) the following output power (Pj) limits for the unit j are obtained:

POj;min  Pj  POj;max

ð7Þ

where POj,min = max(Pj,min, P0j  DRj) and POj,max = min(Pj,max, P0j + URj). (iii) Unit’s j operating zones, considering the prohibited zones are:

8 P 6 Pj 6 Plj;1 > > < j;min u Pj;z1 6 Pj 6 Plj;z ; z ¼ 2; 3; . . . ; NZ j > > : Pu 6 P 6 P j;NZ j

j

ð8Þ

j;max

where NZj is the number of prohibited zones of unit j. P lj;z and Puj;z are the lower and upper boundary of the z prohibited operating zone for the unit j. (iv) Real power balance constraint is:

PG  PL  PD ¼ 0

ð9Þ

where PD is the load demand of the system, in MW. PL represents the transmission line losses, in MW. The transmission line losses PL at the entire system level are quadratic functions in relation to variables Pj and they are calculated using the constant B coefficient formula:

PL ¼

n X n n X X Pi Bij Pj þ B0i Pi þ B00 i¼1 j¼1

ð10Þ

n X Pj

ð11Þ

j¼1

Relation (9) is satisfied with a certain tolerance (TOL), calculated by the relation (12):

TOL ¼ PG  PL  PD

ð12Þ

2.2. The emission dispatch (EmD) problem The EmD problem is similar to the EcD problem, but instead of minimizing the fuel cost it aims minimizing the quantity of emissions released in the environment. The quantity of emissions Ej(Pj), expressed in ton/h, released in the environment by the operation of

ð13Þ

j¼1

where aj, bj, cj, nj, kj are the emission coefficients of the generating unit j. The constraints of the EmD problem are the ones taken into consideration in the EcD problem in the relations (4)–(9). 2.3. The combined economic and emission dispatch (CEED) problem CEED is an optimizing problem with two objectives: minimizing the total fuel cost (F) and the quantity of polluting emissions (E). One way to solving the CEED problem is changing it into an optimizing problem with a single objective, using relation [16]:

FE ¼ w 

n X

Ej ðP j Þ þ ð1  wÞ 

j¼1

n X

ðFðPj;max Þ=EðPj;max ÞÞ  Ej ðPj Þ

ð14Þ

j¼1

where w is a weighting factor with values in the range of [0, 1]. The ratio F(Pj,max)/E(Pj,max) is a scaling factor specific to each generating unit, that makes possible summing up the Fj(Pj) and Ej(Pj) terms. If relation (14) considers w = 1 or w = 0, then two particular forms of the CEED problem are obtained: the EcD problem (when w = 1) and the EmD problem (when w = 0). The two values taken by w determine two points in the objectives space (Fuel cost, Emission), located on the Pareto-front. In order to obtain other points (solutions) located on the Pareto-front, the weight factor varies from 0 value to 1, using a step of 0.1 [16]. In practice, the system’s operator must chose a single operating point as a compromise between the two opposing objectives (Fuel cost and Emission). Following this the question for the solution ranking from Pareto front is put and determine the best compromise solution (BCS), that helps the system’s operator in taking a decision regarding the current functionality of the system (based on a priority list). Next two ways in choosing a BCS are shown, based on a fuzzy mechanism. Considering the Pareto-front was obtained following the previously described process (changing the w factor in relation (14)), the first method consists of the following steps [20]: (i) For each objective i and each solution k from the Paretofront a value using a linear membership function is calculated:

i¼1

where Bij is an element of the loss coefficient matrix of size nxn, B0i is i element of the loss coefficient vector of size n and B00 is the loss coefficient constant. The total generated power (PG) is:

PG ¼

n n X X Ej ðP j Þ ¼ ðcj þ bj Pj þ aj P2j þ nj expðkj Pj ÞÞ j¼1

ð4Þ

where Pj,min and Pj,max represent the minimum and the maximum operating power limits of unit j. (ii) Generator ramp-rate limits:

Pj  P0j þ URj ;

generating unit j, can be modeled through a quadratic function to which an exponential term is added. The total quantity of emissions (E) released in the environment by the operation of n units is described by relation [23]:

li;k ¼

8 > < > :

1 f i;min f i;k f i;max f i;min

f i;k 6 f i;min f i;min < f i;k < f i;max

0

f i;k P f i;max

ð15Þ

where fi,min and fi,max represent the minimum and maximum value of ith objective function. (ii) For each non-dominated solution k, the normalized linear membership function l⁄k is calculated as follows:

P

N li;k lk ¼ PP i¼1 PN k¼1 i¼1 li;k

ð16Þ

where P is the number of non-dominated solutions, and N is the number of the objective functions. (iii) The best compromise solution corresponds to the maximum value of lk : lmax = max(lk , k = 1,2, . . . ,P). The second method is applied in a similar manner with the first one, with the difference that step 2 is replace by the relation (17) for the quantity calculus l1k [19]:

46

D.C. Secui / Energy Conversion and Management 89 (2015) 43–62

l1k ¼ minfli;k ; i ¼ 1; 2; ::; Ng; k ¼ 1; 2; . . . ; P

ð17Þ

Then, BCS is chosen using a similar relation with the one in step 3:

lmax = max(l1k, k = 1,2, . . . ,P). Details regarding the non-dominated solutions and the Paretofront notions can be followed in [19,58]. 3. The original ABC algorithm In the ABC algorithm, the swarm includes three categories of bees (scout bees, employed bees and onlooker bees), each of them carrying on various activities to sustain the hive life. Scout bees perform a random search for food sources in a certain area. Employed bees have the role of exploiting the identified food sources and then sharing various pieces of information (such as orientation and quality of the food source) with other bees waiting in the hive (called onlooker bees) to make an as good as possible decision in choosing a food source. In the ABC algorithm the following associations are made: food sources positions represent possible solutions to the optimizing problem and the amount of food (nectar) in the sources represents the value of fitness function. The number of the employed bees (SN) is equal to the number of the onlooker bees and to the number of food sources from the search area. The ABC algorithm consists of the following phases [35,37]: initialization, employed bee phase, computation of the probabilities for onlooker bees, onlooker bee phase and scout bee phase, as described bellow. In the initialization process SN scout bees randomly search and identify SN food sources in a limited area, and then they become employed bees. The random positions of the food sources, which also represent the initial solutions, can be determined with relation (18):

xij ¼ xj;min þ ðxj;max  xj;min Þ randð0; 1Þ

ð18Þ

where xij is the jth component of i solution vector, at iteration zero; xj,min, xj,max are the minimum and maximum values of j components (j = 1,2 . . . n, where n represents the dimension of solution vector) for the Xi vector; i = 1,2 . . . SN; rand(0, 1) is uniformly distributed random number in the range of (0, 1). After initializing the SN solutions Xi = [xi1, xi2, . . . xij, . . . xin], i = 1,2 . . . SN, they are evaluated using the objective function F(Xi) = Fi. For a minimization problem the fitness function (fitness(Xi) = fitnessi) corresponding to Xi solution is determined by the relation:

fitnessi ¼



1=ð1 þ F i Þ;

if F i P 0

1 þ absðF i Þ; otherwise

ð19Þ

In the case of a maximization problem the fitness function is represented by the objective function (fitnessi = Fi). Also, in this phase, a counter Triali (Triali = 0, i = 1,2 . . . SN) is initialized, memorizing the number of trials in which a solution vector Xi could not be improved. After the initialization phase, an iterative search process is conducted to determine the new positions of food sources, including employed bees, onlooker bees and scout bees. In the ABC algorithm, the employed bees make a change (in their memory) on the positions of food sources, based on the local information. The process of updating the solutions is described by the relation:

v i;j ¼ xi;j þ Ui;j ðxi;j  xk;j Þ

ð20Þ

where, vij, xij represent the new, respective the old position (solution) of food source i for the jth component, i = 1,2 . . . SN, j = 1,2 . . . n; k = 1,2 . . . SN is a randomly chosen index and k – i; Uij is a random number in the range [1, 1].

The new position of the food source Vi = [vi1,vi2, . . . ,vij, . . . vin], i = 1,2 . . . SN, j = 1,2 . . . n, determined through the relation (20), is evaluated and then compared to the old position Xi. If the new food source is better than the old one, then the position of the new source is memorized by the bee. Otherwise, the old position is retained in the memory. After all employed bees SN conducted a search, they share the information regarding the quality of food sources and their positions, with the onlooker bees in a waggle dance. An onlooker bee evaluates the information from all employed bees, and then chooses a food source with a probability based on its quality. The onlooker bee makes the change, in its memory, of the food source’s position, described by relation (20). Then, as in the case of the employed bees, the quality of the new food source (solution) is evaluated and compared with the quality of the old food source. If the new food source is better than the old one, the bee memorizes the new position. Otherwise, the old position is maintained. An onlooker bee chooses a food source with a probability pi, determined by relation (21):

fitnessi pi ¼ PSN u¼1 fitnessu

ð21Þ

Both in the employed bees phase, and in the onlooker bees phase, if the new position Vi established by relation (20) is better than the old position Xi, then Triali counter is set to zero. Otherwise Triali counter is incremented by 1. In the ABC algorithm, if a food source position cannot be improved after a predetermined number of trials (called Limit), then the food source is considered abandoned. The abandoned source is replaced with a new food source discovered by a scout bee. The position of the new source is randomly generated using relation (18). Other details of the ABC algorithm can be seen in [35–38]. The ABC algorithm uses three control parameters: the number of employed bees or onlooker bees (SN), the value of the predetermined parameter (Limit) and the maximum number of iterations (kmax). The main steps of the algorithm are given below [36,38]: Initialization of the population; k=1; Repeat Employed bees phase; Computation of the probabilities for onlooker bees; Onlooker bees phase; Scout bees phase; k = k + 1; (k is the current number of the iterations) Until k = Maximum iterations number (kmax). Memorizing the best solution achieved.

4. Several strategies for solutions update and selection schemes 4.1. Strategies for solutions update In this section three strategies are presented, aiming to improve the performance of ABC algorithm, that are based on modifications of the solution update relations within employed bees phase and onlooker bees phase. These strategies were presented in [42,43,51] and underlay the modified ABC (MABC) proposed in this paper. The relations describing the strategies are:

xi;j ¼ xi;j þ Ui;j ðxr1;j  xr2;j Þ

ð22Þ

47

D.C. Secui / Energy Conversion and Management 89 (2015) 43–62

 xi;j ¼

ð23Þ

xi;j þ Ui;j ðxi;j  xk;j Þ; if r i;j < MR xi;j ;

otherwise

ð24Þ

where xbest,j is the jth component of the best solution; r1 and r2 are random indexes chosen from {1,2, . . . ,SN}, with r1 – r2 – i. MR is a control parameter, and ri,j is a random number in the range [0, 1]. The significance of other parameters was presented in Section 3. The first strategy is defined by (22) and is inspired by the original ABC algorithm [37], but the vectors (xi,j, xk,j) from the difference (xi,j  xk,j) are replaced by two other vectors (xr1,j, xr2,j) randomly selected from the current population. The second strategy is based on the relation (23), proposed in [42], in which the search for new solutions is done around the best solutions (xbest,j). The last strategy uses the relation (24), proposed in [51], in which some parts of the current solution (xij) are altered if the condition rij < MR is satisfied. Basically, the relations (22)–(24) replace the relation (20) of the original ABC algorithm. The algorithms that are based on the relations (23) and (24) may have better performance than the original ABC algorithm [42,43,51] or other population-based algorithms, such as differential evolution, dynamic multi-swarm particle swarm optimizer with local search, covariance matrix adaptation evolution strategy with increasing population size [51]. These aspects recommend the use of relations (22)–(24) in testing new combinations of strategies to update the solutions, which could lead to improve their quality and stability. 4.2. Schemes for solutions selection Besides the classical proportional scheme of solutions selection (called ‘‘Proportional’’) – commonly used the ABC algorithm, the potential of disruptive selection scheme (Disruptive) is analyzed. Proportional and Disruptive schemes are based on the roulette wheel selection method and the probabilistic selection (pi) is determined using the same relation in both schemes – relation (21). The difference between the two selection schemes consists in the computation of the fitness function’s value. Thus, in the Proportional scheme (commonly used in the ABC algorithm) the fitness function value is calculated with relation (19) and in the Disruptive scheme the fitness function value is calculated with relation (25) [63,64]:

fitnessi ¼ jF i  meanðF i Þj; i ¼ f1; 2; . . . ; SNg

ð25Þ

where Fi = Fi(Xi) is the objective function for Xi solution, and mean(Fi) is the average value of the Fi functions from the current iteration. Classical proportional selection (Proportional), presents some disadvantages, one of them being the possibility of creating a ‘‘super-individual’’ towards which the entire population converges gradually, leading to population’s diversity reduction and consequently the algorithm’s stagnation in a local minimum. By its definition, the disruptive selection tends to favor solutions that lead to extreme values (minimum/maximum) for the functions Fi, the other solutions having a lower probability of being selected in the next iteration. Also, in the case of onlooker bees, disruptive selection could improve the diversity of the population and maintain the balance between exploration and exploitation processes [63,64]. In this paper disruptive selection is compared to classical proportional selection of solutions (commonly used the ABC algorithm). 5. Modified ABC (MABC) algorithm

Table 1 Chaotic maps used. Chaotic maps

Equation

Cat map [61]

xk+1 = (xk + yk) mod 1; yk+1 = (xk + 2yk) mod 1; where x mod 1 = x  [x]; x0, y0 2 (0, 1), xk, yk 2 (0, 1). xk+1 = 4xk(1  xk); xk 2 (0, 1); x0 2 (0, 1) x0 – {0.25,0.5,0.75}.

Logistic map [62]

x0, y0 represent the initial values used to generate the chaotic sequences. {xk}k=1. . .1, {yk}sk=1. . .1 represent the numbers generated through a chaotic map, at the kth iteration.

characteristics were modified. These modifications were implemented in employed bees phase and the onlooker bees phase and are presented below. (i) relation (20) of the ABC algorithm for solutions update was replaced with relation (26):

xi;j ¼



xi;j þ ð2cv ij  1Þðxr1;j  xr2;j Þ;

if

ri;j < SF

Relation (26) follows combining two solutions update relations inspired by [37] and respectively [42], through the selection factor (SF). Since the first equation described by relation (26) performs a search around the current solution xij, and the second equation around the best solution of the previous iteration (xbest,j), relation (26) has the capacity of searching for solutions at least as good as the two equations taken individually. The new solution generated by relation (26) differs both from current solution (xij), and the best solution xbest,j. This way of defining relation (26) eliminates the drawback set by using relation (24) [51], in which the new solution generated may coincide with the current solution. (ii) in the MABC algorithm the values of the variables (cvij), from relation (26), are generated using two chaotic maps (Table 1). Chaotic maps have been successfully used in several algorithms, in order to improve the population diversity and to obtain better solutions in the end [42,43,55,57]. With the same purpose, this paper investigates the behavior of two chaotic maps (cat map [61] and logistic map [62]), whose characteristics are shown in Table 1. In Fig. 1 the frequencies of the values generated by the two maps are presented. Logic map 3500 Cat map

3000

Logistic map

2500 2000 1500 1000 500

xk

0

The MABC proposed algorithm has the general structure of the original ABC algorithm, presented in Section 3, but several of its

ð26Þ

xbest;j þ ð2cv ij  1Þðxr1;j  xr2;j Þ; otherwise

where xi,j represents the jth component of the current solution i; xbest,j is the jth component of the best solution; r1 and r2 are random indexes chosen from {1, 2, . . . , SN}, with r1 – r2 – i; cvij are chaotic variables whose values are generated by chaotic maps; ri,j is a uniform random number in the range [0, 1]; SF is a selection factor which takes values in range [0, 1].

Frequence

xi;j ¼ xbest;j þ Ui;j ðxr1;j  xr2;j Þ

0

0.2

0.4

0.6

0.8

1

Fig. 1. The frequency of the values generated by cat and logistic maps.

48

D.C. Secui / Energy Conversion and Management 89 (2015) 43–62

shows higher frequencies at the margins of the (0, 1) range, and cat map shows a relatively uniform frequency of values within the (0, 1) range. The logistic map was commonly used in solving various problems, having good results. The cat map has not been previously used for EcD/EmD problem solving therefore with this paper the aim is to investigate its potential in EcD/EmD problem by MABC algorithm. A detailed procedure for the MABC algorithm is presented below. {Phase 1: Initialization phase} 1: The parameters for the MABC algorithm are set: SN, n, Limit, SF, kmax; 2: Initialize: the iteration counter k (k = 0); the counters Triali = 0, i=1,2 . . . SN; the initial values for the chaotic sequences cx0, cy0 (for Cat map) or cx0 (for Logistic map); 3: Initialize the SN solutions (Xi, i = 1,2, . . . SN) using relation (18) (for each solutions the imposed constraints shown in Section 6, are checked); 4: Assess the value of the functions Fi=F(Xi) assigned to the initial solutions (Xi, i = 1,2, . . . SN), during iteration k = 0; 5: Select the index corresponding to the best solution Xbest, best = Arg(Min(Fi, i=1,2 . . . SN)) and calculate Fbest = F(Xbest). 6: Repeat {Phase 2: Employed bee phase} 7: For i = 1 To SN Do begin 8: Randomly select an index j to be updated, j 2 {1,2, . . . n}; 9: Randomly select the indexes r1 and r2 out of the current population, where r1 – r2 – i, r1, r2 2 {1,2, . . . SN}; 10: Update the values of the chaotic sequences for Cat map (cx, cy) or Logistic map (cx), using the equations displayed in Table 1; 11: Update the xij component for the Xi solution using relation (26); returns a new solution X0 i; 12: Check the imposed constraints on the problem, shown in Section 6; {Compare the old solution Xi with the new solution X0 i}; 13:If F(X0 i) < F(Xi) Then begin Xi = X0 i, Triali = 0 end Else Triali = Triali + 1; {Update the best solution Xbest and function Fbest=F(Xbest)}; 14:If F(Xi) < F(Xbest) Then begin Xbest = Xi; Fbest = F(Xbest); end (end For) {Phase 3: Computation of the probabilities for onlooker bees} 15: Compute the values for the fitnessi functions, i = 1,2 . . . SN using relation (25) (for a Disruptive selection scheme) or relation (19) (for a Proportional selection scheme); 16: For i = 1 To SN Do Compute the probabilities pi using relation (21); {Phase 4: Onlooker bee phase} 17: i = 1, h = 0 18: Repeat 19: If rand(0, 1) < pi Then 20: Randomly select an index j to be updated, j 2 {1,2, . . . n}; 21: Randomly select indexes r1 and r2 out of the current population, where r1 – r2 – i, r1, r2 2 {1,2, . . . SN}; 22: Update the values of the chaotic sequences for Cat map (cx, cy) or Logistic map (cx) using the equations shown in Table 1; 23: Update the xij component of the Xi solution using relation (26); returns a new solution X00 i; 24: Check the imposed constraints on the problem (Section 6);

{Compare the old solution Xi with the new solution X00 i}; 25: If F(X00 i) < F(Xi) Then begin Xi = X00 i, Triali = 0 end Else Triali = Triali + 1; {Update the best solution Xbest and function Fbest = F(Xbest)}; 26: If F(Xi) < F(Xbest) Then begin Xbest = Xi; Fbest = F(Xbest) end 27: h = h + 1 (end If) 28: Else begin i = i + 1; If i = SN + 1 Then i = 1 end 29: Until h 6 SN (end Repeat) {Phase 5: Scout bee phase} 30: The index of the abandoned Xix solution is identified from the condition ix = arg{max(Triali), i=1,2, . . . SN}. 31: If Trialix P Limit Then begin 32: Replace solution Xix with a new solution randomly generated with relation (18); 33: Check the imposed constraints on the problem (Section 6); 34: Triali=0; end (end If) 35: k = k + 1 {increment the iterations counter} 36: Until k 6 kmax (end Repeat) 37: Print the optimal solution Xbest and the best value of the objective function Fbest = F(Xbest).

In the MABC algorithm four parameters are used: the number of employed bees or onlooker bees (SN), the number of cycles required (Limit), the maximum iteration number (kmax) and selection factor (SF). Given the possibility of improving the ABC algorithm characteristics by introducing the relation (26) for the solutions update, the sequences generated by chaotic maps, and also the selection schemes included, the MABC variants could have high performance in terms of solutions’ quality and stability. 6. Constraints handling mechanism (CHM) In order to create the connection between the notations used inside the MABC algorithm (described in Section 5) and the specific notations used by the mathematical model of the stated issue (EcD, EmD or CEED), we assume that all SN solutions, generated during iteration k, are symbolized as follows: Xi  Pi or [xi1, xi2, . . . xij, . . . ,xin]  [Pi1, Pi2, . . . Pij, . . . ,Pin], i = 1,2 . . . SN; j = 1,2 . . . n. In order to meet the inequality and equality constraints, the following mechanism is applied: (i) In order to satisfy relation (4) (if the ramp-rate limits are not imposed) proceed as follows: if Pij < Pj,min then Pij = Pj,min; if Pij > Pj,max then Pij = Pj,max, otherwise Pij remains unchanged. (ii) In order to satisfy relation (7) proceed similarly to case (i): if Pij < POj,min then Pij = POj,min; if Pij > POj,max then Pij = POj,max, otherwise Pij remains unchanged. (iii) In order to satisfy relation (8) proceed as follows: if Pij belongs to prohibition zone z (having the limits Plj,z and Puj,z), then Pij is set to the closest limit (Plj,z or Puj,z); (iv) To ensure the relation (9) complies with the tolerance TOL, an equality constraint handling mechanism similar with the one proposed by [65] is used. 7. Case study The efficiency of the proposed method is shown through the analysis of four different test systems consisting of 6, 13, 40 and 52 generating units. Depending on their basic characteristics (load demand, power losses or a valve-point effect consideration) seven cases are analyzed, presented in Table 2.

49

D.C. Secui / Energy Conversion and Management 89 (2015) 43–62 Table 2 Basic features of the analyzed systems and the parameters setting for MABC variants. Parameters Cases

PD (MW)

kmax

Tested values for kmax

SN

Tested values for SN

Limit

SF

Valve-point effect (VP)

C1:6-units, with power losses (PL) C2:13-units, 1800 MW C3:13-units, 2520 MW C4:13-units, 1800 MW, with PL C5:13-units, 2520 MW with PL C6:52-units (4  13-units) C7:52-units (4  13-units) C8:40-units, (EcD problem) C9:40-units, (EmD problem)

1263 1800 2520 1800 2520 7200 10,080 10,500 10,500

400 18,000 15,000 20,000 15,000 25,000 20,000 15,000 2500

300, 400, 500 1.2  104, 1.5  104, 1.2  104, 1.5  104, 1.2  104, 1.5  104, 1.2  104, 1.5  104, 2.0  104, 2.2  104, 1.5  104, 1.8  104, 1.0  104, 1.2  104, 1500, 2000, 2500

5 12 12 20 12 20 15 15 10

4, 5, 6 10, 12, 14 10, 12, 14 10, 15, 20 10, 12, 15 16, 18, 20 10, 12, 15 10, 12, 15 6, 8, 10

0.02  kmax 0.02  kmax 0.02  kmax 0.02  kmax 0.02  kmax 0.02  kmax 0.02  kmax 0.02  kmax 0.02  kmax

0 0.9 0.4 0.9 0.5 0.9 0.7 0.2 0.2

No Yes Yes Yes Yes Yes Yes Yes Yes

1.8  104 1.8  104 2.0  104 1.8  104 2.5  104 2.0  104 1.5  104

The quality of the solutions is evaluated based on a sequence of 100 trials. For this sequence the values of the following items are kept: best total fuel cost F (B), average total fuel cost F (A), worst total fuel cost F (W), standard deviation (SD) and the success rate (SR). The SR is defined as follows [66]: SR = NS100/NTrials, where NTrials represents the total number of trials (NTrials = 100), and NS is the number of trials for which the condition (F(P)-B)8760 6 dF is satisfied. F(P) represents the fuel total cost F corresponding to a feasible solution P, obtained by a certain method, and B is the best fuel total cost F obtained by the same method. dF quantity is an assumed cost, considered insignificant for the analyzed systems. For this paper an assumption that dF = 1000 $/year (0.1141 $/h) was made. Setting the parameters of the MABC variants (SN, kmax, Limit and SF) was done using the MABC/D/Cat algorithm, following a procedure that is based on the design of experiments. This procedure allows testing the combinations between the different parameters and the values selected for them [8]. The mentioned parameters were set for each individual case (C1–C9) following three stages. In each stage the targeted parameters (SN and kmax for the first stage; Limit for the second stage and SF for the last stage) have been selected in accordance with the best value of item A, and for equal values of items A the tiebreaker was based on the item’s SD value. In the first stage the SN and kmax parameters were set, for which three test values were taken into consideration (pointed out in Table 2). These values were combined among themselves generating 32 pairs of (SN, kmax). For each pair (SN, kmax) 30 trials were performed with MABC/D/Cat algorithm, followed by the computation of items A and SD. In the second stage for the Limit parameter three test values were taken into consideration (0.01kmax, 0.02kmax, 0.03kmax) and the parameters SN and kmax were set with the values obtained during the first stage. With these settings 30 trials were performed for each value of the Limit parameter. In the last stage, to identify the best values of the SF, the MABC/D/Cat algorithm was run for 11 values in the [0, 1] range using a step of 0.1 (the parameters SN, kmax and Limit are set to the values established in the first two stages). For each value of the SF parameter the items A and SD (shown in Table 3) are calculated, this time performing 100 trials. For each case the best value of the SF parameter was recorded in Table 2. Analyzing the values of items A and SD from Table 3 for all the nine cases (C1–C9), the following can be mentioned: (i) for the cases C2, C4, C6, C7 and C8 SF parameter has an influence over A and SD; for the cases C3, C5 and C9 the influence of SF parameter over A and SD is small; for the case C1 the parameter SF has a low influence over the values of the mentioned items; (ii) for most cases (C2–C9) the best values of SF parameter (Table 2) differs from SF = 0 or SF = 1, in relation to items A and SD. This means that relation (26) (proposed for the

solutions update) is better than the individual relations that underlie it. The MABC algorithm was applied in six variants resulting from a combination of three modalities of generating random sequences (cat map, logistic map and ‘‘Random’’) with two selection schemes (Disruptive selection and Proportional selection). The third modality (Random) uses the random number Uij from relation (20). The MABC variants are presented in Table 4 and were applied considering that parameters kmax, SN, Limit and SF have the values shown in Table 2. All case studies were implemented in MathCAD, on a personal computer having a 1.79 GHz processor and 896 MB of RAM. 7.1. Test system 1: 6-unit with power losses (case C1) A six-unit test system that takes into consideration the transmission losses, ramp rate limits and prohibited operating zones of the units is analyzed in solving the EcD problem with the MABC variants. The data related to the cost coefficients (a, b, c), power operating limits, ramp-rate limits, prohibited operating zones of the units, and also the loss coefficient B are taken from [4], and corrected after [67] in regards to the B00 factor. The load demand is PD = 1263 MW. 7.1.1. Solution’s quality and robustness The best solutions obtained by six MABC variants are presented in Table 5. In order to analyze the robustness of MABC variants 100 independent trials were performed. The values obtained for items B, A, W and SD using MABC variants are presented in Table 6, case C1. Following the values for the B, A, W and SD items, a very good stability for all obtained solutions is noticed. Also satisfying the power balance is achieved with a very low tolerance (TOLM < 1010 MW, Table 5). However, from Table 6, it can be noted that the MABC/D/Cat variant has the best value for SD item, the other items (B, A, W) having very close values. 7.1.2. Comparing MABC with other methods Table 7 shows a comparison between the results reached by the MABC/D/Cat variant and the results obtained by other methods recently published. The methods presented in Table 7 are M = {MHSA [16], MSFL [24], IABC, IABC-LS [68], h-PSO [69], DHS [70]}. From Table 7 it is found that TOLM tolerance differs from one method to another. This difference arises because some references [16,24,68,69] used data from [4], but in [70] data from [4] were used, corrected on the basis of observations reached from [67]. In order to compare the MABC/D/Cat variant to methods M, MABC/ D/Cat is applied for each method M, imposing the same tolerance TOLM (indicated in Table 7). The results are then presented in Table 8. MHSA methods [16], IABC and IABC-LS [68] do not provide

D.C. Secui / Energy Conversion and Management 89 (2015) 43–62

15449.8995 1.80  107 17971.2979 3.36 24169.9177 1.95  106 18131.9616 16.58 24514.8756 7.92  107 71812.9504 9.78 96613.8557 13.10 121434.1985 20.16 176682.2647 3.84  105 15449.8995 2.21  107 17969.5088 4.33 24169.9177 6.31  107 18143.6280 16.28 24514.8756 8.37  107 71810.6142 9.80 96610.9146 11.55 121435.2361 21.52 176682.2647 8.81  105 15449.8995 5.67  107 17968.4123 4.51 24169.9177 6.83  107 18137.5887 12.52 24514.8756 3.50  107 71805.5548 11.07 96607.3674 9.52 121436.6106 21.58 176682.2647 1.13  104 15449.8995 4.06  107 17967.1120 4.32 24169.9177 8.19  107 18134.6542 10.43 24514.8756 1.09  106 71801.9949 11.60 96607.2097 9.00 121435.5695 19.57 176682.2648 2.12  104 15449.8995 3.63  107 17964.8207 2.82 24169.9177 1.10  106 18134.1873 10.32 24514.8756 2.34  106 71801.3616 10.69 96604.5720 5.89 121441.9779 22.26 176682.2650 7.35  104 15449.8995 2.37  106 17963.9194 0.90 24169.9177 8.53  106 18134.4222 9.60 24514.8756 3.10  106 71797.9524 10.88 96606.1833 5.67 121444.2366 23.34 176682.2662 3.42  103 15449.8995 3.07  105 17963.8293 2.26  104 24169.9194 2.33  103 18131.9616 6.20 24514.8756 2.28  105 71795.9105 11.25 96608.6093 5.74 121456.6659 28.16 176682.2706 1.61  102 15449.8995 3.73  105 17964.2100 1.65 24169.9444 3.35  102 18133.6854 3.34 24514.8763 1.60  103 71815.8549 15.04 96643.6510 28.26 121474.3266 29.74 176682.3025 9.61  102 C9:40-units, (EmD problem)

C8:40-units, (EcD problem)

C7: 52-units, 10,080 MW

C6: 52-units, 7200 MW

C5: 13-units with PL, 2520 MW

C4: 13-units with PL, 1800 MW

C3: 13-units, 2520 MW

C2: 13-units, 1800 MW

A ($/h) SD ($/h) A ($/h) SD ($/h) A ($/h) SD ($/h) A ($/h) SD ($/h) A ($/h) SD ($/h) A ($/h) SD ($/h) A ($/h) SD ($/h) A ($/h) SD ($/h) A ($/h) SD ($/h) C1: 6-units

Cases

The values in bold correspond to SF factor that determines the lowest average cost (A) or the lowest standard deviation (SD), for a particular case (C1–C9).

15449.8995 1.36  107 17969.7580 4.27 24169.9177 7.48  107 18142.1707 15.26 24514.8756 3.48  107 71809.2706 10.04 96609.6732 11.42 121434.2588 21.92 176682.2647 4.11  105

15449.8995 2.14  107 17971.0172 3.60 24169.9177 9.20  107 18133.6854 15.89 24514.8756 8.98  107 71811.6858 9.78 96611.9507 11.78 121431.7793 19.16 176682.2646 3.82  105

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 SF

Table 3 The values of items A and SD for the analyzed cases (C1–C9) and SF parameter, obtained by MABC/D/Cat algorithm.

Table 4 Different variants of the MABC algorithm implementation.

15449.8995 6.04  108 17971.8415 2.84 24169.9177 3.57  106 18134.4222 15.69 24514.8756 1.06  106 71813.2240 11.35 96611.1548 12.03 121433.6668 22.60 176682.2647 3.85  105

50

MABC variants

Selection scheme

Modalities of generating random sequences (GRS)

M1: M2: M3: M4: M5: M6: M7:

Disruptive (D) Disruptive (D) Disruptive (D) Proportional (P) Proportional (P) Proportional (P) Proportional (P)

Cat map Logistic map Random (Rand) Cat map Logistic map Random (Rand) Random (Rand)

MABC/D/Cat MABC/D/Log MABC/D/Rand MABC/P/Cat MABC/P/Log MABC/P/Rand ABC Classic

sufficient data for conducting the comparisons using statistical tests. However, a qualitative analysis based on items B, A, W and SD was achieved below (Table 7 and Table 8), finding the following: (i) the MABC/D/Cat has a very good stability (item SD has very low values), (ii) the MABC/D/Cat is superior to the five methods, in consideration to the values of items B, A, W and SD. Furthermore, the W item obtained by the MABC/D/Cat (WMABC/D/Cat) is smaller than the B item obtained with any other M method (BM) from Table 7 (WMABC/D/Cat < BM, M = {MHSA, MSFL, IABC, IABC-LS, h-PSO, DHS}); (iii) In order to compare the results obtained by MABC/D/Cat and MSFL/DHS/h-PSO, a parametric statistical test (t-test) was applied, assuming a normal distribution for data collections. In all cases the specific value for t-test (tvalue) is greater than its critical value (tcritic = 2.576 for a significance level of 1%): for MABC/D/Cat – MSFL comparison the obtained tvalue is 13.04, for MABC/D/Cat and h-PSO tvalue is 12.41, and for MABC/D/Cat and DHS tvalue is 13.23. Therefore the null hypothesis (there are not significant statistic differences between the results obtained using the two methods) is rejected for a significance level of 1%. The average computation time for the case C1 is 0.62 s. 7.2. Test system 2: 13-unit with and without power losses (cases C2–C5) The MABC variants are used to solve the EcD problem for a 13 thermal generating unit system with valve-point loading effects, with and without taking into consideration the transmission line losses. The data of the system (a, b, c, e, f coefficients and limits of generated powers) are taken from [60], and the loss coefficient B from [32]. Two values are considered for the load demand: PD = 1800 MW, respectively PD = 2520 MW. 7.2.1. Solution’s quality and robustness The best solutions obtained by MABC variants are shown in Table 9. The solutions are presented for MABC variants that have the best items A and SD. The values of items B, A, W, SD and SR for the cases C2–C5, obtained by MABC variants are presented in Table 6. It can be seen that, generally, the values of the items B, A, W, SD and SR are very close, except for the case C4 (1800 MW with PL) in which there are differences between items. Also, from Table 9 it can be noticed that tolerance values (TOLM) are very small. For cases C2, C3 and C5 the success rate is maximum (SR = 100%) for all MABC variants. This means that the cost of generating achieved through a MABC variant does not differ from the cost obtained by another variant with more than dF = 1000 $/year (assumed limit). Also, for C4, the best MABC variant (MABC/P/ Log) selected based on item A is the same as the one selected based on SR criterion. 7.2.2. Comparing MABC variants with other methods Tables 10 and 11 display the results obtained by various optimization techniques/methods for the case C2 (1800 MW without PL)

51

D.C. Secui / Energy Conversion and Management 89 (2015) 43–62 Table 5 The best solutions obtained through MABC variants, for the case C1 with 6-units (PD = 1263 MW, 100 trials). Algorithms Output

MABC/D/Cat

MABC/D/Log

MABC/D/Rand

MABC/P/Cat

MABC/P/Log

MABC/P/Rand

P1 (MW) P2 (MW) P3 (MW) P4 (MW) P5 (MW) P6 (MW) PG (MW) PL (MW) BM ($/h) TOLM (MW)

447.50320837 173.31776370 263.46315194 139.06507749 165.47352386 87.13552025 1275.95824564 12.95824564 15449.899525 7.208  1011

447.50446888 173.31850746 263.46250373 139.06536865 165.47329094 87.13410047 1275.95824017 12.95824017 15449.899525 3.183  1012

447.50397527 173.31785688 263.46298445 139.06521482 165.47347683 87.13473649 1275.95824476 12.95824476 15449.899525 3.865  1012

447.50362341 173.31843533 263.46344721 139.06581038 165.47291888 87.13399580 1275.95823104 12.95823104 15449.899525 2.092  1011

447.50316334 173.31876029 263.46322108 139.06543237 165.47283280 87.13482349 1275.95823337 12.95823337 15449.899525 7.094  1011

447.5039396 173.3187842 263.4622817 139.0650065 165.4731297 87.13509942 1275.95824111 12.95824111 15449.899525 1.590  1012

Table 6 The values of items B, A, W, SD and SR for different analyzed cases (C1–C9) and MABC variants. Cases

System

Algorithm

SD ($/h)

SR (%)

C1

6-units 1263 MW

MABC/D/Cat (M1) MABC/D/Log (M2) MABC/D/Rand (M3) MABC/P/Cat (M4) MABC/P/Log (M5) MABC/P/Rand (M6) ABC Classic (M7)

B ($/h) 15449.8995 15449.8995 15449.8995 15449.8995 15449.8995 15449.8995 15449.8998

A ($/h) 15449.8995 15449.8995 15449.8995 15449.8995 15449.8995 15449.8995 15449.9091

W ($/h) 15449.8995 15449.8995 15449.8995 15449.8995 15449.8995 15449.8995 15449.9496

6.04  108 1.47  107 1.30  107 6.39  107 2.34  106 1.27  106 1.01  102

100 100 100 100 100 100 100

C2

13-units 1800 MW

MABC/D/Cat (M1) MABC/D/Log (M2) MABC/D/Rand (M3) MABC/P/Cat (M4) MABC/P/Log (M5) MABC/P/Rand (M6) ABC Classic (M7)

17963.8292 17963.8292 17963.8292 17963.8292 17963.8292 17963.8292 17963.8292

17963.8293 17963.8294 17963.8296 17963.8296 17963.8296 17963.8301 17965.2669

17963.8305 17963.8317 17963.8372 17963.8374 17963.8339 17963.8560 17973.2375

2.26  104 3.84  104 9.53  104 9.24  104 6.80  104 3.21  103 3.20

100 100 100 100 100 100 82

C3

13-units 2520 MW

MABC/D/Cat (M1) MABC/D/Log (M2) MABC/D/Rand (M3) MABC/P/Cat (M4) MABC/P/Log (M5) MABC/P/Rand (M6) ABC Classic (M7)

24169.9177 24169.9177 24169.9177 24169.9177 24169.9177 24169.9177 24169.9177

24169.9177 24169.9177 24169.9177 24169.9177 24169.9177 24169.9177 24170.0202

24169.9177 24169.9177 24169.9177 24169.9177 24169.9177 24169.9177 24171.8733

6.31  107 5.77  107 6.33  107 9.19  107 7.13  107 1.00  106 2.49  101

100 100 100 100 100 100 78

C4

13-units 1800 MW with PL

MABC/D/Cat (M1) MABC/D/Log (M2) MABC/D/Rand (M3) MABC/P/Cat (M4) MABC/P/Log (M5) MABC/P/Rand (M6) ABC Classic (M7)

18127.7821 18127.7821 18127.7821 18127.7821 18127.7821 18127.7821 18127.7821

18131.9616 18130.4507 18131.7758 18130.0015 18129.7010 18132.3448 18134.1370

18158.8463 18134.3784 18165.8342 18149.8375 18134.3131 18165.8319 18150.4282

6.20 3.21 6.61 3.59 2.95 7.47 3.74

51 59 56 66 69 52 10

C5

13-units 2520 MW with PL

MABC/D/Cat (M1) MABC/D/Log (M2) MABC/D/Rand (M3) MABC/P/Cat (M4) MABC/P/Log (M5) MABC/P/Rand (M6) ABC Classic (M7)

24514.8756 24514.8756 24514.8756 24514.8756 24514.8756 24514.8756 24514.8756

24514.8756 24514.8756 24514.8756 24514.8756 24514.8756 24514.8756 24515.1370

24514.8756 24514.8756 24514.8756 24514.8756 24514.8756 24514.8756 24540.9439

3.50  107 5.09  107 1.05  106 1.85  106 1.80  106 1.37  106 2.61

100 100 100 100 100 100 99

C6

52-units 7200 MW

MABC/D/Cat (M1) MABC/D/Log (M2) MABC/D/Rand (M3) MABC/P/Cat (M4) MABC/P/Log (M5) MABC/P/Rand (M6) ABC Classic (M7)

71770.3207 71770.2823 71770.2892 71771.0620 71770.2738 71770.2871 71785.6965

71795.9105 71788.5302 71798.8556 71797.4292 71788.5467 71798.6050 71825.8070

71821.3679 71808.8760 71820.8519 71823.8676 71808.5627 71821.7501 71886.9981

11.25 10.52 10.09 11.26 10.61 11.99 18.61

1 8 2 0 16 1 0

C7

52-units 10,080 MW

MABC/D/Cat (M1) MABC/D/Log (M2) MABC/D/Rand (M3) MABC/P/Cat (M4) MABC/P/Log (M5) MABC/P/Rand (M6) ABC Classic (M7)

96601.8504 96601.8517 96601.8509 96601.8556 96601.8549 96601.8529 96601.8531

96604.5720 96605.2397 96606.1761 96607.3925 96605.3712 96605.3877 96632.1958

96629.7631 96628.3065 96630.9067 96631.9366 96629.3396 96628.6297 96706.5650

5.89 6.94 7.32 8.25 6.86 6.07 20.19

20 10 11 10 16 12 9

C8

40-units, 10500 MW (EcD problem, w = 1)

MABC/D/Cat (M1) MABC/D/Log (M2) MABC/D/Rand (M3)

121412.5409 121412.6013 121412.8738

121431.7793 121434.8539 121436.7163

121503.7552 121502.8264 121494.5488

19.16 22.37 19.60

2 2 0

(continued on next page)

52

D.C. Secui / Energy Conversion and Management 89 (2015) 43–62

Table 6 (continued) Cases

C9

System

40-units, 10,500 MW (EmD problem, w = 0)

Algorithm

B ($/h)

A ($/h)

W ($/h)

SD ($/h)

SR (%)

MABC/P/Cat (M4) MABC/P/Log (M5) MABC/P/Rand (M6) ABC Classic (M7)

121412.5796 121412.5918 121412.5390 121414.1619

121432.3026 121431.5763 121435.5004 121474.2957

121502.8316 121493.1885 121502.8826 121558.5607

18.54 18.16 22.78 27.00

1 1 1 0

MABC/D/Cat (M1) MABC/D/Log (M2) MABC/D/Rand (M3) MABC/P/Cat (M4) MABC/P/Log (M5) MABC/P/Rand (M6) ABC Classic (M7)

176682.2646 176682.2646 176682.2646 176682.2646 176682.2646 176682.2646 176682.2646

176682.2647 176682.2647 176682.2647 176682.2647 176682.2647 176682.2647 176682.2723

176682.2648 176682.2655 176682.2650 176682.2654 176682.2673 176682.2649 176682.4185

3.82  105 1.14  104 5.17  105 1.00  104 2.70  104 5.72  105 2.21  102

100 100 100 100 100 100 99

The values in bold represents the best performing MABC variants in relation to items A and SD, for the cases (C1–C9).

Table 7 The results obtained by different optimization techniques for the case C1: 6 units, PD = 1263 MW).

⁄

Method Output

MHSA [16]

MSFL [24]

IABC, IABC-LS [68]

h-PSO [69]

DHS [70]

P1 (MW) P2 (MW) P3 (MW) P4 (MW) P5 (MW) P6 (MW) BM ($/h) AM ($/h) WM ($/h) SDM ($/h) TOLM (MW)⁄

446.7273 173.4889 263.7642 138.8321 165.6496 86.9463 15442.52 15443.594 15444.25 – 0.5450142341

445.0140 175.5156 264.2614 137.3012 162.7899 90.4992 15442.5911 15447.60 15460.29 4.07 0.562145394

451.5204 172.1750 258.4186 140.6441 162.0797 90.3415 15441.108 15441.108 15441.108 – 0.69359931601

447.1045 173.1123 263.6503 139.1516 165.9343 86.5037 15443.1830 15443.2117 15443.3548 0.0291 0.49658658671

447.5285 173.2791 263.4772 139.0291 165.4864 87.1587 15449.8996 15449.9264 15449.9884 0.0204 0.00003723189

TOLM was calculated for each method M with the relation (12), using the solutions indicated in the mentioned references [16,24,68–70].

Table 8 The best solution obtained by MABC/D/Cat algorithm considering the same tolerance TOLM as the methods M shown in Table 7 (the case C1: 6 units, 100 trials). Algorithms Output

MABC/D/Cat for TOLMHSA [16]

MABC/D/Cat for TOLMSFL [24]

MABC/D/Cat for TOLIABC, IABC-LS [68]

MABC/D/Cat for TOLh-PSO [69]

MABC/D/Cat for TOLDHS [70]

P1 (MW) P2 (MW) P3 (MW) P4 (MW) P5 (MW) P6 (MW) BMABC/D/Cat ($/h) AMABC/D/Cat ($/h) WMABC/D/Cat ($/h) SDMABC/D/Cat ($/h) PL (MW) TOLMABC/D/Cat (MW)

447.3913243019 173.2346721727 263.3751759646 138.9724439133 165.3848697672 87.04459336632 15442.51985534 15442.51985542 15442.51985622 1.6152  107 12.94809372002 0.54501423400

447.3871680920 173.2321144459 263.3724113955 138.9692633733 165.3822824272 87.0423939604 15442.28790885 15442.28790891 15442.28790949 9.7913  108 12.9477790883 0.56214539394

447.3599709396 173.2120654879 263.3512011770 138.9470009752 165.3608248035 87.0206659678 15440.50812541 15440.50812547 15440.50812608 8.2998  108 12.9453286670 0.693599316010

447.4005436386 173.2417784405 263.3829635426 138.9809886902 165.3934157108 87.05271960856 15443.17554399 15443.17554406 15443.17554540 1.6725  107 12.94899621794 0.49658658668

447.5037628412 173.3179818758 263.4628031149 139.0650293142 165.4732723149 87.1353562784 15449.89902070 15449.89902075 15449.89902164 9.9131  108 12.9582429712 0.00003723182

and C3 (2520 MW without PL). Among these the following are mentioned: PSO varieties and hybrids (PSO [27,72], CLPSO [80], FAPSO (Fuzzy Adaptive PSO) [72], FAPSO-NM (FAPSO-with Nelder–Mead) [72], PSO-SQP [27], DSPSO-TSA [74], CPSO (Chaotic PSO), CPSO-SQP [31] etc); Harmony Search varieties and hybrids (HS, HHS(Hybrid HS) [6]; evolutionary algorithms (IGA_MU (Improved Genetic Algorithm with Multiplier Updating) [9], GA [74], EP-SQP [27], HGA (Hybrid GA) [26]; Ant Swarm Optimization varieties (ACO [10], CASO – Chaotic Ant Swarm Optimization [73], FCASO-SQP – Fuzzy adaptive CASO [73]); ABC algorithm [11]; Firefly algorithm (FA) [14]; Biogeography-Based Optimization (BBO) [76,78]. A statistical comparison between MABC variants and the most competitive algorithms presented in Table 10 or Table 11 is not possible to be conducted because it lacks certain information, such

as the value of items A and/or SD. For this reason only a qualitative analysis is performed, based on items B, A, W and SD. According with items (B, A, W and SD) it is noticed that the behavior of the MABC variants (Table 6, the case C2: 1800 MW without PL) is better than most of the algorithms shown in Table 10, and is just as good as FAMPSO [7] (algorithms FAPSOVDE [29] and MsEBBO [76] obtain high quality solutions (by item B), but do not respect the imposed constrains). Similarly for the case C3 (2520 MW without PL), where the behavior of MABC variants (Table 6, the case C3) is better than most of the algorithms shown in Table 11, and is just as good as DE [5], FAMPSO [7] or HHS [6]. Needed to be mentioned is that the DE, FAMPSO or HHS algorithms either have a standard deviation (SD) higher than 5.77  107 $/h (obtained by MABC/D/Log), or the SD value was not specified by the authors.

53

D.C. Secui / Energy Conversion and Management 89 (2015) 43–62 Table 9 The best solutions obtained through MABC variants for the cases C2–C5 with 13 units (100 trials). Algorithm Output

MABC/D/Cat (C2:1800 MW, without PL)

MABC/D/Log (C3:2520 MW, without PL)

MABC/P/Log (C4:1800 MW, with PL)

MABC/D/Cat (C5:2520 MW, with PL)

P1 (MW) P2 (MW) P3 (MW) P4 (MW) P5 (MW) P6 (MW) P7 (MW) P8 (MW) P9 (MW) P10 (MW) P11 (MW) P12 (MW) P13 (MW) PL (MW) SD ($/h) TOLM (MW) BM ($/h)

628.3185307138 149.5996501669 222.7490698600 109.8665500538 109.8665500567 109.8665500401 109.8665500472 109.8665500487 60.0000000000 40.0000000000 40.0000000000 55.0000000000 55.0000000000 – 2.26  104 9.87  107 17963.829200

628.3185307163 299.1993002422 299.1992999811 159.7331000954 159.7331001135 159.7331001020 159.7331000836 159.7331000765 159.7331001088 77.3999122700 77.3999120533 87.6845308618 92.3999123012 – 5.77  107 9.94  107 24169.917696

448.7989505145 299.1993059732 224.3994800693 109.8665500881 109.8665501037 109.8665502275 109.8665500941 109.8665500558 109.8665500980 40.0000000000 40.5518593925 55.0000000013 55.0000000000 22.14889661794 2.95 6.02  1011 18127.782085

628.3185307131 299.1993003418 299.1993003372 159.7331001093 159.7331001144 159.7331001111 159.7331001030 159.7331001114 159.7331001144 77.3999124208 113.1111504227 92.3999125176 92.3999124897 40.426619906476 3.50  107 2.91  1011 24514.875629

Table 10 Results of the comparison between MABC/D/Cat and other published algorithms (the case C2 with 13-units, PD = 1800 MW without PL). Items Algorithms

B ($/h)

A ($/h)

W ($/h)

SD ($/h)

DPSO [30] HMAPSO [30] NDS [71] SA [71] PSO_TVAC [18] NEW PSO [18] IFEP [60] DE [5] PSO-SQP [27] EP-SQP [27] PSO [27] ED-DE [81] ST-HDE [25] HDE [25] FAPSO [72] FAPSO-NM [72] PSO [72] FAMPSO [7] FAPSO-VDE [29] HGA [26] CASO [73] FCASO-SQP [73] GSA [13] CLPSO [80] SQP-CLPSO [80] IPSO [82] HHS [33] DEC-SQP [79] BFO [12] GA-PS-SQP [28] TSARGA [17] ABC [11] aBBOmDE [78] IHSWM [21] NPSO [75] MsEBBO/sin [76] MsEBBO [76] h-PSO [69] FA [14] SDE [32] MABC/D/Cat

17976.31 17969.31 17976.9512 – 17963.879 18120.594 17994.07 17963.83 17969.93 17991.03 18030.72 17963.86 17963.89 17975.73 17963.84 17963.84 18030.72 17963.82 17963.82

18084.99 17969.31 17976.9512 18299.2550 18154.562 18227.052 18127.06 17965.48 18029.99 18106.93 18205.78 17972.70 18046.38 18134.80 17969.9187 17963.9577 18205.9247 17963.83041 17963.82484

18310.43 17969.31 17976.9512 – 18358.310 18427.631 18267.42 17975.36 – – – 17975.89 – – 17976.35 17964.21 18401.35 17963.835 17963.832

– 0.00 0.0000 123.8335 – – – – – – – – – – – – – – –

17963.83 17965.15 17964.08 17963.84 17970.67 17973.12 17998.44 17963.8293 17963.9401 17974.48 17964.25 17963.94 17963.86 17963.8521 17963.83 17976.015 17963.8292

17988.04 18022.04 18001.96 18041.21 18019.41 18005.05 18176.95 17972.4822 17973.1339 17997.12 – 17974.31 17987.22 17967.3560 17976.475 – 17967.0705

– – – 18910.31 18071.57 18069.35 – – 17984.8105 18018.75 18,199 18264.23 17995.11 17975.0552 18041.345 – 17972.8105

– – – – 22.67055 23.81023 – 2.4185 1.9735 – – 3.18 – – 25.6079 – 4.1876

17963.8292 17963.8297 17963.83 17963.83 17963.82920

17964.0468 17965.2055 18029.16 – 17963.82933

17969.0323 17980.2030 18168.80 – 17963.83045

1.9215 4.3807 148.542 – 2.26  104

For C2 and C3 cases (1800 MW and 2520 MW without PL) the stability of solutions obtained by MABC variants was tested using a sequence of 100 trials (Table 6). To check if the stability is

maintained another nine sequences were subsequently added, each consisting in 100 trials. The results obtained with MABC/D/ Log (for the case C3) and MABC/D/Log (for the case C3) are shown in Fig. 2 respectively Fig. 3. A very good stability of the solutions (especially for C3 case, Fig. 3) can be noticed. Table 12 shows a comparison of results obtained by MABC/P/ Log, MABC/D/Log and other algorithms (STHDE [25], SDE [32], BBO, DE/BBO, ORCCRO [83]), for C4 and C5 cases (1800 MW and 2520 MW with PL). For PD = 2520 MW, the MABC/D/Cat algorithm is more stable than the algorithms previously presented (SD item has a very small value). It must be mentioned that the value of item BORCCRO = 24513.91 $/h, from [83], is obtained by working with a higher tolerance (|TOLORCCRO| = 0.98 MW > |TOLMABC/D/Cat| = 2.91  1011 MW). TOLORCCRO was calculated based on the solution presented in [83]. In the case of PD = 1800 MW, it can be noticed in Table 12 that MABC/P/Log algorithm obtains better values for B and A items, comparing to SDE algorithm [32] and therefore better quality solutions. The best solutions for cases C4 and C5 are shown in Table 9. The average computation time is 38.2 s. (for the case C2), 32.7 s. (for the case C3), 4.07 min. (for the case C4) and 1.96 min (for the case C5).

7.3. Test system 3: 52-unit without power losses (the cases C6 and C7) This case presents a system consisting of 52 units, resulting by multiplying by four the 13-units system described in paragraph 7.2 (test system 2). The characteristics of this system (cost coefficients, power operating limits) are the same as for the system with 13-units [60]. Two values are considered for the load demand: PD = 7200 MW, respectively PD = 10,080 MW. The transmission losses are neglected in this case.

7.3.1. Solution’s quality and robustness The values of items B, A, W, SD and SR obtained through MABC variants for cases C6 and C7 (52-units with 7200 MW and 10,080 MW) are shown in Table 6. In relation to item A (from Table 6) the best performing variants are MABC/D/Log (for the case C6) and MABC/D/Cat (for the case C7). Also, item B value is the approximately the same for all MABC variants (for both C6 and C7). Following the values for the B, A, W and SD items, a good stability of the obtained solutions using MABC/D/Log, respectively MABC/D/Cat, is noticed.

54

D.C. Secui / Energy Conversion and Management 89 (2015) 43–62

Table 11 Results of the comparison between MABC/D/Log and other published algorithms (the case C3 with 13-units, PD = 2520 MW without PL). Items Algorithms

B ($/h)

A ($/h)

W ($/h)

SD ($/h)

IGA_MU [9] DE [5] FAMPSO [7] FAPSO-NM [72] FAPSO [72] PSO [72] HHS [6] HS [6] ACO [10] PSO-SQP [27] EP-SQP [27] HGA [26] DSPSO-TSA [74] TSA [74] GA [74] PSO [74] SA-PSO [34] CPSO [31] CPSO-SQP [31] CASO [73] FCASO-SQP [73] SDE [32] SDE [77] MABC [22] MABC/D/Log

24169.9790 24169.9177 24169.9176 24169.92 24170.93 24262.73 24169.9 24208.7 24174.39 24261.05 24266.44 24169.92 24169.923 24171.211 24170.804 24170.167 24171.395 24211.56 24190.97 24212.93 24190.63 24169.92 24169.9176968035 24208.8330 24169.917696

24385.4113 – 24169.9176 24170.0017 24173.0069 24271.9231 24169.9 24323.2 24211.09 – – – 24173.137 24184.055 24188.394 24184.849 – – – – – – 24170.7459 – 24169.917698

24754.1450 24169.9180 24169.9176 24170.5 24176.4 24277.81 24169.9 24503.7 24243.90 – – – 24230.803 24392.203 24567.974 24377.890 – – – – – – 24171.4402 – 24169.917700

– 4.45  105 – – – – – – 21.10 – – – 7.72 41 59.53 38.86 – – – – – – – – 5.77  107

‘‘–’’ Data not available.

24169.917701

The parameters specific to each algorithm were obtained by performing several experimental tests. Parameters values which led to the best value for the item A are given in Table 13. In Table 14 we present the items B, W, A and SD obtained through the algorithms: PSO, HS, DE, ABC, MABC/D/Log and MABC/D/Cat, using a sequence of 100 trials. The results of applying the Mann–Whitney U statistic test (MW-U) for comparing MABC/D/Log algorithm for the case C6 (or MABC/D/Cat for the case C7) and competitor algorithms (PSO, HS, DE and ABC) are shown in Table 14. For both cases (C6 and C7) the test MW-U shows that MABC variants are superior to the competitor algorithms (p-values are less than the significance level of 1%). Details regarding the MW-U test are presented in Section 7.5.1. The best solutions obtained through the MABC/D/Log and MABC/D/ Cat algorithms are presented in Table 15. Also, from Table 14 it can be noticed that algorithms HS, PSO, ABC and DE cannot find a solution as good as that obtained by MABC variants (BMABC < BM, where M = {HS, PSO, ABC, DE}). The average computation time is 5.47 min. (for the case C6) and 3.51 min (for the case C7).

24169.917699

7.4. Test system 4: 40-unit without power losses (the cases C8 and C9)

17964.05

Cost F ($/h)

17964.00 17963.95 17963.90 17963.85 No. of trials

17963.80 0

100

200 300 400

500 600 700

800 900 1000

Fig. 2. The best total cost F obtained with MABC/D/Cat, in the case C2 (13-units 1800 MW without PL).

Cost F ($/h)

24169.917703

24169.917697 No. of trials

24169.917695 0

200

400

600

800

1000

Fig. 3. The best total cost F obtained with MABC/D/Log, in the case C3 (13-units 2520 MW without PL).

7.3.2. Comparing MABC variants to other methods The results obtained through the best performing MABC variants (MABC/D/Log and MABC/D/Cat) are compared with the results obtained through other well known algorithms (PSO, HS, DE and ABC), whose performance has been proven for various mathematical functions or technical systems.

In this case a large practical system with 40 thermal units is considered to solving the EcD problem combined with the EmD problem using the MABC variants. The fuel cost coefficients (a, b, c, e, f) for the system’s units are provided by [20,60,85]. The emissions coefficients (a, b, c, f, k) for the system’s units are provided by [20,23,85,86,88]. The transmission losses are neglected in this case. The load demand is PD = 10,500 MW. For this system three situations are analyzed: (i) minimizing only the fuel cost (w = 1); (ii) minimizing only the emissions (w = 0); (iii) minimizing both the fuel cost and the emissions (w varies between 0 and 1 using a step of 0.1). 7.4.1. Solution’s quality and robustness Table 6 (the cases C8 and C9) shows the results obtained using the MABC variants (M1–M6) and ABC Classic algorithm (M7) for

55

D.C. Secui / Energy Conversion and Management 89 (2015) 43–62 Table 12 Results of the comparison between MABC/D/Cat and other published algorithms (the cases C4 and C5, 13-units with PL). Items Algorithms

PD (MW)

B ($/h)

A ($/h)

W ($/h)

SD ($/h)

PL (MW)

SDE [32] STHDE [25] ORCCRO [83] BBO [83] DE/BBO [83] MABC/D/Cat SDE [32] MABC/P/Log

2520 2520 2520 2520 2520 2520 1800 1800

24514.88 24560.08 24513.91 24515.21 24514.97 24514.8756 18134.49 18127.7821

24516.31 – 24513.91 24515.32 24515.05 24514.8756 18138.56 18129.7010

– – 24513.91 24516.09 24515.98 24514.8756 – 18134.3131

– – – – – 3.50  107 – 2.95

40.43 – 39.43 – – 40.426 19.13 22.14

Table 13 The setting of parameters specific to algorithms ABC, PSO, DE and HS. Algorithm Parameters

ABC

PD (MW) kmax NP/Limit

7200 25,000 20/0.02

PSO 10,080 20,000 15/0.02

DE

7200 5000 100/–

10,080 20,000 15/–

7200 10,000 50/–

HS 10,080 3000 100/–

7200 25,000 20/–

10,080 20,000 15/–

For PSO: c1 = 1.5, c2 = 2, wmin = 0.3, wmax = 0.9 (for the cases C6 and C7). For DE: CR = 0.7, AF = 0.1 (for the case C6) and CR = 0.2, AF = 0.1 (for the case C7). For HS: HMCR = 0.99, PAR = 0.1, BW = 0.01 (for the case C6) and HMCR = 0.98, PAR = 0.8, BW = 0.01 (for the case C7).

Table 14 Results of comparison between MABC/D/Log (for C6) and MABC/D/Cat (for C7) with other well-known algorithms HS, PSO, ABC and DE (the system with 52- units, 100 trials). Systems

Algorithms Items

ABC

PSO

DE

HS

MABC/D/Log (for C6) MABC/D/Cat (for C7)

52-units 7200 MW (C6)

B ($/h) A ($/h) W ($/h) SD ($/h) p-value

71785.6965 71825.8070 71886.9981 18.61 0.000

71990.9487 72316.8523 72734.9164 162.2506 0.000

71807.6321 71889.5475 71997.5088 37.3469 0.000

71806.4399 71848.1929 71929.3871 28.0473 0.000

71770.2823 71788.5302 71808.8760 10.5221 –

52-units 10,080 MW (C7)

B ($/h) A ($/h) W ($/h) SD ($/h) p-value

96601.8531 96632.1958 96706.5650 20.19 0.000

96851.7648 97350.5705 97807.2548 192.4468 0.000

96609.7373 96619.0013 96633.6183 5.2046 0.000

96686.8987 96780.1785 96918.8149 74.5547 0.000

96601.8504 96604.5720 96629.7631 5.8908 –

Table 15 Best solutions obtained with MABC/D/Log and MABC/D/Cat (the system with 52-units). Algorithm Output (MW)

MABC/D/Log, PD = 7200 MW (for the case C6)

P1/P14/P27/P40 P2/P15/P28/P41 P3/P16/P29/P42 P4/P17/P30/P43 P5/P18/P31/P44 P6/P19/P32/P45 P7/P20/P33/P46 P8/P21/P34/P47 P9/P22/P35/P48 P10/P23/P36/P49 P11/P24/P37/P50 P12/P25/P38/P51 P13/P26/P39/P52 BM ($/h) TOLM (MW)

448.798950 299.143951 224.379100 109.866550 109.866550 109.866550 109.865942 109.866540 109.866547 40.000000 40.000000 55.000000 55.000000 71770.2823 2.01  1011

538.558741 297.746123 224.383194 109.866480 109.866527 109.866550 109.866543 109.866203 109.866546 40.000000 40.000000 55.000000 55.000000

448.798951 224.390149 224.379956 109.866548 60.000000 109.866532 109.866550 60.000000 109.866550 40.000000 40.000000 55.000000 55.000000

MABC/D/Cat, PD = 10,080 MW (for the case C7) 448.798950 299.180044 224.378884 109.866550 109.866549 109.866550 109.866550 109.866550 109.866549 40.000000 40.000000 55.000000 55.000000

two situations: considering only the total fuel cost (case C8), respectively considering only the total emissions (case C9). In relation to items A and SD, the best variants are MABC/P/Long and MABC/D/Cat (for case C8), respectively MABC/D/Cat and MABC/D/ Rand (for case C9). In C8 case the MABC variants achieve relatively

628.318531 299.199296 299.199300 159.733100 159.733100 109.866550 159.733100 159.733100 159.733100 77.399909 77.399912 92.399823 92.399623 96601.8504 2.00  1011

628.318531 299.199299 299.199292 159.733100 159.733100 159.733100 159.733100 159.733100 159.733100 77.399911 77.399912 92.399869 92.399859

628.318531 299.199298 299.199299 159.733100 159.733100 159.733100 159.733100 159.733100 159.733100 77.399896 108.405810 92.399775 92.399888

628.318531 299.199300 299.199300 159.733100 159.733100 159.733100 159.733100 159.733100 159.733100 77.399905 77.399910 92.399781 92.399859

close results for the items B, A, W and SD. In C9 case all MABC variants achieve very close results and the success rate is SR = 100% (considering dF = 1000 $/year = 0.1141 $/h). The best solutions reached by applying the high performing variants – MABC/P/Log and MABC/D/Cat for EcD problem (w = 1,

56

D.C. Secui / Energy Conversion and Management 89 (2015) 43–62

case C8) and MABC/D/Cat for EmD problem (w = 0, case C9) – are shown in Table 16. Also, in Table 16, the best compromise solutions (BCS) for the MABC/D/Cat and MABC/P/Log variants are shown. The BCS were selected from the Pareto-front solutions using two ranking fuzzy-based mechanism methods described in Section 2.3. The Pareto-front solutions were calculated using relation (14), where the value of the weight factor w was increased from 0 to 1 with a step of 0.1. In order to have a Pareto-front as complete as possible in the promising zone found in the range w 2 (0.6, 0.7), two additional values were considered for the weight factor (w = 0.64 and w = 0.68). For the w values within the range (0.6–0.7) the factor l⁄ has the highest values in the range (0.0828, 0.0856), as can be noticed in Table 17. Thus, the 13 points resulted in the objectives space (Fuel cost, Emission) are shown in Table 17. From Table 17 the BCS is obtained for a w = 0.68 (l⁄w=0.68 = 0.0872) with a total fuel cost of 124491.1616 $/h and an emissions value of 256551.5725 ton/h (with MABC/P/Log algorithm). Also, it can be noticed in Table 17 that both methods applied in order to rank the solutions indicate the same BCS. 7.4.2. Comparing MABC variants to other methods In Table 18 the results obtained with the highest performance MABC variants (MABC/P/Log and MABC/D/Cat) are compared with those obtained by other optimizing algorithms/techniques, for

both Fuel cost, as also for the Emission quantity. For Fuel cost the comparison is done based on the B, A, W, SD items, and for Emission based on Best Emission item. It must be specified that the value of the TOLM item was calculated considering the solutions obtained with the algorithms mentioned in Table 18. The Best Emission achieved with the MABC variants (for example, 176682.26460 ton/h for MABC/P/Log and 176682.264590 ton/ h for MABC/D/Cat) have a better performance than MBFA [20], MODE, NSGA-II [87]. Making the comparison with the results obtained with other algorithms/methods (like HBMO and IHBMO [85], CIABC [86], DE [88], IABC [23] etc) raises a series of issues described below: (i) the HBMO and IHBMO [85] algorithms do not properly calculate the Best emission value. Performing a recalculation, the correct values for the emissions are: 281016.0597 ton/ h instead of 161305.061 ton/h (for HBMO [85]), 280133.0596 ton/h instead of 160807.648 ton/h (for IHBMO [85]). These values are bigger that the values obtained with MABC variants (item B). (ii) in case of the CIABC [86] algorithm the constraint regarding the maximum power operating limit for unit 30 is broken (P30 [86] = 99 MW > P30,max = 97 MW – see the solution from Table 16 for CIABC [86]). If MABC/P/Log algorithm is

Table 16 The solutions corresponding to the Best cost, Best emission and BCS when the MABC/D/Cat, MABC/P/Log and CIABC [86] algorithms are applied (40-units). Algorithm

Best cost

Output

MABC/D/Cat

MABC/P/Log

CIABC [86]

MABC/D/Cat

Best emission CIABC [86]

MABC/D/Cat

Best compromise solution MABC/P/Log

P1 (MW) P2 (MW) P3 (MW) P4 (MW) P5 (MW) P6 (MW) P7 (MW) P8 (MW) P9 (MW) P10 (MW) P11 (MW) P12 (MW) P13 (MW) P14 (MW) P15 (MW) P16 (MW) P17 (MW) P18 (MW) P19 (MW) P20 (MW) P21 (MW) P22 (MW) P23 (MW) P24 (MW) P25 (MW) P26 (MW) P27 (MW) P28 (MW) P29 (MW) P30 (MW) P31 (MW) P32 (MW) P33 (MW) P34 (MW) P35 (MW) P36 (MW) P37 (MW) P38 (MW) P39 (MW) P40 (MW) Cost ($/h) Emission (ton/h) TOL (MW)

110.8000081 110.7998478 97.3999212 179.7330987 87.7999294 139.9999999 259.5996696 284.5996495 284.5996507 130.0000000 94.0000000 94.0000000 214.7597894 394.2793664 394.2793704 394.2793703 489.2793697 489.2793694 511.2793697 511.2793693 523.2793702 523.2793702 523.2793702 523.2793702 523.2793702 523.2793702 10.0000000 10.0000000 10.0000000 87.8004645 189.9999996 190.0000000 189.9999985 164.7999830 200.0000000 194.3968155 110.0000000 109.9999983 109.9999999 511.2793691 121412.540947 359901.397571 9.01  107

110.7998263 110.8029523 97.3999110 179.7330900 87.7999044 139.9999990 259.5996492 284.5996495 284.5996498 130.0000000 94.0000000 94.0000000 214.7597902 394.2793479 394.2793149 394.2793506 489.2793695 489.2793695 511.2793691 511.2793693 523.2793703 523.2793702 523.2793703 523.2793702 523.2793703 523.2793702 10.0000000 10.0000000 10.0000000 87.8067375 190.0000000 189.9999959 189.9999995 164.7998247 194.3884987 199.9999982 109.9994763 109.9999695 109.9999948 511.2793702 121412.591816 359901.503351 7.18  107

110.7704 110.7705 97.3707 179.7038 87.7705 140 259.5699 284.5704 284.5706 130 94 94 214.7298 394.2497 394.2496 394.2497 489.2499 489.2499 511.2497 511.2708 523.2729 523.2737 523.2747 523.2758 523.2758 523.2697 9.9912 9.9915 9.9909 87.7723 189.9824 190 190 164.7731 194.3642 200 110 110 109.938 511.3501 121412.4 359842.0 0.607800

114.0000000 114.0000000 120.0000000 169.3680972 97.0000000 124.2571469 299.7111594 297.9148326 297.2603479 130.0000000 298.4102742 298.0262812 433.5562147 421.7283602 422.7800257 422.7802015 439.4119664 439.4030913 439.4132875 439.4134585 439.4466782 439.4468932 439.7724577 439.7716054 440.1117982 440.1112710 28.9933236 28.9936888 28.9938865 97.0000000 172.3318327 172.3316979 172.3318840 200.0000000 200.0000000 200.0000000 100.8384255 100.8384327 100.8386375 439.4127410 129995.274021 176682.264590 8.97  107

113.9998 114.0000 120.0000 171.0044 96.9995 127.0000 298.2855 298.5914 296.9671 134.0000 297.8076 296.1682 434.7314 421.0079 423.0154 420.5107 441.3810 436.7767 436.4718 439.4165 436.9569 441.8603 438.5185 438.3634 442.1713 441.3550 29.0661 27.5502 31.3780 99.0000 170.6259 173.7795 173.2124 200.0000 200.0000 200.0000 100.0000 100.0000 100.0000 438.0276 129996.6 176633.0 1.81  1012

110.7998216 110.7998237 97.39991059 174.5504493 87.79991926 105.3999112 259.5996498 284.5996451 284.5996474 130 318.1921565 243.5996499 394.2793698 394.2793686 394.2793694 394.2793694 399.5195798 399.5195802 506.1985408 506.1801105 514.1472578 514.1455791 514.5237739 514.5386386 433.5196115 433.5195807 10 10 10 97 159.7330999 159.7330999 159.7330994 200 200 200 89.11413552 89.11413645 89.1141381 506.1879753 124490.903193 256560.266940 9.80  107

110.799823 110.7998243 97.39991167 174.5486362 97 105.3999124 259.5996487 284.5996462 284.5996449 130 318.2129706 243.5996501 394.2793677 394.2793686 394.2793696 394.2793702 399.5195794 399.5195805 506.1716929 506.2206802 514.1105547 514.1472783 514.5664453 514.4868056 433.5195803 433.5196074 10 10 10 87.80420249 159.7330995 159.7331057 159.7330992 200 200 200 89.11413643 89.11413656 89.11413522 506.1951354 124491.161616 256551.572517 7.28  107

57

D.C. Secui / Energy Conversion and Management 89 (2015) 43–62

(iv) the DE [88] algorithm obtains for Best emission a value of 176,680 ton/h, but the exact value is 176683.2612 ton/ h > 176682.2646 ton/h (obtained by MABC/P/Log). Therefore, the MABC variants provide higher quality solutions comparing to DE [88].

Table 17 The fuel cost and the emission obtained by the MABC/P/Log algorithm, for different values of the factor w (40-units). w

0.00 0.10 0.20 0.30 0.40 0.50 0.60 0.64 0.68 0.70 0.80 0.90 1.00

MABC/P/Log Fuel cost ($/h)

Emission (ton/h)

l⁄ (method 1)

l1 (method 2)

129995.2691 129516.6708 129095.4231 128814.0623 128401.8957 127711.2165 126445.8876 125435.8396 124491.1616 124042.7577 123054.0402 122047.9086 121412.5918

176682.2646 203839.7691 204937.2152 205947.6104 208129.6635 213534.9781 226033.5918 238936.8995 256551.5725 270197.0390 306960.4929 331426.7288 359901.5033

0.0724 0.0657 0.0688 0.0708 0.0734 0.0771 0.0828 0.0862 0.0872 (max) 0.0856 0.0794 0.0783 0.0724

0.0000 0.0558 0.1048 0.1376 0.1856 0.2661 0.4136 0.5312 0.5641 (max) 0.4896 0.2889 0.1554 0.0000

Analyzing the results from Table 18, related to the ‘‘Fuel cost’’, the following was found:

The values in bold correspond to BCS.

compared to CIABC [86], enforcing the same limit for unit 30 (P30,max = 99 MW), then the best emission obtained with MABC/P/Log is 176504.613635 ton/h (value that is better than the one obtained with CIABC [86] of 176633.0 ton/h). (iii) IABC, IABC-LS, ABCD, ABCDP-LS [23] algorithms use a relatively high tolerance (|TOL| > 0.0130 MW – see the TOL value from Table 18). If, for example, MABC/P/Log algorithm is applied considering the obtained tolerance by applying the IABC algorithm (TOLIABC [23] = 0.0130 MW), then the Best emission is 176680.8888 ton/h (value which is better than that obtained by IABC [23], so MABC/P/Log algorithm is more efficient than the IABC algorithm [23]).

(i) the MABC/P/Log and MABC/D/Cat variants have better performance than all other algorithms mentioned in Table 18 regarding B, A, W, SD items (exception makes the IABC-LS algorithm [50] which has better values in regards to item W). It must be specified that the CIABC algorithm [86] achieves the best cost given that it breaks two types of constraints: minimum power operating limits for units 27, 28 and 29 (P27 = 9.9912 MW < P27,min = 10 MW, P28 = 9.9915 MW < P28,min = 10 MW and P29 = 9.9909 MW < P29,min = 10 MW – see Table 16) and real power balance constraint (TOLCIABC [86] = 0.6078 MW). If an MABC variant, for example, MABC/P/Log algorithm is applied considering the obtained tolerance by using the CIABC algorithm [86] (TOLCIABC [86] = 0.6078 MW) and for units 27, 28 and 29 the minimum limits from above are imposed, then the Best cost obtained with MABC/P/Log is 121405.7417 $/h (value that is the better than the one obtained with CIABC [86] of 121412.4 $/h); (ii) the MABC/P/Log and MABC/D/Cat algorithms have better performances than the well known and frequently used algorithms in solving different mathematical and engineering problems (like: MBFA [20], DE [88], EP [27], ABC [11,86], GA [10], PSO [10,27,89], ACO [10], CSO [89], BBO

Table 18 The values of the items (B, A, W, SD, TOL) for the Fuel cost, Best emission and TOL for emission.

*

Items

Fuel cost

Algorithm

B ($/h)

A ($/h)

W ($/h)

SD ($/h)

TOLM* (MW)

Emission Best emission

TOLM* (MW)

MODE [87] NSGA-II [87] HBMO [85] IHBMO [85] MBFA [20] DE [88] IABC [23] IABC-LS [23] ABCDP [23] ABCDP-LS [23] CIABC [86] EP [27] CPSO-SQP [31] FCASO-SQP [73] h-PSO [69] PSO [27] PSO [89] PSO [10] ACO [10] GA [10] CSO [89] BBO [90] ABC [11] HDE [25] IABC [50] IABC-LS [50] ABC [86] ABCLogistic [86] ABCTend [86] MABC/P/Log MABC/D/Cat

121836.9839 124963.5028 121416.03 121412.7533 121415.653 121,840 121414.80 121412.72 121412.82 121412.74 121412.4 122624.35 121458.54 121456.98 121420.9027 123930.45 122588.5093 121800.13 121532.41 121996.40 121461.6707 121479.5029 121441.03 121698.51 121412.75 121412.73 121515.1 121440.2 121418.51 121412.591816 121412.540947

– – 122019.65 121875.58 – – – – – – 121623.15 123382.00 122028.16 122026.21 121509.8423 124154.49 123544.8853 121899.57 121606.45 122919.77 121936.1926 121512.0576 121995.82 122304.30 – – 124827.34 123314.12 122831.22 121431.576266 121431.779282

– – – – – – – – – – – – – – 121852.4249 – 124733.6795 122000.80 121679.64 123807.97 122844.5391 121688.6634 122123.77 – 121503.58 121471.61 – – – 121493.188471 121503.755217

– – – – – – – – – – – – – – – – – 84.21 45.58 320.31 – – – – – – – – – 18.16 19.16

9.99  105 1.00  104 2.50  103 0 2.00  104 9.99  105 1.45  102 1.09  102 1.80  102 1.19  102 0.6078 – – – 3.63  1012 – 1.81  1012 – 0.16  102 – 2.99  104 0.277  102 1.82  1012 7.00  105 1.006  104 9.99  105 – – – 7.18  107 9.01  107

176683.2718 176691.9677 161305.061 160807.648 176682.269 176,680 176682.25 176682.25 176682.25 176682.25 176633.0 – – –

1.00  104 1.00  104 4.70  103 0 1.00  104 1.00  104 1.30  102 1.99  102 1.80  102 1.49  102 1.81  1012 – – –

– – – – – – –

– – – – – – –

–  – – – – 176682.264600 176682.264590

– – – – – – 9.42  107 8.97  107

TOLM was calculated for each method M with the relation (12), using the solutions indicated in the mentioned references. ‘‘–‘‘ data not available.

58

D.C. Secui / Energy Conversion and Management 89 (2015) 43–62

[90]) or different variants of the ABC algorithm recently published (IABC, IABC-LS, ABCDP, ABCDP-LS [23], CIABC [86], IABC, IABC-LC [50], ABCLogistic, ABCTend [86]); (iii) The tolerance (TOL) achieved for the solutions obtained with the MABC/P/Log and MABC/D/Cat algorithms are smaller than 106 MW. It must be mentioned that the tolerance TOLM which is relatively high – (0.001, 0.02) MW range – for some methods from Table 18, can be partially explained by the fact that the solutions presented in the corresponding references were rounded to two or three digits. The average computation time is 1.92 min. (for the case C8) and 12.03 s. (for the case C9).

generating very stable and high quality solutions – Table 6). However, from a mathematical point of view, between the MABC variants applied in the five cases (C1, C2, C3, C5 and C9) statistically significant differences can exist. Therefore, these cases will not be taken into consideration in the ranking process of the MABC variants, as not to modify their final standings. (iii) in contrast, in four cases (C4, C6, C7 and C8), there are differences between MABC variants for B, A, W, SD and SR items (these differences may determine values greater than the assumed limit of dF = 1000 $/year). Therefore, only these cases will be taken into consideration in the MABC variants ranking.

7.5. Comparing the results Analyzing the results in Table 6, obtained through the MABC variants for the nine studied cases (C1–C9), the following can be noticed: (i) considering the items A and SD the best performing MABC variants are MABC/D/Cat (for the cases C1, C2, C5, C7 and C9), MABC/D/Log (for the cases C3 and C6) and MABC/P/Log (for the case C4 and C8); (ii) there are five cases (C1, C2, C3, C5 and C9) in which the MABC variants lead to have very close values for B, A, W and SR items. The differences between W and B items, for cases C1, C2, C3, C5 and C9, are lower than the assumed cost of dF = 1000 $/year. Therefore, in such cases, the economic effect obtained through a certain MABC variant, compared to another MABC variant, is considered insignificant (this means that all MABC variants behave just as good,

7.5.1. Statistical test The performance comparison for MABC variants, selection schemes (P and D) and modalities (Cat, Log and Random) of generating random sequences (GRS) is performed using Mann–Whitney U statistical test with Bonferroni–Holm correction. Mann–Whitney U statistical test is a non-parametric test based on ranks, which shows if there are significant differences between two independent samples. Bonferroni–Holm correction is applied when performing repeated comparisons between the same collections of data and is based on the adjustment of p-value obtained by applying Mann–Whitney U test [84]. The null hypothesis (H0) for applying Mann–Whitney U test is: there is no significant difference between the results obtained by two methods Mx and My; the alternate hypothesis (H1) is: there is a significant difference between the results of Mx and My methods. If the adjusted p-value (using Bonferroni–Holm corrections) is less than or equal to the level a, then the null hypothesis is rejected

Table 19 Results of comparison between MABC variants (M1–M6), for the cases (C4, C6, C7, C8), using Mann–Whitney U test. Cases Algorithm

C4

C6

C7

C8

Indicator (RC4,C6,C7,C8)Mx

M1 vs.

M2 M3 M4 M5 M6

Idem Idem Idem Idem Idem

M2 M1 (0.039) Idem M5 Idem

Idem M1 (0.030) M1(0.000) Idem M1 (0.042)

Idem M1 (0.012) Idem Idem Idem

(RC4,C6,C7,C8)M1 = 5 (Rank 3)

M2 vs.

M1 M3 M4 M5 M6

Idem Idem Idem Idem M2 (0.029)

M2 (0.000) M2 (0.000) M2 (0.000) Idem M2 (0.000)

Idem Idem M2 (0.005) Idem Idem

Idem Idem Idem Idem Idem

(RC4,C6,C7,C8)M2 = 6 (Rank 1.5)

M3 vs.

M1 M2 M4 M5 M6

Idem Idem Idem Idem Idem

M1 M2 Idem M5 Idem

M1 Idem Idem Idem Idem

M1 Idem Idem Idem Idem

(RC4,C6,C7,C8)M3 = 0 (Rank 6)

M4 vs.

M1 M2 M3 M5 M6

Idem Idem Idem Idem M4 (0.001)

Idem M2 Idem M5 Idem

M1 M2 Idem M5 M6

Idem Idem Idem Idem Idem

(RC4,C6,C7,C8)M4 = 1 (Rank 4.5)

M5 vs.

M1 M2 M3 M4 M6

Idem Idem Idem Idem M5 (0.005)

M5(0.000) Idem M5 (0.000) M5 (0.000) M5 (0.000)

Idem Idem Idem M5 (0.008) Idem

Idem Idem Idem Idem Idem

(RC4,C6,C7,C8)M5 = 6 (Rank 1.5)

M6 vs.

M1 M2 M3 M4 M5

Idem M2 Idem M4 M5

Idem M2 Idem Idem M5

M1 Idem Idem M6 (0.037) Idem

Idem Idem Idem Idem Idem

(RC4,C6,C7,C8)M6 = 1 (Rank 4.5)

The p-value resulting from Mann–Whitney U test and Bonferroni–Holm correction is presented between round brackets (); Items marked ‘‘Idem’’ indicate that between two compared methods (Mx, My) there is no statistically significant difference; (RC4,C6,C7,C8)Mx represents the number of cases in which the method Mx is superior to other methods, considering the cases which produce more important economic effects (C4, C6, C7 and C8).

59

D.C. Secui / Energy Conversion and Management 89 (2015) 43–62 Table 20 Results of comparison between MABC variants (M1–M6) and ABC algorithm (M7), for each case (C1–C9), using Mann–Whitney U test. Cases Algorithm

C1, C2, C4 C9

C3

M1 M2 M3 M4 M5 M6

M1 M2 M3 M4 M5 M6

M1 M2 M3 M4 M5 M6

vs. vs. vs. vs. vs. vs.

M7 M7 M7 M7 M7 M7

(0.000) (0.000) (0.000) (0.000) (0.000) (0.000)

(RC1–C9)Mx (0.048) (0.048) (0.026) (0.036) (0.014) (0.045)

(RC1–C9)M1 = 9 (RC1–C9)M2 = 9 (RC1–C9)M3 = 9 (RC1–C9)M4 = 9 (RC1–C9)M5 = 9 (RC1–C9)M6 = 9

and the alternative hypothesis is accepted. The significance level for the null hypothesis testing is a = 5%. Mann–Whitney U test was performed using specialized software called SPSSÒ (Statistical Package for the Social Sciences). To test the performance of methods, selection schemes or GRS modalities the data resulting following 100 simulations (the Cost F) is grouped in three ways: (i) by cases; (ii) by both cases and selection schemes; (iii) by both cases and GRS modalities grouping. The results from applying Mann–Whitney U test and Bonferroni–Holm correction are shown in Tables 19–22. When comparing two methods ((Mx vs. My, x – y), in Tables 19–22 the better method was recorded. If between the compared methods there is no significant difference, the table entry is marked with ‘‘Idem’’. Also, alongside Mx or My method, the p-value resulting from Mann–Whitney U test with Bonferroni–Holm correction was specified in round brackets. 7.5.2. Comparison of the results obtained by MABC variants Table 19 shows the pairwise comparisons of the MABC variants (M1–M6), for the cases C4, C6, C7 and C8, using Mann–Whitney U test

and Bonferroni–Holm correction. For each method Mx the number of times when the method was superior to other methods was added up, leading to item (RC4,C6,C7,C8)Mx to be calculated. Then a hierarchy of all MABC variants was set. Analyzing the behavior of the MABC variants (M1–M6) for the analyzed cases (C4, C6, C7, C8), the following can be seen from Table 19: (i) based on RC4,C6,C7,C8 indicator, the MABC variants can be ranked as follows: M2 = M5 > M1 > M4 = M6 > M3; (ii) M2 (MABC/P/Log) and M5 (MABC/P/Log) have a better behavior than other variants, both on large systems (52-units system) and small ones (13 units system); (iii) M3 variant (MABC/D/Rand) has the worst behavior in the cases C4, C6, C7 and C8. 7.5.3. Comparison of the MABC variants and ABC algorithm In Table 20 can be noticed that all the MABC variants (M1–M6) are superior to ABC algorithm (M7) – both when each case was taken under consideration individually, and for all cases (RC1–C9)Mx = 9 analyzed as a whole. Also MABC variants behave better than ABC algorithm (M7) for both large systems (40-units and 52-units systems) and for smaller ones (6-units and 13-units systems). 7.5.4. Comparison of the selection schemes performance This paragraph refers to the results of the comparisons between the performance of the selection schemes (Disruptive(D)/Proportional(P)), with which the MABC variants were endowed. It is noted that when two selection schemes were compared, the GRS modality remains the same (Cat, Log or Random). The results of the comparisons (using Mann–Whitney U test) are shown in Table 21, and it can be seen that the Disruptive scheme (D) is better than the classical Proportional scheme (P), when map cat was used. In all the other situations there are no statistically significant differences between the compared variants.

Table 21 The results of comparing the selection schemes (Disruptive (D)/Proportional (P)), using Mann–Whitney U test, when the random sequence generating is maintained the same (Cat, Logistics or Random). Cases

C4

C6

C7

C8

Indicator RC4,C6,C7,C8

Cat Log Random

Idem Idem Idem

Idem Idem Idem

D (0.000) Idem Idem

Idem Idem Idem

(RC4,C6,C7,C8)D = 1 (Rank 1)

Cat Log Random

Idem Idem Idem

Idem Idem Idem

D Idem Idem

Idem Idem Idem

(RC4,C6,C7,C8)P = 0 (Rank 2)

Selection schemes (D/P)

Map

Disruptive (D) vs. Proportional (P) Disruptive (D) vs. Proportional (P) Disruptive (D) vs. Proportional (P) Proportional (P) vs. Disruptive (D) Proportional (P) vs. Disruptive (D) Proportional (P) vs. Disruptive (D)

Items marked ‘‘Idem’’ indicate that between the compared schemes there is no statistically significant difference.

Table 22 The results of comparing the GRS modalities (Cat/Log/Rand), using Mann–Whitney U test, when the selection scheme is maintained the same (Disruptive (D) or Proportional (P)). Cases

C4

C6

C7

C8

Indicator RC4,C6,C7,C8

Disruptive (D) Disruptive (D) Proportional (P) Proportional (P)

Idem Idem Idem Cat (0.001)

Log Cat (0.037) Log Idem

Idem Cat (0.015) Log Rand

Idem Cat (0.012) Idem Idem

(RC4,C6,C7,C8)Cat = 4 (Rank 2)

Disruptive (D) Disruptive (D) Proportional (P) Proportional (P)

Idem Idem Idem Log (0.005)

Log(0.000) Log(0.000) Log(0.000) Log(0.000)

Idem Idem Log (0.004) Idem

Idem Idem Idem Idem

(RC4,C6,C7,C8)Log = 6 (Rank 1)

Disruptive (D) Disruptive (D) Proportional (P) Proportional (P)

Idem Idem Cat Log

Cat Log Idem Log

Cat Idem Rand(0.037) Idem

Cat Idem Idem Idem

(RC4,C6,C7,C8)Rand = 1 (Rank 3)

GRS modality

Scheme (D/P)

Cat Cat Cat Cat

vs. vs. vs. vs.

Log Rand Log Rand

Log Log Log Log

vs. vs. vs. vs.

Cat Rand Cat Rand

Rand Rand Rand Rand

vs. vs. vs. vs.

Cat Log Cat Log

Items marked ‘‘Idem’’ indicate that between the compared GRS modalities there is no statistically significant difference.

60

D.C. Secui / Energy Conversion and Management 89 (2015) 43–62

0.5

Cost F ($/h)

0.4

C1: 6-units, 1263 MW

C2: 13-units, 1800 MW

C3: 13-units, 2520 MW

C4: 13-units, 1800 MW, with PL

C5: 13-units, 2520 MW, with PL

C6: 52units, 4x1800 MW

C7: 52-units, 4x2520 MW

C8:40-units, 10500 MW

0.1

0.3

0.05

0.2 0 0

200

400

600

800

1000

0.1 No. of iterations

0

0

5000

10000

15000

20000

25000

Fig. 4. The convergence process for MABC/D/Cat algorithm, cases C1–C8.

7.5.5. Comparison of the random sequences The comparison of the performance for GRS modalities (Cat, Log and Random), with which the MABC variants were endowed, are presented in Table 22. It needs to be mentioned that when two GRS modalities were compared, the selection scheme (D or P) remains the same. Based on Table 22 the following can be seen: (i) for all the cases C4, C6, C7 and C8, logistic map ((RC4,C6,C7,C8)Log = 6) is superior to the other two (cat map or Random). The hierarchy of the GRS modalities is Logistic > Cat > Random; (ii) logistic map has a good behavior in combination with Proportional scheme ((RC4,C6,C7,C8)Log = 4 for the scheme Proportional); (iii) logistic and cat maps have a close performance if Disruptive scheme is used; (iv) ‘‘Random’’ has poor performance for both D scheme and P scheme ((RC4,C6,C7,C8)Rand = 1). 7.6. Convergence process Fig. 4 shows the Cost F variation throughout the optimizing process by MABC/D/Cat algorithm for the eight cases (C1–C8). It is noticed that convergence is good in all the situations. For a single graphical representation of Cost F values, these are normalized within the [0, 1] range for each case. The value zero (0) corresponds to the best value of the Cost F (item B), and the value one (1) corresponds to Cost F maximum (item W). The number of iterations for the analyzed systems is different. Therefore, to get a clearer picture on the convergence, the optimization process for the four cases (C1: 6-units, C2 and C3: 13-units with 1800 MW and 2520 MW, C4: 13-units, 1800 MW with PL and C8: 40-units, 10,500 MW) is detailed in Fig. 4.

8. Conclusions This paper proposes a modified ABC algorithm (MABC), which is successfully applied for the economic dispatch problem solving, taking into account inequality and equality constraints, the valve-point effects, the emission pollutions and the transmission line loses. The MABC algorithm proposes a new relation for the solutions update, applied both in employed bee phase, and in onlooker bees phase. Four test systems having 6, 13, 40 and 52 units were analyzed. The results indicate that the proposed relation is effective, especially for 13 units, 40 units and 52 units systems. The behavior of MABC algorithm is also assessed when it is endowed with three modalities for generating random sequences (Cat, Log and Random) and two selection schemes of the solutions (disruptive selection and classical proportional selection). The combination between the GRS modalities and the selection schemes determine six variants of MABC, which are applied to nine cases (C1–C9) and then compared in pairs through Mann–Whitney U test. Mann–Whitney U test shows that the MABC

variants are superior to classical ABC algorithm for all the analyzed cases. The best performing MABC variants – MABC/D/Cat (for cases C1, C2, C5, C7 and C9), MABC/D/Log (for cases C3, C6) and MABC/P/Log (for the case C4 and C8) – behave better than the vast majority of optimization methods with which they were compared (methods relatively recently reported in the scientific literature and also well known methods (PSO, HS, ABC, DE) applied by the author). The best solutions obtained meet all of the equality and inequality constraints with high accuracy. The MABC variants have the ability to offer high-quality solutions for both small systems (6-units and 13 units) and large systems (40-units and 52-units). Also, the solutions obtained have very good stability for the 6-units and 13- units systems, and good stability for the 40-units and 52-units systems. Following the MABC algorithm structure it is noticed that it retains the positive characteristics of the original ABC algorithm: simple in concept, easy to implement, allows the hybridization with other classical or artificial intelligence based algorithms, flexible, can include the chaos element either in the solution update relations, or as a local search procedures, uses a moderate number of specific parameters. Compared to other metaheuristic algorithms (ABC, PSO, HS, DE, BBO, FA, GA etc) the MABC algorithm keeps the exploration–exploitation balance better, having the capacity to obtain stable and high quality solutions.

References [1] Nanda J, Hari L, Kothari ML. Economic emission load dispatch with line flow constraints using a classical technique. IEE Proc Gener Trans Distrib 1994;141(1):1–10. [2] Papageorgiou L, Fraga E. A mixed integer quadratic programming formulation for the economic dispatch of generators with prohibited operating zones. Electr Power Syst Res 2007;77(10):1292–6. [3] Chen SD, Chen JF. A direct Newton–Raphson economic emission dispatch. Int J Electr Power Energy Syst 2003;25(5):411–7. [4] Gaing Z-L. Particle swarm optimization to solving the economic dispatch considering the generator constraints. IEEE Trans Power Syst 2003;18(3):1187–95. [5] Noman N, Iba H. Differential evolution for economic load dispatch problems. Electr Power Syst Res 2008;78(3):1322–31. [6] Fesanghary M, Ardehali MM. A novel meta-heuristic optimization methodology for solving various types of economic dispatch problem. Energy 2009;34(6):757–66. [7] Niknam T, Mojarrad HD, Nayeripour M. A new fuzzy adaptive particle swarm optimization for non-smooth economic dispatch. Energy 2010;35(4):1764–78. [8] Liang Y-C, Juarez JRC. A normalization method for solving the combined economic and emission dispatch problem with meta-heuristic algorithms. Int J Electr Power Energy Syst 2014;54(1):163–86. [9] Chiang C-L. Genetic-based algorithm for power economic load dispatch. IET Gener Trans Distrib 2007;1(2):261–9. [10] Pothiya S, Ngamroo I, Kongprawechnon W. Ant colony optimisation for economic dispatch problem with non-smooth cost functions. Int J Electr Power Energy Syst 2010;32(5):478–87. [11] Hemamalini S, Simon SP. Artificial bee colony algorithm for economic load dispatch problem with non-smooth cost functions. Electr Power Compo Syst 2010;38(7):786–803. [12] Vijay R. Intelligent bacterial foraging optimization technique to economic load dispatch problem. Int J Soft Comput Eng 2012;2(2):55–9. [13] Swain RK, Sahu NC, Hota PK. Gravitational search algorithm for optimal economic dispatch. Proc Technol 2012;6:411–9. [14] Yang X-S, Hosseini SSS, Gandomi AH. Firefly algorithm for solving non-convex economic dispatch problems with valve loading effect. Appl Soft Comput 2012;12(3):1180–6. [15] Bijami E, Jadidoleslam M, Ebrahimi A, Askari J, Farsangi MM. Implementation of imperialist competitive algorithm to solve non-convex economic dispatch problem. J Chin Inst Eng 2014;37(2):232–42. [16] Jeddi B, Vahidinasab V. A modified harmony search method for environmental/ economic load dispatch of real-world power systems. Energy Convers Manage 2014;78(2):661–75. [17] Subbaraj P, Rengaraj R, Salivahanan S. Enhancement of self-adaptive real coded genetic algorithm using Taguchi method for economic dispatch problem. Appl Soft Comput 2011;11(1):83–92. [18] Chaturvedi KT, Pandit M, Srivastava L. Particle swarm optimization with time varying acceleration coefficients for non-convex economic power dispatch. Int J Electr Power Energy Syst 2009;31(6):249–57.

D.C. Secui / Energy Conversion and Management 89 (2015) 43–62 [19] Pandit N, Tripathi A, Tapaswi S, Pandit M. An improved bacterial foraging algorithm for combined static/dynamic environmental economic dispatch. Appl Soft Comput 2012;12(11):3500–13. [20] Hota PK, Barisal AK, Chakrabarti R. Economic emission load dispatch through fuzzy based bacterial foraging algorithm. Int J Electr Power Energy Syst 2010;32(7):794–803. [21] Pandi VR, Panigrahi BK, Mohapatra A, Mallick MK. Economic load dispatch solution by improved harmony search with wavelet mutation. Int J Comput Sci Eng 2011;6(1):122–31. [22] Hardiansyah. Solving economic dispatch problem with valve-point effect using a modified ABC algorithm. Int J Electr Comput Eng 2013;3(3):377–85. [23] Aydin D, Ozyon S, Yasar C, Liao T. Artificial bee colony algorithm with dynamic population size to combined economic and emission dispatch problem. Int J Electr Power Energy Syst 2014;54(1):144–53. [24] Sayedi E, Farsangi MM, Barati M, Lee KY. A modified shuffled frog leaping algorithm for nonconvex economic dispatch problem. Power Energy Soc General Meeting, 2012 IEEE, San Diego, CA, p. 1–8. [25] Wang S-K, Chiou J-P, Liu C-W. Non-smooth/non-convex economic dispatch by a novel hybrid differential evolution algorithm. IET Gener Trans Distrib 2007;1(5):793–803. [26] He Da-K, Wang FL, Mao ZZ. A hybrid genetic algorithm approach based on differential evolution for economic dispatch with valve-point effect. Int J Electr Power Energy Syst 2008;30(1):31–8. [27] Victoire TAA, Jeyakumar AE. Hybrid PSO–SQP for economic dispatch with valve-point effect. Electr Power Syst Res 2004;71(1):51–9. [28] Alsumait JS, Sykulski JK, Al-Othman AK. A hybrid GA–PS–SQP method to solve power system valve-point economic dispatch problems. Appl Energy 2010;87(5):1773–81. [29] Niknam T, Mojarrad HD, Meymand HZ. A novel hybrid particle swarm optimization for economic dispatch with valve-point loading effects. Energy Convers Manage 2011;52(4):1800–9. [30] Kumar R, Sharma D, Sadu A. A hybrid multi-agent based particle swarm optimization algorithm for economic power dispatch. Int J Electr Power Energy Syst 2011;33(1):115–23. [31] Cai J, Li Q, Li L, Peng H, Yang Y. A hybrid CPSO–SQP method for economic dispatch considering the valve-point effects. Energy Convers Manage 2012;53(1):175–81. [32] Reddy AS, Vaisakh K. Shuffled differential evolution for large scale economic dispatch. Electr Power Syst Res 2013;96(3):237–45. [33] Pandi RV, Panigrahi BK, Bansal RC, Das S, Mohapatra A. Economic load dispatch using hybrid swarm intelligence based harmony search algorithm. Electr Power Compo Syst 2011;39(8):751–67. [34] Kuo CC. A novel coding scheme for practical economic dispatch by modified particle swarm approach. IEEE Trans Power Syst 2008;23(4):1825–35. [35] Karaboga D, An idea based on honey bee swarm for numerical optimization. Erciyes University, Kayseri, Turkey, Technical Report-TR06; 2005. [36] Karaboga D, Akay B. A comparative study of artificial bee colony algorithm. Appl Math Comput 2009;214(1):108–32. [37] Karaboga D, Basturk B. A powerful and efficient algorithm for numerical function optimization: artificial bee colony (ABC) algorithm. J Global Optim 2007;39(3):459–71. [38] Karaboga D, Akay B. A modified Artificial Bee Colony (ABC) algorithm for constrained optimization problems. Appl Soft Comput 2011;11(3): 3021–31. [39] Luo J, Wang Q, Xiao X. A modified artificial bee colony algorithm based on converge-onlookers approach for global optimization. Appl Math Comput 2013;219(20):10253–62. [40] Banharnsakun A, Achalakul T, Sirinaovakul B. The best-so-far selection in Artificial Bee Colony algorithm. Appl Soft Comput 2011;11(2):2888–901. [41] Zhu G, Kwong S. Gbest-guided artificial bee colony algorithm for numerical function optimization. Appl Math Comput 2010;217(7):3166–73. [42] Gao WF, Liu SY, Huang LL. A global best artificial bee colony algorithm for global optimization. J Comput Appl Math 2012;236(11):2741–53. [43] Gao WF, Liu SY. A modified artificial bee colony algorithm. Comput Oper Res 2012;9(3):687–97. [44] Sonmez M. Artificial Bee Colony algorithm for optimization of truss structures. Appl Soft Comput 2011;11(2):2406–18. [45] Akay B, Dervis Karaboga D. Artificial bee colony algorithm for large-scale problems and engineering design optimization. J Intell Manuf 2012;23(4):1001–14. [46] Haluk G, Cengiz MT. Comparative performance analysis of artificial bee colony algorithm for automatic voltage regulator (AVR) system. J Franklin Inst 2011;348(8):1927–46. [47] Sßahin AS, Kılıç B, Kılıç U. Design and economic optimization of shell and tube heat exchangers using Artificial Bee Colony (ABC) algorithm. Energy Convers Manage 2011;52(11):3356–62. [48] Huang SJ, Liu XZ. Application of artificial bee colony-based optimization for fault section estimation in power systems. Int J Electr Power Energy Syst 2013;44(1):210–8. [49] Ayan K, Kiliç U. Artificial bee colony algorithm solution for optimal reactive power flow. Appl Soft Comput 2012;12(5):1477–82. [50] Aydın D, Ozyon S. Solution to non-convex economic dispatch problem with valve point effects by incremental artificial bee colony with local search. Appl Soft Comput 2013;13(5):2456–66. [51] Akay B, Karaboga D. A modified Artificial Bee Colony algorithm for realparameter optimization. Inform Sci 2012;192(1):120–42.

61

[52] Eslami M, Shareef H, Khajehzadeh M. Optimal design of damping controllers using a new hybrid artificial bee colony algorithm. Int J Electr Power Energy Syst 2013;52(11):42–54. [53] Li Y, Wang Y, Bin L. A hybrid artificial bee colony assisted differential evolution algorithm for optimal reactive power flow. Int J Electr Power Energy Syst 2013;52(11):25–33. [54] Kıran MS, Gündüz M. A recombination-based hybridization of particle swarm optimization and artificial bee colony algorithm for continuous optimization problems. Appl Soft Comput 2013;13(4):2188–203. [55] Liao X, Zhou J, Ouyang S, Zhang R, Zhang Y. An adaptive chaotic artificial bee colony algorithm for short-term hydrothermal generation scheduling. Int J Electr Power Energy Syst 2013;53(12):34–42. [56] Khorsandi A, Hosseinian SH, Ghazanfari A. Modified artificial bee colony algorithm based on fuzzy multi-objective technique for optimal power flow problem. Electr Power Syst Res 2013;95(2):206–13. [57] Gao F, Fei F-X, Xu Q, Deng Y-F, Qi Y-B, Balasingham I. A novel artificial bee colony algorithm with space contraction for unknown parameters identification and time-delays of chaotic systems. Appl Math Comput 2012;219:552–68. [58] Akbari R, Hedayatzadeh R, Ziarati K, Hassanizadeh B. A multi-objective artificial bee colony algorithm. Swarm Evol Comput 2012;2(2):39–52. [59] Gu W, Yin M, Wang C. Self adaptive artificial bee colony for global numerical optimization. IERI Proc 2012;1:59–65. [60] Sinha N, Chakrabarti R, Chattopadhyay PK. Evolutionary programming techniques for economic load dispatch. IEEE Trans Evol Comput 2003;7(1):83–94. [61] Cat Map, Wikipedia Web Site, [retrieved 25.02.14]. [62] Logistic Map, Wikipedia Web Site, [retrieved 25.02.14]. [63] Alzaqebah M, Abdullah S. Hybrid artificial bee colony search algorithm based on disruptive selection for examination timetabling problems, Lecture Notes in Computer Science (including subseries Lecture Notes in Artificial Intelligence and Lecture Notes in Bioinformatics) 6831 LNCS 2011; p. 31–45. [64] Kuo T, Hwang S-Y. Using disruptive selection to maintain diversity in genetic algorithms. Appl Intel 1997;7(3):257–67. [65] Park J-B, Jeong Y-W, Shin J-R, Lee KY. An improved particle swarm optimization for nonconvex economic dispatch problems. IEEE Trans Power Syst 2010;25(1):155–66. [66] Alatas B. Chaotic harmony search algorithms. Appl Math Comput 2010;216(9):2687–99. [67] Gaing Z-L. Closure to discussion of ‘‘particle swarm optimization to solving the economic dispatch considering the generator constraints’’. IEEE Trans Power Syst 2004;19(4):2122–3. [68] Özyön S, Aydin D. Incremental artificial bee colony with local search to economic dispatch problem with ramp rate limits and prohibited operating zones. Energy Convers Manage 2013;65(1):397–407. [69] Hosseinnezhad V, Babaei E. Economic load dispatch using h-PSO. Int J Electr Power Energy Syst 2013;49(7):160–9. [70] Wang L, Li LP. An effective differential harmony search algorithm for the solving non-convex economic load dispatch problems. Int J Electr Power Energy Syst 2013;44(1):832–43. [71] Lin W-M, Gow H-J, Tsai M-T. Combining of direct search and signal-to-noise ratio for economic dispatch optimization. Energy Convers Manage 2011;52(1):487–93. [72] Niknam T. A new fuzzy adaptive hybrid particle swarm optimization algorithm for non-linear, non-smooth and non-convex economic dispatch problem. Appl Energy 2009;87(1):327–39. [73] Cai J, Li Q, Li L, Peng H, Yang Y. A hybrid FCASO-SQP method for solving the economic dispatch problems with valve-point effects. Energy 2012;38(1):346–53. [74] Khamsawang S, Jiriwibhakorn S. DSPSO-TSA for economic dispatch problem with nonsmooth and noncontinuous cost functions. Energy Convers Manage 2010;51(2):365–75. [75] Tsai MT, Yen CW. The influence of carbon dioxide trading scheme on economic dispatch of generators. Appl Energy 2011;88(12):4811–6. [76] Xiong G, Shi D, Duan X. Multi-strategy ensemble biogeography-based optimization for economic dispatch problems. Appl Energy 2013;111(11):801–11. [77] Reddy AS, Vaisakh K. Shuffled differential evolution for economic dispatch with valve point loading effects. Int J Electr Power Energy Syst 2013;46(3):342–52. [78] Lohokare MR, Panigrahi BK, Pattnaik SS, Devi S, Mohapatra A. Neighborhood search-driven accelerated biogeography-based optimization for optimal load dispatch. IEEE Trans Syst Man Cybern 2012;42(5):641–52. [79] Coelho LS, Mariani VC. Correction to ‘‘combining of chaotic differential evolution and quadratic programming for economic dispatch optimization with valve-point effect’’. IEEE Trans Power Syst 2006;21(3):1465. [80] Wang Y, Li B, Yuan B. Hybrid of comprehensive learning particle swarm optimization and SQP algorithm for large scale economic load dispatch optimization of power system. Sci China-Inform Sci 2010;53(8):1566–73. [81] Wang Y, Li B, Weise T. Estimation of distribution and differential evolution cooperation for large scale economic load dispatch optimization of power systems. Inform Sci 2010;180(12):2405–20. [82] Abdullah MN, Bakar AHA, Rahim NA, Jamian JJ, Aman MM. Economic dispatch with valve point effect using iteration particle swarm optimization. In: 47th

62

[83]

[84] [85]

[86]

D.C. Secui / Energy Conversion and Management 89 (2015) 43–62 international universities power engineering conference (UPEC), London, UK, Sept. 2012. p. 1–6. Bhattacharjee K, Bhattacharya A, Dev SHn. Oppositional real coded chemical reaction optimization for different economic dispatch problems. Int J Electr Power Energy Syst 2014;55(2):378–91. Mann–Whitney U test, Wikipedia Web Site, [retrieved 25.02.14]. Ghasemi A. A fuzzified multi objective interactive honey bee mating optimization for environmental/economic power dispatch with valve point effect. Int J Electr Power Energy Syst 2013;49(7):308–21. Shayeghi H, Ghasemi A. A modified artificial bee colony based on chaos theory for solving non-convex emission/economic dispatch. Energy Convers Manage 2014;79(3):344–54.

[87] Sharma R, Samantaray P, Mohanty DP, Rout PK. Environmental economic load dispatch using multiobjective differential evolution algorithm. In: Proceeding of international conference on energy, automation, and signal (ICEAS); 2011. p. 1–7. 28–30 Dec. 2011, DOI: 10.1109/ICEAS.2011.6147132. [88] Basu M. Economic environmental dispatch using multi-objective differential evolution. Appl Soft Comput 2011;11(2):2845–53. [89] Selvakumar AI, Thanushkodi K. Optimization using civilized swarm: solution to economic dispatch with multiple minima. Electr Power Syst Res 2009;79(1):8–16. [90] Bhattacharya A, Chattopadhyay PK. Solving complex economic load dispatch problems using biogeography-based optimization. Expert Syst Appl 2010;37(5):3605–15.