A Modified Strength Pareto Evolutionary Algorithm 2 based Environmental / Economic Power Dispatch

A Strength Pareto Evolutionary Algorithm 2 (SPEA 2) approach for solving the multi-objective Environmental / Economic Power Dispatch (EEPD) problem is presented in this paper. In the past fuel cost consumption minimization was the aim (a single objective function) of economic power dispatch problem. Since the clean air act amendments have been applied to reduce SO2 and NOX emissions from power plants, the utilities change their strategies in order to reduce pollution and atmospheric emission as well, adding emission minimization as other objective function made economic power dispatch (EPD) a multi-objective problem having conflicting objectives. SPEA2 is the improved version of SPEA with better fitness assignment, density estimation, and modified archive truncation. In addition fuzzy set theory is employed to extract the best compromise solution. Several optimization run of the proposed method are carried out on 3-units system and 6-units standard IEEE 30-bus test system. The results demonstrate the capabilities of the proposed method to generate well-distributed Pareto-optimal non-dominated feasible solutions in single run. The comparison with other multi-objective methods demonstrates the superiority of the proposed method.


INTRODUCTION
The objective of Economic Power Dispatch EPD of electrical power system is to schedule the committed generating unit outputs so as to meet the load demand plus real power transmission loss at minimum operating cost while satisfying all units and system equality and inequality constraints.
The increasing public awareness of environmental protection and the passage of the U.S clean air act amendments of 1990 have forced utilities to modify their design or operational strategies to reduce pollution and atmospheric emission of thermal power plants, Abido, 2001.Strategies as switching to lower emission fuels or installation of pollutant clearing equipment requires considerable capital outlay, while considering emission as a constraint to be satisfied in EPD problem solution, or minimization of emission side by side with fuel cost requires modifying existing dispatching programs to include emission ,Talaq et al, 1994.During recent years the last idea received much attention due to development of a number of multiobjective techniques.
Linearly combined fuel cost and the amount of emission as a weighted sum convert the multiobjective EEPD problem to a single-objective optimization problem, Perez-Guerrero, 2005.By varying the weights a set of potential solutions Paretooptimal set were found, unfortunately this requires many runs as many as number of Pareto-front individuals to form the Pareto-optimal solutions as well as the diversity of Pareto-optimal set along Pareto-front depend on the diversity of the weights.Alternatively, many attempts to solve the EEPD were done by handling emission minimization as another objective function using stochastic multi-objective approaches.M. Abido used NSGA , Abido, 2001, NPGA ,Abido, 2003aand SPEA Abido, 2003 to solve EEPD problem and obtained the better results using SPEA, although T.F.Robert, King et al, 2004, and S. Agrawal ,Agrawal et al, 2008 obtains a better results less fuel cost than SPEA using NSGA-II, FCPSO respectively but the corresponding emission amount increased due to conflicting objectives.An extensive evaluation may be done by inserting a real penalty price factor the tax legislated to calculate the overall cost and helping to state the better solution.SPEA2 which is used in this paper is the improved version of SPEA to solve the multi-objective EEPD problem, SPEA2 uses an archive to store the nondominated individuals, fitness assignment which takes into account both dominated and dominating other individuals and archive truncation to maintain a constant archive size with good diversity ,Zitzler et al, 2001.A fuzzy-based mechanism is used to extract a Pareto-optimal solution as the best compromise solution.Two test systems were used to compare and state the superiority SPEA2 with other multi-objective approaches.

PROBLEM FORMULATION
The EEPD, involves the simultaneous optimization of fuel cost and emission amount as multi-objective conflicting problem, is generally formulated as follows.

Problem's Objectives
1. Minimization of fuel cost: The generator consumption fuel cost curves are represented by quadratic functions where the total consumed fuel cost F(PG) in ($/h) can be expressed as Where n is the number of generators, , , and are cost coefficients of the i th generator, and is the real power output of the i th generator , Abido, 2001.2. Minimization of emission amount: The total emission E(PG) in (ton/h) of atmospheric pollutants such as sulfur oxides SO x or nitrogen oxides NO x caused by fossil-fueled thermal generating units may be expressed as ,Perez-Guerrero, 2005.
Where , , , and are i th generator emission coefficients and P G is the vector of real power outputs of system generators so,

Problem's Constraints
1. Generator capacity constraints: The real power output of each generator is restricted by the lower and upper limits as follows:

The EEPD Problem Formulations
The EEPD problem is formulated as: Minimize [ F(P G ) , E(P G ) ] Subject to: The real power loss in transmission lines PL can be considered as another objective function which needs to be minimized, Abido, 2001.

MULTI-OBJECTIVE OPTIMIZATION 3.1 Basic Concepts
Multi-objective optimization (MOP) to several objective functions (which are often competing and conflicting objectives) simultaneously obtains a set of optimal solutions (instead of one solution) since none of them are considered better with respect to all objectives.The MOP to N obj objective function can be formulated as ,Abraham et al, 2005: Subject to the m equality constraints: And the p inequality constraints: x n ] T is the vector of decision variables.The Pareto-optimality is explained as follows: a vector x* is Pareto-optimal if for every other vector x for all i=1,2,…,N obj and f j (x*) < f j (x) for at least one j .
Where is the feasible set (its elements satisfy Eqs.(7.a) and (7.b), the vector x* is called non-dominated since there is no such x which dominate it, all nondominated solutions (vectors) forms Pareto-optimal set.

The Strength Pareto Evolutionary
Algorithm 2 (SPEA2) E. Zitzler, M. Laumanns, and L. Thiele, had developed SPEA in 1999 , Zitzler and Thiele, 1999, yet in 2001they published SPEA2 ,Zitzler et al, 2001 as enhancement version, by fixing (improving) the potential weakness in fitness assignment, density estimation and archive truncation.Like the earlier SPEA, SPEA2 has external archive to store nondominated individual and only these individuals form the mating pool.The SPEA2 has an overall algorithm as follow: Step1: Initialization: Generate an initial population P 0 and create the empty archive (external set) ̅ 0 =Ø; Set t = 0.
Step 2: Fitness assignment: Calculate fitness values of individuals in P t and ̅ t , each individual i in the archive ̅ t and the population P t is assigned a strength value S(i), representing the number of solutions (individuals) which it dominates: Where | | denotes the cardinality of a set (the number of elements in a set), + stands for multi-set union, corresponds to Pareto dominance relation ( refers to that individual i which dominates individual j) and i,j ̅ t + P t .On the basis of the S values, the raw fitness is determined by the strengths of its dominators in both archive and population, The raw fitness R(i) of an individual i is calculated: Since the non-dominated individuals would have the same raw fitness value R(i)=0, while a high R(i) value means that individual i is dominated by many individuals.Therefore additional density information is incorporated to discriminate between individuals having identical raw fitness values.The density estimation technique used in SPEA2 is an adaptation of the k th nearest neighbor method ,Zitzler et al, 2001, for each individual i the distances (in objective space) to all other individuals j in archive and population are calculated and stored in a list.After sorting the list in ascending order, the k th element gives the distance sought.As a common setting, k equal to the square root of the entire population size √ ̅ .Where N is a population (P t ) size and ̅ is an archive ( ̅ ) size.
The distance between the individuals i an k is denoted as .Density D(i) corresponding to i is defined by By adding D(i) to the raw fitness value R(i) of an individual i yields its fitness F(i): Step

Real-Valued (Coded) Genetic Algorithm
The calculations of objective functions many times in the process of the simulation computer program make coding and decoding individuals from and to binary-coded time consuming as well as coding one individual of 3-unit system having the same accuracy obtained in real-valued string needs 3×20 bits, since there are 3 variables of P G and each variable needs at least 4 digits of 0, 1, assuming the population size is 5. Therefore the genetic algorithm string is represented in a vector of real-valued of power outputs of system generators as: 1. Recombination (crossover): A blending crossover operator has been employed.This operator recombines the i th parameter (gene) values of individuals x, y (selected for recombination), the offspring appear as follows: Where the offspring are: And β is a randomly generated number between 0 and 1. 2. Mutation: a non-uniform mutation operator is employed in this study ,Michalewicz, 1996, the new value x i " of the parameter x i after mutation at generation t is given as: Where rand is a random number generator between (0,1), t is iteration index, T is maximum number of iterations and b=5 ,Michalewicz, 1996.

Best Compromise Solution
The solutions obtained for best (minimal) fuel cost and for best (minimal) emission amount were giving an image about the optimized objective functions, but in order to adopt one solution as the best compromise solution to the decision maker's judgment, the proposed approach presents a fuzzy-based mechanism to extract a Pareto-optimal solution as the best compromise solution , Abido, 2003.Due to the imprecise nature of the decision maker's judgment, the i th objective function of a solution in the Pareto-optimal set Fi is represented by a membership function which is a Z-function asymmetrical polynomial curve, defining by: Where and are the maximum and minimum values of the i th objective function, respectively.For each non-dominated solution h, the normalized membership function µ h is calculated as: Whereas M is the number of non-dominated solutions.The best compromise solution is the one having the maximum value of µ h .By arranging all solutions in Pareto-optimal set in descending order according to their membership function will provide the decision maker with a priority list of nondominated solutions.This will guide the decision maker in view of the current operating conditions.The implementation flow chart of the proposed approach is shown in Fig. 1 4

. IMPLEMENTATION OF THE PROPOSED METHOD
To satisfy the problem constraints, the following steps had been made: (a) The initial population is generated within a capacity limit of each generator as well as the recombination and mutated elements are as follows: Whereas rand is a random number generator between (0, 1).

(b)
The power balance constraint is satisfied as follows, the traditional B-matrix loss formula is used to calculate the real power transmission loss ∑ ∑ ∑ (17) By choosing the r th unit randomly, it's assumed that the r th reference unit power output is responsible of bucking up the remaining load after the other (n-1) units have been assigned theirs output power, ∑ Rewriting eq. ( 18) in order to form a polynomial with P Gr is the variable as follows: (19) Where a ,b, and c represent parameters depending on the B-matrix coefficients of the power loss equation of the test system used and on the powers of the (n-1) generators.
Substituting PL from Eq. ( 19) in Eq. ( 18), yields The roots of the eq.( 20) represent the value of P Gr satisfying equality constraint, if neither roots located within unit power capacity limit, other generator is randomly chosen as r th reference unit also, if all units were filed to backup the remaining load another individual is generated randomly to replace this one.After recombination and mutation process each individual is checked for equality constraint violation, it is worth to mention that this process works as another mutation operator.Since, it may vary the value of one gene (generating unit) to achieve equality constraint satisfaction, or it may be the reason of destroying the recombination operator process (if the recombination gene is chosen to be the first or the last in the individual string).(c) Archive truncation modification By formulating (a distance matrix) with size [ ̅ t+1 × ̅ t+1 ], this matrix contains each possible pair of two individuals in archive set ̅ t+1 (before truncation process), reducing the dependency on the k th nearest neighbor method, the process is as follows: Step1: form the distance matrix, calculating the distance between each pair of non-dominated individuals (i,j) in multi-objective space.
Step2: searching the matrix for smallest element (represents the minimum distance between any two individuals i,j) .
Step3: If individual i would be eliminated from archive and its row and column would be eliminated too, else the individual j is chosen for elimination and the size of ̅ t+1 and distance matrix reduced by one.

RESULTS AND DISCUSSION
In this research two systems were adopted in order to investigate the effectiveness and applicability of the proposed method, several simulations runs were done for each test and an identical population and archive sizes were used ,Zitzler et al, 2001, the parameters used for all cases are as follows: 200 individuals for each population size and archive size for case (1) and 100 individuals for each population size and archive size for case (2) and the simulations were run for 1000 generations, crossover and mutation probabilities were 1(100%) and 0.01 respectively.

Case (1): Three Generating Units System
The three generators test system whose data are given in Tables 1, 2 and 3 ,King, 2003, the system demand is 850 MW, and the system transmission losses are calculated using a simplified loss expression ,King, 2003: +0.00009 +0.00012 ( 21) The coefficients of eq. ( 19) are a=B rr , b=0 and c =∑ ∑ It is important to maintain that the cross-point in the proposed approach applied to solve this test system is fixed in the second gene (second generating unit position) to grantee new individuals generated next population.

Test (1): Fuel Cost and SO 2 Emission Objective Functions:
In this test fuel cost with SO 2 emission were taken as objective functions to be minimized, Tables 4, 5, and 6 show the simulation results obtained in one run for best (minimum) fuel cost, minimum emission and best compromise solution in the Pareto front respectively, as compared toTABU search and NSGA-II, while the Pareto-front were plotted in Fig. 2. Table 4 shows a reduction in the consumption fuel cost by more than 100 $ per year than the results obtained by NSGA-II approach.While, Table 5 shows a reduction in SO 2 harmful emission by 5.6 ton per year than the results obtained by NSGA-II approach for this small system.
The minimum fuel cost and minimum emission solutions were drawn against generations (iterations) in Fig. 3.The average simulation run time for this test is 550 second.The convergence of non-dominated solutions to the true Pareto-optimal front region was done in early stages of the search process (20-30 % of the maximum generations limit), as shown in Fig. 3, the later stages is for convergence to the exact solutions (fine tuning).

Test (2): Fuel Cost and NOx Emission Objective Functions:
In this test fuel cost with NOx were taken as objective functions to be minimized, Tables 7, 8, and 9 shows the simulation results obtained in one run as compared to TABU search and NSGA-II, while the Pareto-front was plotted in Fig. 4. The average simulation run time for this test is 550 second.

Test (3): Fuel Cost, SO2 Emission and NOx Emission objective functions:
In this test fuel cost with SO 2 emission and NOx emission were taken as objective functions to be minimized Tables 10, and 11 shows the simulation results obtained in one run as compared to NSGA-II, while the Pareto-front was plotted in Fig. 5.The average simulation run time for the test is 780 second.
Although NSGA-II is a new powerful multi-objective technique, the results in the previous tables show that SPEA 2 gave not only better results but also with less population size and less maximum generations number, NSGA-II has population size 500 individuals and was run for 20000 generation ,King, 2003.

Case (2): Six Generating Units System
The standard 30-bus IEEE test system with 6generating units, King, 2004 is used with load demand 2.834 p.u. (100 MW power base), since this system has been already solved and validated by several multi-objective optimization techniques the comparison of SPEA2 with such techniques can show the potential of the proposed method, system data listed in Tables 12, and 13.Test (1): Fuel Cost and Emission Objective Functions (Real Power Transmission Loss is neglected): In this test, fuel cost with harm emission were taken as objective functions to be minimized, the system is considered as lossless and the equality constraint is as follows: ∑ Tables 14, 15, and 16 show the simulation results obtained in one run as compared to other approaches, while the Pareto-Optimal front was plotted in Figure 6.The average simulation run time for the test is 70 second.The minimum fuel cost and minimum emission solutions were drawn against generations (iterations) in Fig. 7.
The results obtained in Tables 14 and 15 are clearly demonstrated the superior of SPEA 2 over other multiobjective GA methods and also over the new multiobjective-PSO approach (FCPSO) with reduction in consumption fuel cost more than 190 $ per year than FCPSO approach (table IXV), also the fuel cost corresponding to minimum emission in table XV is less than the fuel cost corresponding to minimum emission FCPSO approach.

Test (2): Fuel Cost and Emission objective functions (transmission loss is included):
The exact value of the system losses can only be determined by means of a power flow solution, but in this research the B-coefficient matrix from, Perez-Guerrero, 2005 is used: B io =[(-0.01070.0060 -0.0017 0.0009 0.0002 0.0030)] the coefficient of eq. ( 19) is

∑ ∑
and c =∑ ∑ ∑ Tables 17, 18, and 19 show the simulation results obtained in one run as compared to other approaches, while the Pareto-front was plotted in Fig. 8.The average simulation run time for the test is 70 second.Other researches use load flow solution approach to determine the transmission loss instead of Bcoefficient matrix approach, the transmission loss obtained by load flow solution depends on the accuracy of the load flow solution such as in FCPSO approach for example, while the accuracy of the transmission loss obtained by the B-coefficient matrix approach is within 0.25 % of the load demand, which may be accepted when is taking into consideration the computational time which is spent in load flow solution approach compared to those which is using Bcoefficient matrix approach.

Test (3) Fuel Cost, Emission and Transmission Loss
Objective Functions: In this test real power transmission loss, fuel cost and harm emission were taken as objective functions to be minimized, Table 20 shows the three objective functions optimization results as compared with reference ,Wang, 2008, results.While, Table 21 shows the best compromise solution.The Pareto-front was plotted in Fig. 9.The average simulation run time for the test is 90 second.
It is important to mention that the research in ,Wang, 2008, neglects the linear coefficient Bio and the constant coefficient Boo of the B-matrix used to calculate the power loss PL then, in order to obtain a good comparison results these coefficients are neglected in the obtained results in this test only.All the simulation results obtained in this research were implemented on personal computer Pentium 4, 3.59 GHz with 1GB RAM using MATLAB version 7 programming language.

CONCLUSIONS
This paper presents multi-objective environmental/ economic power dispatch (EEPD) solutions using the proposed Strength Pareto Evolutionary Algorithm 2 (SPEA2).The proposed method has a diversitypreserving mechanism to find widely different Paretooptimal solutions.A distance matrix technique is implemented to provide the decision maker with diverse and manageable Pareto-front without destroying the shape and the boundary of the trade-off front.Moreover, a Fuzzy-based mechanism is employed to extract the best compromise solution over the trade-off curve as in Fig. 1.The Fuzzy system contains a fuzzy inference system, fuzzy controller with a rule base, and a defuzzifier.A triangular membership function is used.The results show that the proposed method is efficient for solving multiobjectives optimization whereas multiple Paretooptimal solutions can be found in one simulation run.The simulation results of all the tests which were done reveal that SPEA2 presents low computational time (less population size and less generations) and is suitable for on-line EEPD solutions.The proposed method has reliable convergence, high accuracy of solution, and better results than other multi-objective Figure2.Pareto -optimal front for case (1) test (1).
Figure3.Convergence of min.fuel cost and min.emission of the proposed method for system (1), test (1).
Figure7.Convergence of min.fuel cost and min.emission of the proposed method for system 2, test (1).

Hassan Abdullah Kubba Saif Sabah Sami A Modified Strength Pareto Evolutionary Algorithm 2 based Environmental /Economic Power Dispatch 108
Copy all nondominated individuals in P t and ̅ t to ̅ t+1 .If size of P t+1 exceeds ̅ then reduce ̅ t+1 by means of the modified truncation operator, so, at each iteration the individual i which has the minimum distance to another individual j and < is chosen for removal; otherwise if size of ̅ t+1 is less than ̅ then fill ̅ t+1 with dominated individuals in P t and ̅ t , by sorting the multi-set P t + ̅ t according to the fitness values and copy the first N − | | individuals i with F(i) ≥ 1 from the resulting ordered list to ̅ t+1 .

Table 4 .
Best fuel cost.

Table 7 .
Best fuel cost.

Table 10 .
Comparison of results for three objective functions.