Performance Assessment of a Set of Multi-Objective Optimization Algorithms for Solution of Economic Emission Dispatch Problem

,


Introduction
The Economic Load Dispatch (ELD) problem deals with the estimation of the scheduled real power generation from the committed units for best economic operation. Over the years the problem has become more complex due to the increasing effects of emissions from fossil fuel based power plants on the environment. The emission and fuel cost of each unit depend on the quantity of power to be generated. Both of them are nonlinear functions of power output. Minimum operating cost does not ensure minimum emission. Each operating condition must satisfy the power balance criterion and should obey the generating limits of the committed units. These can be considered as constraints. Generally, better quality fuel ensures less emission but it can be further reduced by proper scheduling of generation from different units. The cost coefficients and emission coefficients of these generating units do not match. Hence, achieving these two objectives, i.e. less cost and less emission is contradictory in nature. Thus, the EED problem has evolved as a modification of the ELD problem. Therefore, the EED problem is a multi-objective optimization problem with nonlinear constraints.
In [1][2], the power engineers solved the ELD problem by scheduling of the generation of multi-unit systems using the derivative based Gauss-Siedel and Newton-Raphson algorithms along with the Lagrangian multiplier. These conventional methods suffer from the problem of getting trapped in local minima and also fail for system discontinuities due to prohibited zones. These techniques are inadequate to solve multi-objective problems with nonlinear constraints. Chang et al. [3] rehabilitated the inherently multi-objective EED problem to a single objective one by assigning weights to the operating cost and emission. This weighted sum approach requires many runs of the same algorithm to find the Pareto optimal front. The solutions arrived at by this method do not ensure a uniform Pareto front. The trade-off information is lost when the function is concave. To avoid this bottom-hole different evolutionary based heuristic approaches have been introduced by many researchers [4][5]. These evolutionary algorithms have considered the two objectives simultaneously and are shown to perform better as compared to the conventional ones. Chiang et al. [6] made a further refinement and proposed an improved genetic algorithm to speed up the search process. He used the ϵ-constraint technique for efficient constraint handling and proposed a multiplier updating mechanism for better exploration of the search space.
Deb et al. [7] proposed the non-dominated sorting genetic algorithm which utilized rank and crowding distance as parameters to arrive at a compromise between the two conflicting objectives. This was applied to the multi-objective environmental economic load dispatch problem in [8]. The Pareto optimal front could be obtained by a single run of the algorithm. But, this population based genetic algorithm depends upon biologically inspired factors like mutation and crossover parameters. It needs further improvement in terms of exploring a wider area in the search space. Brar et al. [9] made improvements in the search space by adding the fuzzy inference system. Muthuswamy et al. [10] modified the non-dominated sorting technique by incorporating a dynamic crowding distance to improve the diversity of solutions in the search space. These algorithms fail when there are discontinuities in the cost function.
Nayak et al. [11] implemented another evolutionary algorithm, the artificial bee colony (ABC) optimization, and improved the convergence rate and reliability under the presence of the prohibited zones and ramp rate limits. Liang et al. [12] modified the ABC algorithm to form an improved artificial bee colony (IABC) by addition of a new skill called chaos ques in the search process. Mori et al. [13] made an excellent improvement in the exploration of search space through the implementation of the particle swarm optimization (PSO) for this multimodal problem. They also used adaptive parameter adjustment to improve the results. A significant improvement in search space exploration was made by Hadji et al. [14]. They incorporated a time varying acceleration of the particles to improve the robustness of the algorithm. Recently, a differential evolution (DE) algorithm came up which generates the next set of population of new particles by the addition of a differential vector obtained from the difference of the position vectors of two different particles other than the particle undergoing evolution [15]. This algorithm is still dependent on the bio-inspired parameters but is able to avoid premature convergence. Meza et al. [16] improved the algorithm by incorporating spherical pruning for better exploitation of the search space. Di et al. [17] introduced a marginal analysis correction operator to improve the constraint handling.
In [18], the particle swarm optimization algorithm has been developed which is based on the intelligence of flock of birds. The same has been improved and tested for multi-objective problems in [19][20][21]. The EED problem has been solved to decide the unit commitment of the power system by considering operational power flow and environmental constraints in [22]. But, it again utilized the method of conversion of the multi-objective problem to a single objective one. A new approach to optimization is proposed in [23] which hybridized adaptive PSO and DE for improvement of the search space. An improvement over ABC called as multiobjective global best artificial bee colony (MOGABC) optimization is suggested in [24] for better constraint handling in EED problem. The EED problem has been further modified and applied to the micro-grid containing renewable sources along with the conventional thermal power stations in [25]. It also converts the problem to a single objective one by incorporating a h-index.
In this paper, a new constraint handling mechanism has been implemented, and a new multi-objective optimization (MOP) algorithm, namely the multiobjective differential evolution with recursive distributed constraint handling (MODE-RDC) has been proposed. The constraint handling mechanism is suitably incorporated in three multi-objective optimization (MOP) algorithms, and the effectiveness of the algorithms has been tested under various load conditions.

Multi-objective optimization: a review
The main aim of the multi-objective optimization technique is to optimize two or more conflicting objectives simultaneously. The MOP is denoted by a decision variable vector, each element of which represents the objective functions [21]. The solution to the MOP is the optimum value of the vector function by considering all the constraints. A multi-objective minimization problem can be generalized as follows: h k (x ⃗ ) = 0; k =1, 2, …, K Solutions that are not dominated by other solutions within the given solution space are said to be non-dominated solutions. The front obtained by mapping such points onto the objective space is said to be the Pareto optimal front (POF) where, p is the set of non-dominated particles.

Economic emission dispatch problem
The generation schedule for minimum operating cost is called ELD. This schedule is obtained when the committed units of the power system are able to supply the load demanded and the associated transmission losses by satisfying the generator limits. The thermal generating units are associated with emissions which are highly polluting in nature. Therefore, it is essential to minimize the emissions along with the fuel cost. The problem has been transformed to EED problem. It is considered as a multi-objective optimization problem as minimum operating cost does not ensure minimum emission [5]. The operating cost of thermal power systems depends mostly on the cost of fuel used. The quantity of fuel used by each unit depends on the generated power, efficiency of turbine etc. The fuel cost characteristics of the generating units are normally of the second order polynomial of the generated power. Thus, the operating fuel cost of the ith generator supplying a real power P G i is given by where, ai, bi and ci are the coefficients of cost function.
The emission from the generator can be approximated as E i = α i + β i P Gi + γ i P Gi 2 (8) where, αi, βi, γi are the coefficients of emission function.
The aim of the ELD problem is to determine generation schedule for the minimum total fuel cost given by (9) subject to the constraints.
The total real power generation must be equal to the demand plus transmission losses ∑ P Gi N i=1 = P D + P L (10) where, PD is the load demand on the system and PL is the transmission loss. It is given by Kron's formula P L = ∑ ∑ P Gi B ij P Gj The constants Bij, Bi0 and B00 are dependent on the line parameters. The generated real power of each unit must be within the feasible lower and upper bounds. P Gi(min) ≤ P Gi ≤ P Gi(max) ; i =1, 2, …, N.
Thus, the sole objective of the EED problem is to optimize both the fuel cost and emission simultaneously. Hence, it is inherently a multi-objective optimization problem where these two objectives which are conflicting in nature need to be optimized. The set of all the potential compromised solutions is represented by the Pareto optimal front.
The problem can be stated as: minimize Subject to: g(P G )=0; h(P G )≤0 where, the equality constraint is represented by equation (14) and inequality constraint by equation (15).

Multi-objective optimization algorithms
The multi-objective evolutionary algorithms can be categorized mostly into four types in accordance with the algorithmic framework, such as indicator based, convolution based, memetic based and non-dominated sorting based. In this paper, we have solved the economic emission dispatch problem using three different multiobjective optimization algorithms that are based on nondominated sorting.

Non-dominated Sorting Genetic Algorithm-II (NSGA-II)
This algorithm was formulated by Deb et al. [7], [8] for solving the multi-objective optimization problems. It is initialized with a random population, and used some operators for covering the objective space uniformly on the Pareto set. For multi-criteria optimization it uses three strategies: non-domination sorting, ranking based on density and crowding comparison. The individuals are classified into several layers based on their rank and crowding distance. The diversity in the solutions is maintained by rejecting the solutions with lower crowding distance. The quality of solutions is ensured by selecting the individuals with lower ranks. The advantage of this algorithm is that the complexity of computation is lowered and elitism is maintained.

Non-dominated Sorting Multi-Objective Particle Swarm Optimization (NS-MOPSO)
Kennedy and Eberhart proposed that any optimization problem can be solved by mimicking the movement of a flock of birds and school of fish [18]. The social behavior of the swarm is to change their position and velocity to maximize their chance of getting food and follow the best successful neighbor. This lead to the formulation of particle swarm optimization (PSO). In this method of optimization, a local best and a global best solution are identified. The i th particle in the population having the best position (pbest) may be represented by p i , that gives the best fitness value represented as p i = (p i1 , p i2 , …, p iN ) (16) The old and new velocity of the particles will be shown in equation (17) and (18) respectively.
and the new position of the particle will be where, d = 1, 2, … , D is the dimension of the decision variables and i = 1,2, … , N, andχ is the constriction factor which constricts and controls the velocity magnitude. w, c 1 and c 2 are weight parameters and r 1 , r 2 are random numbers known as acceleration constants in the range [0, 1]. This method of improvement of position and velocity is applied to the non-dominating vectors to solve the multi-objective problem [14], [19] & [20].

Multi-Objective Differential Evolution with Recursive Distributed Constraint handling (MODE-RDC)
The differential evolution (DE) algorithm as developed by Stern and Price [15] is less dependent on bio-inspired mechanisms, and serves better for multivariable problems. The multi-objective optimization using differential evolution (MODE) proposed by Meza et al. [16] is an improvement of DE to solve multiple number of conflicting objectives simultaneously. The MODE is an evolutionary multi-objective optimization algorithm that retains the diversity of solutions on the Pareto front. This real coded stochastic algorithm uses an initial population to explore the search space by avoiding convergence to local optimal points. It uses two main operators: mutation and cross over. Each initial particle of the population is improved using these two operators.
The mutation operator uses a differential vector selected from the particles other than the target particle. Three vectors x r0,g , x r1,g and x r2,g are selected randomly from the population and the first one is updated with the difference of the other two. It is done as follows: where, v i,g is the mutant vector created from the target vector x r0,g and x r1,g , x r2,g are two other vectors; F∈(0,1+) is the scale factor that controls the rate of evolution.
A new trial vector (child) is created from the mutant vector and the target vector after cross over.
where, the cross over probability C r ∈[0,1]. The child and the parent are tested for their fitness. The one with the best fitness is selected for participation in the next generation.
This is done for all i =1, 2, ⋯, n where n = population size. The steps of the proposed MODE-RDC algorithm are summarized as follows: I. Initialize the number of individuals N and the population P(0) by random selection within the limits of the search space. Set the fitness function, constraints, maximum number of generations, mutation factor, cross over rate and tolerance.

Constraint handling
The main problem in finding the solution to the EED problem is that every new population evolved must satisfy the upper and lower bounds along with the nonlinear power balance constraint. The power balance constraint, being a polynomial of the individual solutions, makes it a complicated task to ensure convergence. Therefore, a continuous effort has to be made to restrict the solutions in the feasible area. Here, we propose a recursive distributed constraint handling method. The constraint c is the mismatch of power defined as where, ε is the maximum allowable constraint violation. The constraint violation is evaluated and distributed over the decision variables. The values of the N decision variables will change until the constraint violation is restricted within the permissible limits. Thus the solutions are confined within the feasible work space.
The proposed algorithm is implemented using sequences shown in the flowchart ( fig. 1). In the first block, the parameters including the limits are set. The second block initializes the population randomly between the upper and lower limits of the generating units. The population created is evaluated for the constraint violation; if the violation is high, then the population is updated using the constraint handler. The non-dominated sorting is performed, and Pareto optimal front is plotted from the updated population. This population represents the first iteration; subsequent iterations are performed to get further modified populations by applying equations (20), (21) and (22). This process is continued till convergence or the maximum number of iterations performed.

Performance measures
The following three performance metrics [21] have been applied to investigate the performance quality of the nondominated solutions obtained in the form of Pareto fronts using different algorithms.

Generation Distance (GD)
It is the estimated distance between the non-dominated solution vectors from the standard efficient front. It is mathematically expressed as (24) where, n is the number of solution vectors and di is the minimum Euclidean distance between each of them. GD=0 indicates that all the solution vectors are present in the standard Pareto front. A smaller positive value of GD means the Pareto front obtained from the proposed algorithm is closer to the standard Pareto front.

Spacing (S)
Spacing is the measure of the spread of the solution vectors. It is expressed as for all i, j=1, 2, …, n and d ̅ =mean of all d i and n is the number of non-dominated solution vectors found so far. The lower the value of the better is the Pareto solution.

Diversity metric (Δ)
It measures how evenly the solution vectors are distributed in the search space, i.e. extent of the spread on the Pareto front. It is found from the Euclidean distances as follows where, di is the distance between the consecutive solution vectors in the non-dominated solution set. The average of these distances is d ̅ .
Here,d f and d l are the first and last Euclidean distances. A low value of Δ indicates better diversity, and ∆=0 means the non-dominated solution vectors are uniformly spread on the Pareto front.

Simulation study and results
The EED problem was simulated for two different standard test cases i.e., IEEE 14 bus and IEEE 30 bus. The system data of these two test cases were obtained from the website www.ee.washington.edu/research/pstca. The cost and emission coefficients were also recorded from standard sources [22] and are presented in the Appendix along with the B coefficients. Each test case is solved for three different load demands. The scaling factor of differential evolution is set as 0.5 and crossover rate as 0.5 with a population size of 100.

Test case I: IEEE-14 bus system
The algorithms were applied to this test case for three different load demands; i.e., for 200 MW, 259 MW and 300 MW. The performance was compared with respect to Pareto optimal front, computational time, fuel cost, transmission loss involved and statistical performance metrics.

IEEE-14 bus system: load demand 200 MW
The Pareto optimal fronts obtained by applying the three algorithms are compared and shown in Fig.2 for a load demand of 200MW.
The comparative generation schedule, fuel cost, emission and transmission loss are presented in Table 1. The Pareto front obtained by applying the proposed MODE-RDC algorithm and other two algorithms for different load conditions are depicted from Fig.2 to Fig.  7. The fuel cost value (518.399$/hour) obtained is lowest among all the three competitive algorithms. The emission obtained from the proposed algorithm is also comparable with other two. Similarly, the algorithms were run for 25 times, and the performance metrics were calculated which are presented in Table 2. The average spacing obtained is 0.37959 which is the lowest among all the three algorithms. The average values of other two performance matrices show improved performance of the proposed algorithm. This indicates that the Pareto solutions obtained by the proposed algorithm are superior to the competing algorithms. Figure 3 shows the Pareto optimal fronts obtained for the test case with the load demand of 259 MW.      The generation schedule for this load is presented in Table 3, and the statistical performance metric values are shown in Table 4. The values of fuel cost and emission obtained by the proposed algorithm are better as compared to the other two algorithms. The generation schedule obtained by the former leads to less transmission loss. The proposed MODE-RDC algorithm shows lower average values of GD, S and ∆ as compared to the other two competitive algorithms.

IEEE-14 bus system: load demand 300 MW
With a higher load of 300 MW the Pareto optimal fronts obtained using the algorithms are shown in Fig. 4.
The generation schedule and performance metrics of the solution points are presented in Tables 5 and 6 respectively. The fuel cost and emission obtained using proposed MODE-RDC algorithm is better in terms of emission and PL in MW. The transmission losses involved due to the schedule obtained by the proposed algorithm is also lower. The average value of GD, S and ∆ metrics of the solutions arrived from the proposed algorithm are less than that of the other two. Thus the quality of solutions is maintained.

IEEE-14 bus system: summary of results
The summary of generation schedules obtained for the three load conditions using the three algorithms as presented in Tables 1, 3 and 5 show that the fuel cost for the proposed algorithm provides improved performance for the load conditions. For the IEEE 14 bus test case, the two-tailed Sign tests [21] is conducted for the pair-wise comparison of the algorithms by considering three different performance metrics, and are presented in Tables 7, 8 and 9. The tests were carried out for all the three load conditions and by repeated run of the three algorithms for 25 times. It is observed from Tables 7, 8 and 9 that the proposed algorithm wins over the other two for all loads in terms of all the three performance metrics i.e. the generation distance, spacing and diversity metric as winning parameters. It can be pointed in Table-10 that in Sign test, if more than 17 wins are recorded, then the algorithm is better with a detected difference α=0.1; and if more than 18 wins are recorded then, α=0.05.      It is evident from the Table 7 and 8 that the proposed MODE-RDC algorithm shows better performance as compared to other two competitive algorithms. The algorithm is better with respect to NSGA-II a detected difference α=0.1 for load of 200 MW, and a detected difference of α=0.05 for load of 259MW and 300 MW. However, the proposed algorithm does not show much improvement as compared to MOPSO algorithm.

Test case II: IEEE-30 bus system
The algorithms were applied to this test case for three different load demands; i.e., for 200 MW, 283.4 MW and 350 MW. The performance is compared with respect to Pareto optimal front, computational time, fuel cost, transmission loss involved and statistical performance metrics.

IEEE-30 bus system: load demand 200 MW
The Pareto optimal fronts obtained by applying the three algorithms for the load demand of 200 MW are presented in Figure 5. From the simulation output it reveals that the Pareto curve obtained by applying the proposed algorithm covers a wider area of the search space. The generation       schedule for an acceptable emission level is achieved, and is presented in Table 11. The quality of solutions is assessed by evaluating all the three performance metrics, and is presented in Table  12. The cost of fuel obtained by employing the proposed algorithm is lower, and the corresponding schedule causes lower transmission loss. It surpasses the performance of NSGA-II in terms of computational time. The lower average values of all the three performance matrices imply superior performance of the proposed algorithm over other two.

IEEE-30 bus system: load demand 283.4 MW
The test case is subjected to the rated load of 283.4 MW, and performance of the algorithms is assessed. The Pareto fronts obtained for this demand are presented in Figure 6. The Pareto curve corresponding to the proposed algorithm covers a wide area of the search space. The generation schedules obtained by applying the algorithms are shown in Table 13. The proposed algorithm yields lower values of fuel cost and emission. The transmission losses involved with this schedule is less than the other two algorithms. The quality of solutions is assessed by running the algorithms for 25 times; statistical values of the performance metrics thus obtained are presented in Table 14. The average diversity of the solutions using the proposed algorithm is less than that of the other two algorithms. It also requires less computational time as compared to the NSGA-II.

IEEE-30 bus system: load demand 350 MW
The test case is further subjected to a higher load of 350 MW, and performance of the algorithms is assessed. The Pareto optimal fronts are shown in Figure 7. The generating schedule for the system obtained for this demand is presented in Table 15.
The statistical behavior of the performance metrics obtained for the solutions are shown in Table 16. It is observed that the proposed algorithm performs better than NSGA-II in terms of diversity of solutions, computational time, fuel cost and emission values.

IEEE-30 bus system: summary of results
The summary of results for this test case with the three load demands by applying all the three algorithms is presented in Tables 11, 13 and 15. It shows that the fuel cost obtained using the proposed algorithm provides improved performance as compared to other two algorithms. The transmission loss involved with the generation schedule thus arrived is also lower. The quality of solutions arrived using the algorithms is further estimated from pair wise sign test [21] Tables 17 and 18. This is clear in these tables as it does not have a significant detected difference that needs at least 17 wins out of 25 runs. This is due to the reduction in the number of non-dominated solutions in the successive iterations in the proposed algorithm.
Here, it can be pointed out that in most cases the average value of S metric has increased with rise in load demand for all the three algorithms this can be verified from Tables 2, 4 and 6 for test case I and Tables 12, 14 and 16 for test case II. This is due to the fact that when the load demand on the system rises, the size of the feasible space decreases due to the generator limits. The situation is further limited by the power balance constraint as the transmission losses increase with increase in power generation. Thus, the feasible solutions fall apart causing increase in the spacing parameter.

Conclusion
A set of three multi-objective optimization algorithms have been applied to solve the EED problem for two test cases on three different load demands. The performance of the proposed MODE-RDC algorithm along with other two is assessed by considering three different performance metrics. The performances of these algorithms have been critically analyzed. The Pareto optimal fronts obtained by all the three algorithms incorporating the proposed recursive distributed constraint handling technique have sufficient diversity by exploiting the entire available range of search space. In particular, the Pareto front obtained by the multiobjective differential evolution with the recursive distributed constraint handling (MODE-RDC) approach has a better diversity in most cases. The spacing between the Pareto solutions has found to be increased with rise in the load demand on the system for all the three     algorithms. Moreover, the time requirement to achieve the Pareto front by applying the proposed recursive distributed constraint handling based technique is satisfactory.
Further research on this topic may include the inclusion of different evolutionary local search mechanisms into the approaches. It is expected to obtain lower computational speed and exploitation of the multidimensional search space. There is a need for further investigation to explore the strengths and weaknesses of the proposed algorithm, so that it can be applied to other multi-objective problems in power systems, such as management of voltage profiles, reactive power compensation etc. The performance of the proposed algorithm can also be investigated by considering other real world constraints like ramp rate limits, power loss etc.

Appendix
The standard test cases taken for the solution are IEEE 14 bus and IEEE 30 bus power systems. The parameters of the test cases have been adopted from standard sources. The cost and emission coefficients used for solution of the problem are shown in tables 20 and 21 below.
The values of B coefficients used for the IEEE 14 bus test case are