Abstract

Efficient and reliable power production is necessary to meet both the profitability of power systems operations and the electricity demand, taking also into account the environmental concerns about the emissions produced by fossil-fuelled power plants. The economic emission load dispatch problem has been defined and applied in order to deal with the optimization of these two conflicting objectives, that is, the minimization of both fuel cost and emission of generating units. This paper introduces and describes a solution to this famous problem using a new metaheuristic nature-inspired algorithm, called firefly algorithm, which was developed by Dr. Xin-She Yang at Cambridge University in 2007. A general formulation of this algorithm is presented together with an analytical mathematical modeling to solve this problem by a single equivalent objective function. The results are compared with those obtained by alternative techniques proposed by the literature in order to show that it is capable of yielding good optimal solutions with proper selection of control parameters.

1. Introduction

Biology-inspired metaheuristic algorithms have recently become the forefront of the current research as an efficient way to deal with many NP-hard combinatorial optimization problems and non-linear optimization constrained problems in general. These algorithms are based on a particular successful mechanism of a biological phenomenon of Mother Nature in order to achieve optimization, such as the family of honey-bee algorithms, where the finding of an optimal solution is based on the foraging and storing the maximum amount of flowers’ nectar [1]. A new algorithm that belongs in this category of the so-called nature inspired algorithms is the firefly algorithm which is based on the flashing light of fireflies.

Although the real purpose and the details of this complex biochemical process of producing this flashing light is still a debating issue in the scientific community, many researchers believe that it helps fireflies for finding mates, protecting themselves from their predators and attracting their potential prey [14]. In the firefly algorithm, the objective function of a given optimization problem is associated with this flashing light or light intensity which helps the swarm of fireflies to move to brighter and more attractive locations in order to obtain efficient optimal solutions.

In this research paper we will show how the recently developed firefly algorithm can be used to solve the famous economic emissions load dispatch optimization problem. This hard optimization problem constitutes one of the key problems in power system operation and planning in which a direct solution cannot be found and therefore metaheuristic approaches, such as the firefly algorithm, have to be used to find the near optimal solutions. This optimization problem deals with allocating loads to power generators of a plant for minimum total fuel cost and emissions while meeting the power demand and transmission losses constraints. There are numerous variations of this problem which model the two objective functions and the constraints in many different ways.

Moreover, we will demonstrate how the firefly algorithm works and how this method can be easily adapted in order to solve this multiobjective optimization problem. Therefore, we will discuss why this method is sufficiently accurate and easy to implement for real-time operation and control of power systems. For the efficiency and validation of this algorithm, we will use, as an example, a sample realistic test system having six power generators. We will also compare the solutions obtained with the ones obtained by alternative optimization techniques that have been successfully applied by many scientists in order to solve these types of problems, such as the goal attainment SQP method, ant colony optimization, and particle swarm optimization [510].

The remainder of this paper is organized as follows: Section 2 gives a brief description of the multiobjective optimization and why this is important in our case. In Section 3, a mathematical formulation and description of the economic emissions load dispatch problem is given. Section 4 gives a brief description of the goal attainment SQP method which was used as an alternative way to solve an example test system of the economic emission load dispatch problem, while Section 5 briefly describes the Firefly algorithm in general. Section 6 presents the computational results simulated in Matlab and acquired when we applied both the firefly algorithm and the goal attainment SQP method to solve the economic emission load dispatch problem. Moreover, the efficiency of the firefly algorithm is measured by comparing its results with those obtained by other stochastic algorithms proposed by the researchers for different test systems of this problem. Finally, Section 7 provides some conclusions concerning the solutions obtained by the firefly algorithm and some suggestions and ideas for further research.

2. Multiobjective Optimization and Problem Formulation

Multiple conflicting objectives arise naturally in most real-world combinatorial optimization problems, such as the economic emissions load dispatch problem. Several principles and strategies have been developed and proposed for over a decade in order to solve these problems, some of which will be discussed in this section. In the multiobjective optimization problem, as its name implies, we have multiple objective functions with a possibility of conflict with each other. The aim is to find a vector of decision variables that satisfies constraints and optimizes (minimizes or maximizes) these functions. In such cases, we have to construct an overall objective function as a linear combination of the conflicting multiple objective functions using a weighting factor for each function. In a more precise mathematical way, the multiobjective optimization problem can be defined as follows [1, 7, 1014]:𝑥Findvector𝑋=1,𝑥2,,𝑥𝑛𝑇𝑓Ωwhichoptimizes𝑓(𝑥)=1(𝑥),𝑓2(𝑥),,𝑓𝑖(𝑥)subjectto𝑔𝑗(𝑥)0or𝑔𝑗(𝑥)0,𝑗=1,2,𝑚,(2.1) where 𝑓1,𝑓2,,𝑓𝑖 denote the objective functions to be optimized simultaneously, 𝑋 is the vector of discrete decision variables or search/decision space, Ω is the finite set of feasible solutions and 𝑔𝑗(𝑥) denotes the inequality constraints. The functions 𝑓𝑖(𝑥) and 𝑔𝑖(𝑥) may be linear or non-linear. The multiobjective optimization problem is sometimes called vector minimization problem.

2.1. Pareto Optimal Solutions

For any given problem having more than one objective functions, any two solutions 𝑥 and 𝑦 of this problem can have one of two distinct possibilities: one solution may dominate over the other or none of them dominates over the other, since there can be no solution vector 𝑋 that minimizes all the 𝑖 objective functions simultaneously. Therefore, we introduce the well known Pareto optimum solution concept in multiobjective optimization problems. A feasible solution 𝑋 is called Pareto optimal solution if there exists no other feasible solution 𝑌 such that 𝑓𝑖(𝑌)𝑓𝑖(𝑋) for all 𝑖=1,2,,𝑘with 𝑓𝑗(𝑌)<𝑓𝑖(𝑋) for at least one 𝑗, which means that the solution 𝑋 is no worse than 𝑌 in all objective functions, and the solution 𝑋 is strictly better than 𝑌 in at least one objective function [10, 11, 13, 14]. In other words, the solutions that are nondominated within the entire search space are denoted as Pareto optimal solutions and constitute the Pareto optimal set or Pareto optimal frontier (i.e., the image of the Pareto set in objective space). The knowledge of this set is crucial in many multiobjective optimization problems, as this enables the decision maker choose the best compromise solution [11, 13, 14].

As it is very difficult to effectively handle with all the conflicting objective functions, several methods have been developed for this purpose, such as the utility function method and the goal attainment method. In most of these methods, the multiobjective problem is transformed into a single-objective problem, then a set of Pareto optimal solutions is generated, and some additional criterion or rule to select one particular Pareto optimal solution is used as a solution of the multiobjective problem.

2.2. The Utility Function Method

In this paper, in order to apply and solve the economic emissions load dispatch problem to the firefly algorithm, we use the utility function method or weighting function method as an efficient way to deal with conflicting discrete goals and combine the two conflicting objective functions of the problem into one objective function. In this method, the problem is transformed into a single-objective function problem by using a utility function 𝑈𝑖(𝑓𝑖) for each objective function based on its importance compared to the other objective functions. The overall utility function of the problem is defined as [1214]𝑈=𝑘𝑖=1U𝑖𝑓𝑖.(2.2)

Thus, the solution vector 𝑋 is found by maximizing the total utility function 𝑈 subject to the constraints 𝑔𝑗(𝑥)0, 𝑗=1,2,𝑚. Then, the previous equation can be rewritten as𝑈=𝑘𝑖=1𝑈𝑖=𝑘𝑖=1𝑤𝑖𝑓𝑖(𝑋),(2.3) where 𝑤𝑖 is a scalar weighting factor associated with the 𝑖th objective function and 𝑤1+𝑤2++𝑤𝑘=1. In case of 𝑔𝑖(𝑥)0, then 𝑤1+𝑤2++𝑤𝑘=1.

Since the economic emissions load dispatch problem has two conflicting objective functions, the problem can be defined as 𝑓(𝑥)=𝑎𝑓1(𝑥)+(1𝑎)𝑓2(𝑥), where 𝑓1(𝑥) and 𝑓2(𝑥) represent these two conflicting objectives, and a is the weighting coefficient, as we will see in Section 6. However, in all the tables of this paper, we will not use the weights in the presentation of the values of the two objectives.

3. The Economic Emissions Load Dispatch Problem

The economic emission load dispatch problem is considered a hard optimization problem which normally renders classical exact optimization methods ineffective. For this reason, many approximate algorithms and particularly evolutionary methods such as artificial immune systems, genetic algorithms, cultural algorithms, particle swarm optimization, and differential evolution, have gained a lot of popularity for efficiently solving this problem [5, 6, 9, 10, 13, 1518]. Furthermore, many hybrid methods in which two or more methodologies are combined together in order to enhance the final model have been introduced and attracted much attention by many researchers in the last few years [6, 810, 13, 1517]. In fact, some of these proposed hybrid approaches are considered very effective, robust, and well promising so as to deal with other hard combinatorial optimization problems.

This problem is a multiple, conflicting objective function problem. Its primary objective is to determine the optimal quantity of both energy and emission objectives to reliably serve consumers. In particular, the main objective is the minimization of both the fuel cost of operation and emission of each power unit, subject to satisfying the equality and inequality constraints concerning operational limits of generations, the load demand, and the transmission losses of the facilities. Notably, in the recent years, many countries in the world have agreed to decrease the amount of pollutants fossil fuel power generating units taking into account the environmental concerns and protection policies about the emissions produced by these plants. However, this decrease of the total emission from a system plant implies an increase in the system operating cost, justifying the high applicability and paramount importance of this optimization problem. Therefore, the objective of this problem is to find out an operating point, which keeps a balance between cost and emission. This can be formulated mathematically as a minimization multiobjective problem as follows [5, 6, 810, 13, 1519].

3.1. The Fuel Cost Objective

The aim is to minimize the total fuel cost (operating cost) of all committed plants can be stated as follows:min𝑓1(𝑋)=𝑛𝑖=1𝐶𝑖𝑃𝐺𝑖=𝑛𝑖=1𝛼i+𝑏𝑖𝑃𝐺𝑖+𝑐𝑖𝑃2𝐺𝑖,(3.1) where 𝑛 is the number of units power generators of a power plant, 𝐶𝑖 is the fuel cost of the 𝑖th generator, 𝑃𝐺𝑖 is the out power of generator 𝑖 and 𝑎𝑖, and 𝑏𝑖 and 𝑐𝑖 are the fuel cost coefficients of the 𝑖th generator. Normally, the fuel cost equation 𝑓1(𝑋) is expressed as continuous quadratic (higher order) equation, as we have here, but sometimes it can be expressed in linear form, when the coefficient 𝑐𝑖 is equal to zero. However, in both cases, the equation expresses the variation of fuel cost ($ or Rs) with generated power or time (MW or hr). In our paper we use the $/hr as the only unit of measurement of the fuel cost function.

3.2. The Emissions Objective

The aim is to minimize the sum of all types of emissions from fossil power plants, namely, the gaseous particle pollutants, such as NOx, SOx, and CO2, thermal and any other chemical emissions with suitable emission weights/coefficients on each pollutant power generator. In our research, without loss of generality, we will consider only one type of emission, the amount of NOx, which is formulated in the following quadratic equation:min𝑓2(𝑋)=𝑛𝑖=1102𝛼𝑖+𝛽𝑖𝑃𝐺𝑖+𝛾𝑖𝑃2𝐺𝑖+𝜁𝑖𝜆exp𝑖𝑃𝐺𝑖,(3.2) where, 𝛼𝑖, 𝑏𝑖, 𝛾𝑖, 𝜁𝑖, and 𝜆𝑖 are the emission coefficients of NOx for the 𝑖th power generator. In our research, we use the unit of measurement of ton/hr for the emissions function (some research papers use kg/hr).

3.3. The Necessary Constraints of the Problem

The total power generation must satisfy the total required demand (power balance) and transmission losses. This can be formulated as follows:𝑛𝑖=1𝑃𝐺𝑖=𝐷+𝑃loss,(3.3) where 𝐷 is the real total load demand of the system, 𝑃𝐺𝑖 is the 𝑖th generator’s power, and 𝑃loss is the transmission losses. These can be determined from either the load/power flow or the matrix 𝐵𝑖𝑗 of coefficients. In this paper, only the 𝐵𝑖𝑗 coefficients are considered𝑃loss=𝑛𝑖𝑛𝑗𝐵𝑖𝑗𝑃𝑖𝑃𝑗,(3.4) where, 𝐵𝑖𝑗 are the elements of the loss coefficient matrix B and 𝑃𝑖 and 𝑃𝑗 are the out powers of the 𝑖th and 𝑗th generator; respectively. In our paper, we use the MW as the only unit of measurement of the power balance constrain.

Apart from the total demand and transmission loss constrain, there is also the generator capacity constrain in which the power limits of each generator are formulated in order to have a stable operation of a plant. The upper and lower limits are defined as follows:𝑃MIN𝐺𝑖𝑃𝐺𝑖𝑃MAX𝐺𝑖,for𝑖=1𝑛,(3.5) where 𝑃MIN𝐺𝑖and 𝑃MAX𝐺𝑖 are the lower and upper limit of the 𝑖th generator’s out power 𝑃𝐺𝑖, respectively. The power load of each generator unit is measured in MW. As we will see later, in order incorporate these conditions into one equation, we will use an auxiliary (Lagrange) function 𝑓3(𝑥).

4. The Goal Attainment SQP Method

The goal attainment SQP method is actually a hybrid method which combines the Goal Attainment (GA) method with the Sequential Quadratic Programming (SQP). In particular, a multiobjective non-linear optimization problem is initially reformulated by the goal attainment method so as it can be used by the Sequential Quadratic Programming method (SQP) in order to optimize it. This hybrid method has been successfully applied by many scientists as an effective way to solve many optimization problems. It constitutes a highly effective and state-of-the-art technique to obtain the best compromise solution in a multiobjective problem [7]. This method is already implemented in Matlab 2008 and belongs to the standard optimization toolbox multiobjective solver of Matlab. We will apply this method in a realistic test system, described in Section 6, in order to compare and contrast the solutions obtained by this alternative method and those obtained by the firefly algorithm. By this way, we will prove the effectiveness and robustness of the proposed algorithm.

4.1. The Goal Attainment Method

The main idea behind this method is to reduce a set of non-linear functions 𝑓𝑖(𝑋) below a set of goals of 𝑓𝑖. However, since all goals cannot be achieved, this set of functions is reformulated in an appropriate, well-defined way and enabling the designer to be relatively imprecise about the initial goals [6, 8, 14, 16]. In particular, in this method, we define a set of desired design goals 𝑏=[𝑏1,𝑏2,,𝑏𝑛]𝑇 which is associated with a set of objective functions 𝑓(𝑥)={𝑓1(𝑥),𝑓2(𝑥),,𝑓𝑛(𝑥)}. The relative degree of underachievement or overachievement of the desired goals is controlled by a vector of weighting coefficients, 𝑤=[𝑤1,𝑤2,,𝑤𝑛] for each objective function, which express the importance of the 𝑖th objective function relative to the other objective functions in meeting the respective 𝑖th goal. The goal 𝑏𝑖 is found by first solving the single-objective function optimization problem [14]minimize𝑓𝑖(𝑋),𝑖=1,2,,𝑛subjectto𝑔𝑗(𝑋)0or𝑔𝑗(𝑋)0,𝑗=1,2,𝑚,(4.1) where 𝑋=(𝑥1,𝑥2,,𝑥𝑛)𝑇Ω are the design variables. In this case, the goal 𝑏𝑖 is taken as the optimum value of the objective 𝑓𝑖(𝑋), that is, 𝑏𝑖=𝑓𝑖=𝑓(𝑋𝑖), where 𝑋𝑖 denotes the solution of the problem defined in (4.1).

In order to find the best compromise solution of the problem, we introduce a scalar variable 𝛾 as a design variable in addition to 𝑛 design variables 𝑥𝑖, for 𝑖=1,2,,𝑛. Then, the problem is solved as a standard optimization problem using the following equation [68, 14]:𝑥Find1,𝑥2,,𝑥𝑛𝑥Ω,𝛾𝑅tominimize𝑓1,𝑥2,,𝑥𝑛,𝛾=𝛾subjectto𝑔𝑗(𝑋)0or𝑔𝑗(𝑋)0,𝑗=1,2,𝑚,𝑓𝑖(𝑋)𝛾𝑤𝑖𝑏𝑖for𝑔𝑗(𝑋)0or𝑓𝑖(𝑋)+𝛾𝑤𝑖𝑏𝑖for𝑔𝑗(𝑋)0,𝑖=1,2,,𝑛,(4.2) where 𝑋=(𝑥1,𝑥2,,𝑥𝑛)𝑇Ω are the design variables (parameters), 𝛾 is a scalar variable unrestricted in sign, and the weights 𝑤=[𝑤1,𝑤2,,𝑤𝑛]𝑅𝑛 satisfy the necessary normalization condition which in case of 𝑔𝑗(𝑋)0 the condition is 𝑤1+𝑤2++𝑤𝑘=1 or in case of 𝑔𝑖(𝑥)0 the condition is 𝑤1+𝑤2++𝑤𝑘=1.

Moreover, the term 𝛾𝑤𝑖 introduces the degree of slackness into the problem, which otherwise imposes the goals to be rigidly met. The weighting vector 𝑤 enables the decision maker to quantitatively express a measure of the relative tradeoffs between the objectives [68]. For example, if we set the weighting vector equal to the initial goals, we assume equal underachievement and overachievement of the goals 𝑏𝑖. We can also include hard constraints into the design by setting a particular weighting factor to zero. In particular, if a component of vector w is equal to zero (i.e., some 𝑤𝑖=0), this means that the maximum limit of an objective 𝑓𝑖(𝑋) is equal to 𝑏𝑖. We can easily see that a set of solutions can be generated by varying 𝑤 over 𝑅𝑛 even for nonconvex problems, and therefore the multiobjective optimization involves generating and selecting solutions that characterize the objectives, where the improvement of one objective causes a declination in another objective.

During this optimization procedure, the scalar variable 𝛾 is varied changing the size of the feasible region. It is interesting to mention that the optimum value of 𝛾, which is also called attainment factor, will inform the decision maker whether the decision goals are attainable or not. A negative value of 𝛾 implies that the goal defined by the decision maker is attainable and an improved solution can be obtained. If the value of 𝛾 is positive, then the goal defined by the decision maker is unattainable, and no improvement of the solution is possible.

4.2. The Sequential Quadratic Programming Method

The Sequential Quadratic Programming (SQP) method constitutes one of the most popular, efficient, and accurate algorithms for non-linear continuous optimization which has been implemented and tested over a large number of several combinatorial optimization problems [7]. This method is very closely related to Newton’s method for unconstrained optimization. In particular, at each iteration, an approximation local model of the problem is constructed and solved providing the necessary direction for finding the solution of the originally defined problem. In the case of an unconstrained minimization, only the objective function must be approximated, then the local model is quadratic and the algorithm reduces to the famous Newton's method.

Given a problem defined as in (4.1), both the objective function and constraints have to be modeled. The objective function is formulated as a quadratic model while constraints are formulated as a linear model. Both of these models are based on the quadratic approximation of the Langragian function of the problem which is updated using a quasi-Newton updating method. The Langragian function of the problem defined in (4.1) is given as follows [1, 7, 14]:𝐿(𝑋,𝜆)=𝑓𝑖(𝑋)𝑚𝑗=1𝜆𝑗𝑔𝑗(𝑋),(4.3) where, 𝜆𝑗 is a Langrage multiplier, 𝑔𝑗(𝑋)0 and for 𝑖=1,2,,𝑛. In case of 𝑔𝑗(𝑋)0, we have a plus sign instead of minus before the sum.

Now, we can define the quadratic programming model as follows [1, 7, 14]:min𝑑𝐿𝑥𝑘,𝜆𝑘𝑥+𝐿𝑘,𝜆𝑘𝑇1𝑑+2𝑑𝑇2𝐿𝑥𝑘,𝜆𝑘𝑥𝑑,s.t.𝑔𝑘𝑥+𝑔𝑘𝑇𝑥𝑑0for𝑔𝑘0,(4.4) where, 𝑥𝑘 is the current estimate of a solution 𝑥, and 𝑑 is the solution of the quadratic programming model which is used as a search direction for finding the appropriate solution of the given problem.

5. The Firefly Algorithm

5.1. Description

The Firefly Algorithm (FA) is a metaheuristic, nature-inspired, optimization algorithm which is based on the social (flashing) behavior of fireflies, or lighting bugs, in the summer sky in the tropical temperature regions [13, 20]. It was developed by Dr. Xin-She Yang at Cambridge University in 2007, and it is based on the swarm behavior such as fish, insects, or bird schooling in nature. In particular, although the firefly algorithm has many similarities with other algorithms which are based on the so-called swarm intelligence, such as the famous Particle Swarm Optimization (PSO), Artificial Bee Colony optimization (ABC), and Bacterial Foraging (BFA) algorithms, it is indeed much simpler both in concept and implementation [24, 20]. Furthermore, according to recent bibliography, the algorithm is very efficient and can outperform other conventional algorithms, such as genetic algorithms, for solving many optimization problems; a fact that has been justified in a recent research, where the statistical performance of the firefly algorithm was measured against other well-known optimization algorithms using various standard stochastic test functions [13, 20]. Its main advantage is the fact that it uses mainly real random numbers, and it is based on the global communication among the swarming particles (i.e., the fireflies), and as a result, it seems more effective in multiobjective optimization such as the economic emissions load dispatch problem in our case.

The firefly algorithm has three particular idealized rules which are based on some of the major flashing characteristics of real fireflies [24, 20]. These are the following: (1) all fireflies are unisex, and they will move towards more attractive and brighter ones regardless their sex. (2) The degree of attractiveness of a firefly is proportional to its brightness which decreases as the distance from the other firefly increases due to the fact that the air absorbs light. If there is not a brighter or more attractive firefly than a particular one, it will then move randomly. (3) The brightness or light intensity of a firefly is determined by the value of the objective function of a given problem. For maximization problems, the light intensity is proportional to the value of the objective function.

5.2. Attractiveness

In the firefly algorithm, the form of attractiveness function of a firefly is the following monotonically decreasing function [2, 3, 20]:𝛽(𝑟)=𝛽0exp(𝛾𝑟𝑚),with𝑚1,(5.1) where, 𝑟 is the distance between any two fireflies, 𝛽0 is the initial attractiveness at 𝑟=0, and 𝛾 is an absorption coefficient which controls the decrease of the light intensity.

5.3. Distance

The distance between any two fireflies 𝑖 and 𝑗, at positions 𝑥𝑖 and 𝑥𝑗, respectively, can be defined as a Cartesian or Euclidean distance as follows [2, 3, 20]:𝑟𝑖𝑗=𝑥𝑖𝑥𝑗=𝑑𝑘=1𝑥𝑖,𝑘𝑥𝑗,𝑘2,(5.2) where 𝑥𝑖,𝑘 is the 𝑘th component of the spatial coordinate 𝑥𝑖 of the 𝑖th firefly and 𝑑 is the number of dimensions we have, for 𝑑=2, we have𝑟𝑖𝑗=𝑥𝑖𝑥𝑗2+𝑦𝑖𝑦𝑗2.(5.3)

However, the calculation of distance 𝑟 can also be defined using other distance metrics, based on the nature of the problem, such as Manhattan distance or Mahalanobis distance.

5.4. Movement

The movement of a firefly 𝑖 which is attracted by a more attractive (i.e., brighter) firefly 𝑗 is given by the following equation [2, 3, 20]:𝑥𝑖=𝑥𝑖+𝛽0exp𝛾𝑟2𝑖𝑗𝑥𝑗𝑥𝑖1+𝑎rand2,(5.4) where the first term is the current position of a firefly, the second term is used for considering a firefly’s attractiveness to light intensity seen by adjacent fireflies, and the third term is used for the random movement of a firefly in case there are not any brighter ones. The coefficient 𝛼 is a randomization parameter determined by the problem of interest, while rand is a random number generator uniformly distributed in the space [0,1]. As we will see in this implementation of the algorithm, we will use 𝛽0=1.0, 𝛼[0,1] and the attractiveness or absorption coefficient 𝛾=1.0, which guarantees a quick convergence of the algorithm to the optimal solution.

5.5. Convergence and Asymptotic Behavior

The convergence of the algorithm is achieved for any large number of fireflies (𝑛) if 𝑛𝑚, where 𝑚 is the number of local optima of an optimization problem [1, 3]. In this case, the initial location of 𝑛 fireflies is distributed uniformly in the entire search space. The convergence of the algorithm into all the local and global optima is achieved, as the iterations of the algorithm continue, by comparing the best solutions of each iteration with these optima. However, it is under research a formal proof of the convergence of the algorithm and particularly that the algorithm will approach global optima when 𝑛 and 𝑡1 [3]. In practice, the algorithm converges very quickly in less than 80 iterations and less than 50 fireflies, as it is demonstrated in several research papers using some standard test functions [13, 20]. Indeed, the appropriate choice of the number of iterations together with the 𝛾, 𝛽, 𝛼, and 𝑛 parameters highly depends on the nature of the given optimization problem as this affects the convergence of the algorithm and the efficient find of both local and global optima. Note that the firefly algorithm has computational complexity of 𝑂(𝑛2), where 𝑛 is the population of fireflies. The larger population size becomes the greater the computational time is [13].

5.6. Special Cases

There are two important special cases of the firefly algorithm based on the absorption coefficient 𝛾; that is, when 𝛾0 and 𝛾 [1, 3, 20]. When 𝛾0, the attractiveness coefficient is constant 𝛽=𝛽0, and the light intensity does not decrease as the distance 𝑟 between two fireflies increases. Therefore, as the light of a firefly can be seen anywhere, a single local or global optimum can be easily reached. This limiting case corresponds to the standard Particle Swarm Optimization (PSO) algorithm. On the other hand, when 𝛾, the attractiveness coefficient is the Dirac delta function 𝛽(𝑟)𝛿(𝑟). In this limiting case, the attractiveness to light intensity is almost zero, and as a result, the fireflies cannot see each other, and they move completely randomly in a foggy place. Therefore, this method corresponds to a random search method.

5.7. Hybridization

In a recent bibliography, a new metaheuristic algorithm has been developed and formulated based on the concept of hybridizing the firefly algorithm. In particular, the new Levy flight Firefly algorithm was developed by Dr. Xin-She Yang at Cambridge University in 2010 and it combines the firefly algorithm with the Levy flights as an efficient search strategy [4]. It combines the three idealized rules of the firefly algorithm together with the characteristics of Levy flights which simulate the flight behavior of many animals and insects. In this algorithm, the form of the attractiveness function and the calculation of distance between two fireflies are the same as in firefly algorithm, but in the movement function, the random step length is a combination of the randomization parameter together with a Levy flight. In particular, the movement of a firefly is a random walk, where the step length is drawn by the Levy distribution [4].

6. Case Study

6.1. Application Example

In order to apply the firefly algorithm to the economic emissions load dispatch problem, we need to effectively deal with the necessary constraints of the problem and the fact that the firefly algorithm can directly solve only maximum optimization problems, not minimization problems. This is expected as the light intensity or brightness of a firefly at a particular location 𝑥 can be chosen as analogous or equal to the value of the objective function at location 𝑥.

For this reason, in order to avoid the violation of constraints, which could cause infeasible solutions, we have converted the constrained optimization problem into an unconstrained problem by penalizing infeasible solutions of the objective function, instead of repairing them, which is very computationally expensive in our case [1, 11, 14]. Therefore, we used the method of Lagrange multipliers as an efficient strategy for finding the maxima (orminima) of the given function subject to the constraints of total required demand and transmission losses (Section 3.3). The Lagrange multiplier we used as a weighting factor in the Lagrange (constraints) function 𝑓3(𝑋) is 1000. This value was chosen experimentally as a parameter that maximizes 𝑓3(𝑋). Furthermore, without loss of generality, in an optimization problem, the minimum of a function can be found by seeking the maximum of the negative of the same function, or alternatively, we can subtract from that function a large positive number which is usually the global optimum of this function [1, 11, 12, 14]. In this problem we adopt the first choice, where we maximize the negative form of the objective function. Finally, we use 𝑤1=𝑤2=𝑎=0.5 as a weighting factor for each objective function, since both of the objectives (i.e., fuel cost and emissions) have equal importance in the given problem, and therefore, they have to contribute equally to the formulation of the objective function. Moreover, the choice of equal weights is also justified by some experiments we conducted with different weights, in which we did not achieve better results, as the ranges of the values of the two objective functions 𝑓1(𝑋) and 𝑓2(𝑋) are very different (Figures 3 and 6). Now the problem can be formulated as follows:max𝑓(𝑋)=max0.5𝑓1(𝑋)0.5𝑓2(𝑋)𝑓3(𝑋),(6.1) where 𝑓1(𝑋) is the fuel cost objective function defined in the following equation (in $/hr):𝑓1(𝑋)=6𝑖=1𝛼i+𝑏𝑖𝑃𝐺𝑖+𝑐𝑖𝑃2𝐺𝑖,(6.2) where 𝑓2(𝑋) is the emissions objective function defined in the following equation (in ton/hr):𝑓2(𝑋)=6𝑖=1102𝛼𝑖+𝛽𝑖𝑃𝐺𝑖+𝛾𝑖𝑃2𝐺𝑖+𝜁𝑖𝜆exp𝑖𝑃𝐺𝑖,(6.3) and 𝑓3(𝑋) is the constraints equation which is defined as follows (in MW):𝑓3(𝑋)=1000abs6𝑖=1𝑃𝐺𝑖𝐷6𝑖=16𝑗=1𝐵𝑖𝑗𝑃𝐺𝑖𝑃𝐺𝑗,(6.4)

It can be easily seen that the generalized form of the objective function given in (6.1) is non-linear, and the number of equality and inequality constraints of the problem increase with the size of the power test system we use. Therefore, the application of conventional straightforward optimization techniques, such as the gradient-based methods, is highly unsuccessful as the solution of the non-linear objective functions with many inequality constraints of the problem is very computationally expensive and the search space inefficiently large to explore.

In this test problem, we consider a plant with 6 power generators. The test system was derived from [10]. The generation data are given in Tables 1 and 2. The total power system load demand 𝐷 is 1 MW. In particular, the values of the fuel cost and emission coefficients together with the power limits of each generator are given in Table 1 of the next page.

In order to test and demonstrate the efficiency of the algorithm and how it deals with different problem complexities and tradeoff surfaces we have considered two different cases. In case A the test problem is considered as lossless; that is, the transmission power loss matrix 𝐵 is not considered. In case B, the power loss is taken into account and the transmission power loss coefficients matrix 𝐵 is defined in Table 2 of the following page. Notice that the power loss coefficients matrix 𝐵 is not symmetric (i.e., 𝐵46𝐵64), which depicts a realistic sample system for application. The reported results of both of these cases are discussed and analyzed in the following Section 6.3.

6.2. Application of the Proposed Firefly Algorithm

In order to solve the economic emissions load dispatch problem, we have implemented the firefly algorithm in Matlab 2008 and it was run on a portable computer with an Intel Core2 Duo (1.8 GHz) processor, 2 GB RAM memory and MS Windows 7 as an operating system. Mathematical calculations and comparisons can be done very quickly and effectively with Matlab and that is the reason that the proposed Firefly algorithm was implemented in Matlab 2008 programming environment. In this proposed method, we represent and associate each firefly with a valid power output (i.e., potential solution) encoded as a real number for each power generator unit, while the combined fuel cost and emission objective (i.e., the objective function of the problem) is associated and represented by the light intensity of the fireflies. In this simulation, the values of the control parameters are: 𝛼=0.2, 𝛾=1.0, 𝛽0=1.0, and 𝑛=12, and the maximum generation of fireflies (iterations) is 50. The values of the fuel cost, the emission coefficients, the power limits of each generator, the power loss coefficients, and the total power load demand are supplied as inputs to the firefly algorithm. The power output of each generator, the total system power, the fuel cost, and the emissions with/without transmission losses are considered as outputs of the proposed Firefly algorithm. The pseudocode of this algorithm can be summarized in Pseudocode 1 of the next page [1, 3, 20].

Pseudocode 1 gives a brief description and a practical application of the proposed algorithm for the economic emission load dispatch problem. Initially, the objective function of the given problem is formulated as defined in (6.1) and it is associated with the light intensity of the swarm of the fireflies. The initial solution of the given problem is generated based on the mathematical formulation given below:𝑥𝑗=rand(upper-rangelower-range)+lower-range,(6.5) where 𝑥𝑗 is the new solution of 𝑗th firefly, that is, created, rand is a random number generator uniformly distributed in the space [0,1], while upper range and lower range are the upper range and lower range of the 𝑗th firefly (variable), respectively.

After the evaluation of the initial population/generation (i.e., solution), the firefly algorithm enters its main loop which represents the maximum number of generations of the fireflies. This is actually the termination criterion that needs to be satisfied for the termination of the loop.

The generation of a new solution (i.e., the movement of a firefly) of the given problem is made based on the following mathematical formulation:𝑥𝑖=𝑥𝑖+𝛽0exp𝛾6𝑖,𝑗=1𝑥𝑖𝑥𝑗2𝑥𝑗𝑥𝑖1+𝑎rand2,(6.6) where 𝑥𝑖 is the current solution of the 𝑖th firefly and 𝑥𝑗 is the current (optimal) solution of 𝑗th firefly. The values of the algorithm’s control parameters is 𝛼=0.2, 𝛾=1.0, 𝛽0=1.0, and rand is a random number which is uniformly distributed in the space [0,1]. As we can see the distance between two fireflies is calculated using the Euclidean distance (Section 5.3) and the generation of a new solution is actually a sum of the current solution (𝑥𝑖), the metric of the evaluation of the current solution based on the current optimal solution (Euclidian metric), and a random step/move of the algorithm (Section 5.4).

After the generation of the new solutions, we have to apply the generator capacity constraints so as the new solutions are within the given operational power ranges. To avoid such violation, a repair process is applied to each solution (firefly) in order to guarantee that the generated power outputs are feasible. 𝑃𝐺𝑖, 𝑃𝑖min and 𝑃𝑖max denote the current, the minimum, and the maximum power outputs of the 𝑖th unit, which is associated with the 𝑖th firefly. The implemented procedure is shown in Pseudocode 2 of the previous page.

Finally, it is notable that for each generation (iteration), the swarm of 12 fireflies is ranked based on their light intensity, and the firefly with the maximum light intensity (i.e., the solution with the higher objective function value) is chosen as the brighter one (i.e., it is a potential optimal solution), while the others are updated based on (6.6). In the final iteration, the firefly with the brighter light intensity among the swarm of 12 fireflies is chosen as the brightest one which represents the optimal solution of the problem.

6.3. Results and Discussion

The main characteristic feature of the firefly algorithm is the fact that it simulates a parallel independent run strategy, where in every iteration, a swarm of n fireflies (12 in our case) has generated n solutions. Each firefly works almost independently, and as a result the algorithm, will converge very quickly with the fireflies aggregating closely to the optimal solution [1, 3, 20]. In our case, we can see that the total number of 600 function evaluations (12 fireflies * 50 generations) is sufficient, and as a result, the algorithm stably converges to the optimal solution very quickly (approximately from the 10th generation/iteration). This algorithm also differs from other alternative approaches in the selection procedure in which each firefly constructs its own solution based on a weighted sum of the objective functions, where the weight attached to multiobjective criteria is constant. It is also observed that the proposed implementation of the firefly algorithm is very fast and predicts accurate results while satisfying inequality generator capacity constraints at various power load levels. It also offers a considerable saving in computer memory. The algorithm generated the optimal solution in less than 5 seconds, on the personal computer we used with specifications discussed in the previous section. From the results we obtained, it is clear that the firefly algorithm gives a global optimum in a computation time of less than 3 seconds, which is much lower than the time of other alternative techniques, such as the goal attainment SQP method, where the computation time was approximately 5 seconds. From the simulation results, as we show in Table 3 of the next page, we can conclude that the algorithm is very fast, efficient, and computational inexpensive in finding the Pareto optimal solutions of the given combinatorial optimization problem.

6.3.1. Application to the Test System

The results in terms of values of the power losses, the best fuel cost, and the best emission objectives of both the proposed firefly algorithm and the goal attainment SQP method can be summarized in the following Table 3 of the next page. Two different cases (A and B) have been considered as discussed earlier, and in both cases, the proposed algorithm has been demonstrated through a sample test system consisting of six power generators with load demand 𝐷 of 1 MW and equal weights of the two objectives (𝑤1=𝑤2=𝑎=0.5) with and without transmission losses B.

The Obj. function in Table 3 represents the value of the single-objective function of the problem, as defined in (6.1), which is used by the firefly algorithm in order to find the true Pareto optimal solution. The values of the fuel cost and emission objectives have been calculated directly from (6.2) and (6.3) without using weights.

According to the analysis of the simulation results, the firefly algorithm generated a good quality of Pareto set of solutions in contrast to the goal attainment SQP method. Using this algorithm, we can find the global minimum of the problem in about 10 iterations for a population of 12 fireflies. Thus, the total number of function evaluations is just 120. This is indeed very efficient and robust.

It is also worth to mention that in Table 3, the amount of transmission losses PLoss is very low, and as a result, it does not seem to affect the finding of the Pareto optimal solution using the proposed firefly algorithm since both the Pareto optimal solutions and the fuel cost and emission objectives remain the same. The small increase in the objective function of the problem (Obj. function) is expected as we have penalized infeasible solutions, based on the problem constraints, into the single-objective function of the problem. It is notable, however, that the transmission losses slightly affect the Pareto optimal solution in the goal attainment method.

In the goal attainment SQP method, we gave the following input parameters: 𝑋0=[0.1,0.1,0.1,0.1,0.1,0.1] as initial values of power generators, weight = [0.5 0.5] as a weight of the two objectives, and goal = [603 0.2] as an initial solution of the two objectives for the algorithm. We should consider that the goals are usually found empirically by first solving every single objective of the optimization problem individually, which constitutes higher initial values than those given here. However, in this case, we give as goals the solutions we obtained by the firefly algorithm as an effective way to test the robustness of these solutions and whether these could be further improved. It is also notable that the goal attainment method gave positive attainment factor (𝛾=0.0644 and 0.0645 for transmission losses) which implies that the goal defined is unattainable, and as a result, no further improvement of the solution is possible. Therefore, the solution obtained by the firefly algorithm is both robust and efficient.

A graphical representation of the convergence of the objective function of the problem (with equal weights 𝑤1=𝑤2=𝑎=0.5) given in (6.1) with time (iterations) is shown above in Figure 1. Figures 2 and 3 show the convergence of fuel cost and emission objectives when optimized and their convergence through time (iterations) respectively. The red point in Figure 2 represents the Pareto optimal solution found by the firefly algorithm in Table 3. Moreover, in Figure 3 we can see that the decrease of fuel cost causes the increase of emission and vice versa, since these are conflicting objectives and cannot be minimized simultaneously. The aim of this optimization problem is to find a good balance between them by converting the biobjective function into a single-objective function with different weights for each objective.

6.3.2. Comparison with Other Optimization Algorithms

Recently, Abd El-Wahed et al. [10] used an efficient hybrid ant colony optimization and modified simplex method in order to solve the economic emissions load dispatch problem with random weighting factors [10]. The test system they used is the same as those described in this paper, with the same emission/cost coefficients as defined in Tables 1 and 2 and total load demand 𝐷=2.8 MW. The best Pareto optimal solution found is at 𝑥 (0.135, 0.411, 0.467, 1.00, 0.462, 0.357), which corresponds to 𝑓1891.1966 and 𝑓20.2192. This means that the lowest fuel cost is about 891.1966 $/hr, while the emissions are 0.2192 ton/hr. The transmission losses are PLoss = 0.0035. However, it is notable a minor mistake which gives erroneously that the lowest fuel cost is about 603.3 $/hr instead of 891.1966 $/hr which is the correct value. The fuel cost and emission values can be easily verified if we apply the given solutions to (6.2) and (6.3) or if we use the same equations given in [10]. Applying the firefly algorithm to the same problem, we have found a slightly better solution but with higher power losses. The Pareto optimal solution we obtained with 12 fireflies after 50 iterations is at 𝑥 (0.050000, 0.511433, 0.215548, 0.935979, 0.581323, 0.507463), which corresponds to 𝑓1 881.533506 and 𝑓20.222158 and transmission losses PLoss = 0.007614. This clearly shows that the firefly algorithm is very efficient and effective. However, obviously, further applications are highly needed to see how it may behave for solving various similar tough engineering optimization problems, such as the economic, rapid and environmental power dispatch problem [13]. The following three Figures show the stable convergence of the objective function of the problem given in (6.1) with time (Figure 4), the convergence of fuel cost and emission objectives when optimized (Figure 5) and their quick convergence through time (Figure 6). Notice again the red point in Figure 5, which represents the Pareto optimal solution found and the fact that the two objectives are conflicting and they cannot be minimized simultaneously (Figure 6).

In another research article, Kumar et al. [5] used an improved version of particle swarm optimization algorithm in order to solve the economic emissions load dispatch problem for a test system of 6 power generators, for several load demands, with transmission losses, using a price penalty factor [5]. The equations used in that article for the problem formulation constitute a slight variation of those defined in this paper.

In order to apply the firefly algorithm to this example problem, we have slightly changed the firefly algorithm. In particular, there is a change of the generation of the initial solution since the difference of the ranges of each solution (power generator) is relatively large. For the case of load demand equal to 900 MW, the initial solution is generated using the same equation (6.5), while in the other two cases the initial solution is created based on the following two equations: for𝐷=500MW;𝑥𝑗=randup-rangelow-range4+low-range,(6.7)for𝐷=700MW;𝑥𝑗=rand(up-rangelow-range)2+(up-range)3,(6.8) where 𝑥𝑗 is the new solution of 𝑗th firefly that is created, rand is a random number generator uniformly distributed in the space [0,1], while up range and low range are the upper range and lower range of the 𝑗th firefly (power generator), respectively.

The results of that research paper including the generation cost, emission level, and power losses can be summarized in Table 4, while Table 5 presents the results obtained when we applied the firefly algorithm to the same test system.

From Tables 4, 5 (next page) we can see that the firefly algorithm gives good quality of Pareto optimal solutions for the given test problem and in some cases it outperforms the improved version of particle swarm optimization. In particular, in the case of 𝐷=500 MW the Firefly algorithm obtains lower fuel cost, emissions and power losses than those found by the particle swarm optimization. In the other two cases the firefly algorithm managed to obtain a lower value in one of the two objectives, while in all three cases the algorithm achieved lower power losses.

It is also notable that the firefly algorithm solves this test problem by converting the multiobjective problem to a single-objective problem by a linear combination of different objectives as a weighted sum, while the particle swarm optimization introduces a price penalty factor for the same purpose. Moreover, by using a population of solutions (fireflies) in its search, multiple Pareto optimal solutions can be found more quickly, even in one run, in contrast to particle swarm optimization algorithm where each agent (particle) corresponds to one single solution of the problem. Finally, it is important to point out that the firefly algorithm converges in an acceptable time, which for this test system was approximately 3 seconds.

In general, the analysis of the experimental results has demonstrated that the firefly algorithm performs better than other methods used for the same problem, or at least it obtains good quality Pareto optimal solutions in significantly low computing times. It is characterized by a stable and fast convergence compared to other conventional methods and good computation efficiency, as it has been demonstrated by its application. This much improved speed of computation allows for additional searches and improvements that could be made in order to increase the confidence and efficiency of the generated solutions. As we have seen, one such improvement is the generation of the initial solution. Overall, the firefly algorithm is shown to be very helpful and promising in solving optimization problems in power systems although more example problems are necessary in order to prove its efficiency and superiority towards other alternative optimization algorithms used in the same domain, such as the particle swarm optimization and ant colony optimization algorithms.

7. Conclusions and Future Work

There is no doubt that the firefly algorithm, developed by Dr. Xin-She Yang in 2007, is a very powerful novel population-based method for solving constrained optimization problems and particularly NP-hard problems. The idea behind this algorithm is that the social behavior and especially the flashing light of fireflies can be easily formulated and associated with the objective function of a given optimization problem [1, 3, 20].

In this paper, we have proposed, presented, and tested the recently developed Firefly algorithm for application to the multiobjective minimization problem of economic emissions load dispatch. Focus was on the reduction of a single pollutant nitrogen oxide, NOx, while the transmission losses, the equality constraint of power balance, and the inequality generator capacity constraints are also considered. The results of the implementation of this proposed method clearly showed the efficiency and effectiveness of the firefly algorithm for solving the particular optimization problem and finding the true Pareto optimal solutions. The algorithm achieved good results comparable to those achieved by other stochastic nature-inspired algorithms as these reported in the literature.

For various load demands of two different test systems consisting of 6 power generators the simulation results of the algorithm compared favorably with other competitive state-of-the-art algorithms, and as a result, we can say that the firefly algorithm is very efficient and accurate in obtaining global optima with high success rates for the given constrained optimization problem. In particular, the results indicate that the algorithm obtained good quality Pareto optimal solutions which are better or at least equal to the solutions obtained by other stochastic alternative optimization algorithms such as the goal attainment SQP method, the particle swart optimization, and the Genetic algorithms, in terms of efficiency and success rates.

In all cases, the algorithm converged to the optimal solution very quickly (mainly from the 10th iteration), and its computation time in the simulations we run was less than 3 seconds. Moreover, the fact that the proposed algorithm is very simple in concept and easy to implement clearly implies that it could also be effectively applied to other multiobjective optimization problems and especially to NP-hard combinatorial optimization problems, such as the vehicle-routing problem.

However, from the simulation results, it seems that the proper selection of the population size, the number of generations (iterations), and the absorption coefficient is of paramount importance for the convergence of the algorithm as this heavily depends on the nature of the applied problem. Moreover, a refinement and improvement of the generation of the initial solution and the decreasing random step size seem very promising and beneficial to further increase the algorithm’s performance, while it might also be possible to hybridize the algorithm together with other heuristic search methods for better results, such as the Levy flight firefly algorithm [4]. Furthermore, future comparison studies with other methods are necessary so as to identify the strengths and weaknesses of the current metaheuristic algorithm and prove its superiority and effectiveness. In future works, we intend to analyze the behavior of the proposed algorithm in other variations of the economic emissions load dispatch problem, such as the environmental power load dispatch problem with more security constraints [13]. Ultimately, even better optimazation algorithms may emerge as a promising way to solve NP-hard optimization problems.

We hope that this paper will provide the necessary background for further improvement of optimization algorithms and will be the milestone conducive to further research and development in the domain of analyzing and applying evolutionary optimization algorithms to NP-hard combinatorial optimization problems.

Input: 𝛼 , 𝛾 , 𝛽 0 , 𝑛 , MaxGeneration, 𝐵 , cost-emission-capacities coefficients
Output: 𝑃 G i for 𝑖 = 1 , 6 , 𝑓 ( 𝑋 ) , 𝑓 1 ( 𝑋 ) , 𝑓 2 ( 𝑋 )
Begin of algorithm
Define the multiobjective function: max 𝑓 ( 𝑃 G i ) , with 𝑖 = 1 , , 6
Generate initial population of fireflies 𝑛 = 1 , , 1 2   %generate 𝑛 = 1 2 initial solutions
Light Intensity of firefly 𝑛 is determined by objective function, 𝐼 𝑛 𝑓 ( 𝑃 G i )
Define 𝛼 = 0 . 2 , 𝛽 0 = 1 . 0 and 𝛾 = 1 . 0   %necessary algorithm’s parameters
While ( 𝑡 ≤ MaxGeneration = 50)
 For 𝑖 = 1 1 2   %for all fireflies (solutions)
  For 𝑗 = 1 1 2   %for all fireflies (solutions)
   If ( 𝐼 𝑖 < 𝐼 𝑗 )
   Then move firefly 𝑖 towards firefly 𝑗 (move towards brighter one)
   Attractiveness varies with distance 𝑟 𝑖 𝑗 via exp ( 𝛾 𝑟 𝑖 𝑗 )
   Generate and evaluate new solutions and update Light Intensity
  End for 𝑗 loop
 End for 𝑖 loop
Check the ranges of the given solutions and update them as appropriate
Rank the fireflies, find and display the current best %max solution for each iteration
End of while loop
% Post-process results and visualization
Find the firefly with the highest Light Intensity among all fireflies %optimal solution
Plot the increase of the Light Intensity with timeiterations
Plot the two objectives with time %best solution with time
End of algorithm

Repeat for each firefly 𝑖
If 𝑃 𝐺 𝑖 < 𝑃 𝑖 m i n
 Then 𝑃 𝐺 𝑖 = 𝑃 𝑖 m i n
If 𝑃 𝐺 𝑖 > 𝑃 𝑖 m a x
 Then 𝑃 𝐺 𝑖 = 𝑃 𝑖 m a x
End of Repeat